class LatSync_McRobie(): """ Class used to implement the multiple modes + multiple TMDs analysis method used to consider the phenomenon of Lateral Synchronous crowd loading response, per the following [paper](../references/The Lateral Dynamic Stablity of Stockton Infinity Footbridge using Complex Modes.pdf) """ def __init__(self, modalsys, cp_func=None, mp_func=None, store_mtrxs=False, makePlots=False, verbose=True, **kwargs): """ Initialise analysis *** Required: * `modalsys`, instance of `ModalSys` class, defining the modal system to be analysed. For situations where TMDs are present, it is assumed that TMDs have already been appended to the parent modal system, e.g. by use of the `append_TMDs()` method in `damper.py`. *** Optional: * `cp_func`, function of form g(fn), where fn is frequency in Hz, to define how _negative_ damping per pedestrian varies with frequency """ if verbose: print("Initialising 'LatSync_McRobie' analysis...") self.modalsys = modalsys """ Modal system that form the basis of this analysis. May have subsystems pre-appended (e.g. TMDs) """ if verbose: print("Defining pedestrian added damping function...") if cp_func is None: if verbose: print("Default cp(f) function will be used") cp_func = self.init_cp_func() elif isinstance(cp_func,float) or isinstance(cp_func,int): if verbose: print("cp = %.1f will be used for all f" % cp_func) cp_val = cp_func cp_func = scipy.interpolate.interp1d([0,2.0],[cp_val,cp_val], bounds_error=False, fill_value=cp_val) elif isfunction(cp_func): if verbose: print("User-supplied cp(f) function will be used") else: raise ValueError("Unexpected `cp_func` argument") self.cp_func = cp_func """ Function to define how added damping per pedestrian varies with frequency. Function with single argument 'fn' expected, i.e. cp = g(fn), cp in [Ns/m] """ if verbose: print("Defining pedestrian added mass function...") if mp_func is None: if verbose: print("Default mp(f) function will be used") mp_func = self.init_mp_func() elif isinstance(mp_func,float): if verbose: print("mp = %.1f will be used for all f" % mp_func) mp_val = mp_func mp_func = scipy.interpolate.interp1d([0,2.0],[mp_val,mp_val], bounds_error=False, fill_value=mp_val) elif isfunction(mp_func): if verbose: print("User-supplied mp(f) function will be used") else: raise ValueError("Unexpected `mp_func` argument") self.mp_func = mp_func """ Function to define how added mass per pedestrian varies with frequency. Function with single argument 'fn' expected, i.e. mp = g(fn), cp in [kg] """ self.store_mtrxs = store_mtrxs """ Boolean option, if True then system damping and stiffness matrices will be stored via nested list for all Np values and modes considered (this can result in a large amount of data being stored!) """ # Class variables defined later on self.eigenvalues = None """ Ndarray of complex-valued eigenvalues of shape (N_Np,N_modes), where: * _N_Np_ is number of pedestrian cases considered * _N_modes_ is number of system matrix eigenmodes """ self.damping_ratios = None """ Ndarray of effective damping ratio for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.damped_freqs = None """ Ndarray of damped natural frequencies for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.eigenvectors = None """ Ndarray of complex-valued eigenvalues of shape _(N_Np,N_modes,N_dof)_, where: * _N_Np_ is number of pedestrian cases considered * _N_modes_ is number of system matrix eigenmodes * _N_dof_ is the total number of degrees of freedom / state variables for the system analysed """ self.N_pedestrians = None """ Array listing the set of pedestrian numbers analysed """ self.cp_vals = None """ Ndarray of effective damping rate per pedestrian for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.mp_vals = None """ Ndarray of effective added mass per pedestrian for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.ped_effect_mtrx = None """ Matrix to represent the smeared effect of one pedestrian on the deck Used to calculate additive system matrices due to pedestrians """ self.Np_crit = None """ Critical number of pedestrians to give zero net damping """ self.crit_mode = None """ Index of (complex) eigenmode that loses stability at Np=Np_crit """ self.fd_crit = None """ Damped natural frequency of mode that loses stability at Np=Np_crit """ self.eta_crit = None """ Effective damping ratio of mode that loses stability at Np=Np_crit (damping ratio = 0.00 expected) """ self.s_crit = None """ Complex-valued pole / eigenvalue at Np=Np_crit """ self.X_crit = None """ Complex-valued eigenvector of mode that loses stability at Np=Np_crit """ self.C_mtrx_crit = None """ Damping matrix (of 'parent' modal system) applicable to mode that loses stability at Np=Np_crit. Will in general include off-diagonal terms due to mode coupling effect of pedestrians, determined per method presented in McRobie's paper """ self.M_mtrx_crit = None """ Mass matrix (of 'parent' modal system) applicable to mode that loses stability at Np=Np_crit. Will in general include off-diagonal terms due to mode coupling effect of pedestrians, determined per method presented in McRobie's paper """ if makePlots: self.plot_cp_func() def init_cp_func(self,conservative=False): """ Initialise default function to use to define frequency-dependence of effective added damping per pedestrian (Ns/m per pedestrian). Note damping taken in conventional sense, such that positive values define added damping i.e. increased attenuation of vibrations. Correspondingly the effect of pedestrians will tend to be represented by a negative damping rate. This default profile has been defined with regards to the work done by Dallard et al (2006) and John MacDonald (2008). The default curve defined within this function is taken to be a conservative representation of the damping values presented in Figure 9(a) of the latter's paper _"Lateral excitation of bridges by balancing pedestrians"_ If `conservative=True` is used then the negative 300Ns/m value proposed by Dallard et al will be used. This value is taken to be independent of frequency. """ if conservative: f_vals = [0.0,2.0] cp_vals = [-300,-300] else: f_vals = [0.0,0.5,1.0,1.5,2.0] cp_vals = [0,-300,-300,0,0] cp_func = scipy.interpolate.interp1d(f_vals,cp_vals,kind='linear', bounds_error=False, fill_value=0.0) return cp_func def init_mp_func(self): """ Initialise default function to use to define frequency-dependence of effective added mass per pedestrian (kg per pedestrian) The default curve defined within this function is taken to be a simplified representation of the added mass values presented in Figure 9(b) of John MacDonald's 2008 Proc. of the Royal Society paper _"Lateral excitation of bridges by balancing pedestrians"_ """ f_vals = [0.0,0.5,0.8,1.25,2.0] deltaM_vals = [+70,-45,-45,0,0] deltaM_func = scipy.interpolate.interp1d(f_vals,deltaM_vals, kind='linear', bounds_error=False, fill_value=0.0) return deltaM_func def calc_pedestrian_effect(self,f=1.0,num=1000): """ Calculates modal damping matrix due to pedestrians, per the method presented in McRobie paper """ # Get system matrix to represent distributed action of pedestrians # along deck modeshapes ped_effect_mtrx = self.ped_effect_mtrx if ped_effect_mtrx is None: modalsys = self.modalsys modeshape_func = modalsys.modeshapeFunc L = modalsys.Ltrack # Evaluate modeshape ordinates at unif dL = L/num x = numpy.linspace(0,L,num,endpoint=False) + dL/2 phi = modeshape_func(x) # Evaluate matrix product of modeshape ordinate matrices phi_product = phi.T @ phi # Return mode-generalised damping matrix ped_effect_mtrx = dL / L * phi_product self.ped_effect_mtrx = ped_effect_mtrx # Evaluate cp and mp given mode natural frequency provided f = numpy.abs(f) # positive frequencies to be used cp = self.cp_func(f) mp = self.mp_func(f) C_pa = cp * ped_effect_mtrx M_pa = mp * ped_effect_mtrx # Return results as dict rslts = {} rslts['C_pa']=C_pa rslts['M_pa']=M_pa rslts['cp']=cp rslts['mp']=mp return rslts def run(self,Np_vals,verbose=True, calc_Np_crit=True,append_rslts=False, **kwargs): """ Run analysis to explore eigenvalues of system state matrix for various assumed pedestrian crowd densities *** Required: * `Np_vals`, _array-like_: each value defines total number of pedestrians on bridge; a range of values will typically be provided, to allow exploration of how system eigenproperties (in particular effective damping) varies with pedestrian numbers """ if verbose: print("Running lat sync eigenvalues analysis...") # Run analysis self._run_analysis(Np_vals,append_rslts,**kwargs) # Calculate critical number of pedestrians for onset on instability if calc_Np_crit: Np_crit = self.calc_Np_crit(verbose=verbose) else: Np_crit = None if verbose: print("Analysis complete!") return Np_crit def _run_analysis(self,Np_vals, append_rslts=True, **kwargs): """ Run analysis for various pedestrian numbers, as provided to function To cater for frequency-dependent nature of pedestrian damping/mass effect an iterative procedure is adopted, similar to the p-k method devised by Theodorsen for solving aeroelastic problems """ modalsys = self.modalsys store_mtrxs = self.store_mtrxs # Take copy of system damping matrix with no pedestrians C0 = deepcopy(modalsys._C_mtrx) M0 = deepcopy(modalsys._M_mtrx) # Define function to iterate with def calc_modal_properties(f,mode_index,Np,return_rslts=False): # Calculate change in model damping matrix # due to smeared effect of N=1 pedestrian rslts = self.calc_pedestrian_effect(f=f,**kwargs) C_pa = rslts['C_pa'] M_pa = rslts['M_pa'] # Adjust bridge damping and mass matrices modalsys._C_mtrx = C0 + Np * C_pa modalsys._M_mtrx = M0 + Np * M_pa # Carry out eigevalue analysis using updated system matrices eig_props = modalsys.CalcEigenproperties() f_new = eig_props['f_d'][mode_index] f_error = f_new - f #print("Calculating: Np=%d, mode=%d, f=%.3f, f_error=%.3f" # % (Np,mode_index,f,f_error)) if return_rslts: d = {} d['eig_props']=eig_props d = {**d, **rslts} d['C_mtrx'] = modalsys._C_mtrx d['M_mtrx'] = modalsys._M_mtrx return f_error, d else: return f_error # Carry out eigevalue analysis of system with no pedestrians eig_props = modalsys.CalcEigenproperties() fd_vals_last = eig_props['f_d'] # Loop over all pedestrian numbers provided to function cp_vals = [] mp_vals = [] s_vals = [] eta_vals = [] fd_vals = [] X_vals = [] if store_mtrxs: C_mtrx_list = [] M_mtrx_list = [] for _Np in Np_vals: cp_vals_inner = [] mp_vals_inner = [] s_vals_inner = [] eta_vals_inner = [] fd_vals_inner = [] X_vals_inner = [] if store_mtrxs: C_mtrx_inner_list = [] M_mtrx_inner_list = [] for _mode_index, _fd in enumerate(fd_vals_last): # Solve for consistent frequencies, allowing for # frequency-dependence of pedestrian-related mass and damping fd_sol = scipy.optimize.newton(func=calc_modal_properties, x0=_fd, args=(_mode_index,_Np)) # Rerun for converged frequency of given mode ferr, rslts = calc_modal_properties(fd_sol, _mode_index,_Np, return_rslts=True) eig_props = rslts['eig_props'] cp = rslts['cp'] mp = rslts['mp'] C_mtrx = rslts['C_mtrx'] M_mtrx = rslts['M_mtrx'] # Unpack results for this mode and append to inner lists cp_vals_inner.append(cp) mp_vals_inner.append(mp) s_vals_inner.append(eig_props['s'][_mode_index]) eta_vals_inner.append(eig_props['eta'][_mode_index]) fd_vals_inner.append(eig_props['f_d'][_mode_index]) X_vals_inner.append(numpy.ravel(eig_props['X'][:,_mode_index])) if store_mtrxs: C_mtrx_inner_list.append(C_mtrx) M_mtrx_inner_list.append(M_mtrx) # Update last frequencies fd_vals_last = fd_vals_inner # Append inner lists to outer lists cp_vals.append(cp_vals_inner) mp_vals.append(mp_vals_inner) s_vals.append(s_vals_inner) eta_vals.append(eta_vals_inner) fd_vals.append(fd_vals_inner) X_vals.append(X_vals_inner) if store_mtrxs: C_mtrx_list.append(C_mtrx_inner_list) M_mtrx_list.append(M_mtrx_inner_list) # Convert nested lists to numpy ndarray type cp_vals = numpy.array(cp_vals) mp_vals = numpy.array(mp_vals) s_vals = numpy.array(s_vals) eta_vals = numpy.array(eta_vals) fd_vals = numpy.array(fd_vals) X_vals = numpy.array(X_vals) # Restore original no-pedestrians system matrices modalsys._C_mtrx = C0 modalsys._M_mtrx = M0 # Check to see if previous results exist if not hasattr(self,'eigenvalues'): append_rslts = False # no previous results avaliable # Record key results as attributes, or append to previous results if not append_rslts: self.eigenvalues = s_vals self.damping_ratios = eta_vals self.damped_freqs = fd_vals self.eigenvectors = X_vals self.N_pedestrians = Np_vals self.cp_vals = cp_vals self.mp_vals = mp_vals if store_mtrxs: self.C_mtrx_list = C_mtrx_list self.M_mtrx_list = M_mtrx_list if append_rslts: self.eigenvalues = numpy.vstack((self.eigenvalues,s_vals)) self.damping_ratios = numpy.vstack((self.damping_ratios,eta_vals)) self.damped_freqs = numpy.vstack((self.damped_freqs,fd_vals)) self.eigenvectors = numpy.vstack((self.eigenvectors,X_vals)) self.N_pedestrians = numpy.hstack((self.N_pedestrians,Np_vals)) self.cp_vals = numpy.vstack((self.cp_vals,cp_vals)) self.mp_vals = numpy.vstack((self.mp_vals,mp_vals)) if store_mtrxs: self.C_mtrx_list.append(C_mtrx_list) self.M_mtrx_list.append(M_mtrx_list) def plot_results(self): """ Plots results from the above analysis N.b: Figures intentionally emulate the format of figures in McRobie's paper, as these have been used to validate the above routine """ fig = plt.figure() fig.set_size_inches((14,9)) gs = gridspec.GridSpec(3, 2, height_ratios=[1.5, 1, 1]) gs.update(hspace=0.7) # Prepare subplots ax = fig.add_subplot(gs[0, 0]) self.plot_damping_vs_freq(ax=ax) ax = fig.add_subplot(gs[0, 1]) self.plot_poles(ax=ax) ax = fig.add_subplot(gs[1:, 0]) self.plot_damping_vs_pedestrians(ax=ax) ax1 = fig.add_subplot(gs[1, 1]) self.plot_cp_func(ax=ax1) ax = fig.add_subplot(gs[2, 1],sharex = ax1) self.plot_mp_func(ax=ax) return fig def plot_damping_vs_freq(self,ax=None): """ Plot damping ratio against damped natural frequency (per Figure 5 in McRobie's paper) """ if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() eta_vals = self.damping_ratios fd_vals = self.damped_freqs ax.plot(eta_vals,fd_vals,'k.',markersize=0.5) ax.plot(eta_vals[::10],fd_vals[::10],'ko',markersize=1.5) ax.plot(eta_vals[0],fd_vals[0],'bo',markersize=3.0) ax.axvline(x=0.0,color='r',alpha=0.3) # denotes stability limit ax.set_ylim([0.0,ax.get_ylim()[1]]) ax.set_xlabel("Damping ratio") ax.set_ylabel("Damped natural frequency (Hz)") ax.set_title("Frequency vs Effective Damping\n") if self.eta_crit is not None: ax.plot(self.eta_crit,self.fd_crit,'r.') return fig def plot_poles(self,ax=None): """ Plot eigenvalues (poles) on complex plane (per Figure 4 in paper) """ if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() s_vals = self.eigenvalues ax.plot(real(s_vals),imag(s_vals),'k.',markersize=0.5) ax.plot(real(s_vals[::10]),imag(s_vals[::10]),'ko',markersize=1.5) ax.plot(real(s_vals[0]),imag(s_vals[0]),'bo',markersize=3.0) ax.axvline(x=0.0,color='r',alpha=0.3) # stability limit ax.axhline(y=0.0,color='k',linewidth=0.5) # illustrates Im(z)=0 axis ax.set_xlabel("Real(s)") ax.set_ylabel("Imag(s)") ax.set_title("Eigenvalues of system state matrix") s_crit = self.s_crit if s_crit is not None: ax.plot(real(s_crit),imag(s_crit),'r.') return fig def plot_damping_vs_pedestrians(self,ax=None): """ Plots effective damping ratio of poles against number of pedestrians """ Np_vals = self.N_pedestrians eta_vals = self.damping_ratios ax.plot(Np_vals,eta_vals,'k.',markersize=0.5) ax.axhline(y=0.0,color='r',alpha=0.3) ax.set_xlim([0,ax.get_xlim()[1]]) ax.set_xlabel("Number of pedestrians") ax.set_ylabel("Effective damping ratio") ax.set_title("Effect of pedestrians on effective damping") Np_crit = self.Np_crit if Np_crit is not None: ax.axvline(x=Np_crit,color='r',alpha=0.3) ax.plot(Np_crit,self.eta_crit,'r.') def plot_cp_func(self,ax=None,f_vals=None): """ Plots function defining how negative damping per pedestrian varies with frequency """ cp_func = self.cp_func if f_vals is None: if hasattr(cp_func,'x'): f_vals = cp_func.x else: raise ValueError("`f_vals` must be provided!") cp_vals = cp_func(f_vals) if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() ax.plot(f_vals,cp_vals,label='Defined $c_p(f)$ function') ax.set_xlim([0,ax.get_xlim()[1]]) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("$c_p$ (Ns/m)") ax.set_title("Damping rate per pedestrian, $c_p$") if hasattr(self,'cp_vals'): fd_vals = self.damped_freqs cp_vals = self.cp_vals lines = ax.plot(fd_vals,cp_vals,'r|') for m, line in enumerate(lines): if m==0: line.set_label('$c_p$ values for modes') ax.legend() return fig def plot_mp_func(self,ax=None,f_vals=None): """ Plots function defining how negative damping per pedestrian varies with frequency """ mp_func = self.mp_func if f_vals is None: if hasattr(mp_func,'x'): f_vals = mp_func.x else: raise ValueError("`f_vals` must be provided!") mp_vals = mp_func(f_vals) if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() ax.plot(f_vals,mp_vals,label='Defined $m_p(f)$ function') ax.set_xlim([0,ax.get_xlim()[1]]) ax.set_ylim([-100,+100]) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("$\Delta M$ (kg)") ax.set_title("Added mass per pedestrian, $\Delta M$") if hasattr(self,'mp_vals'): fd_vals = self.damped_freqs mp_vals = self.mp_vals lines = ax.plot(fd_vals,mp_vals,'r|') for m, line in enumerate(lines): if m==0: line.set_label('$\Delta M$ values for modes') ax.legend() return fig def calc_Np_crit(self,rerun_factor=2.0,verbose=True): """ Calculates from results the critical number of pedestrians for the onset of pedestrian-induced lateral instability """ if verbose: print("Determining critical number of pedestrians " + "for onset on instability...") # Check that analysis has been run if not hasattr(self,"N_pedestrians"): raise ValueError("Analysis does not appear to have been run!" + "Cannot calculate Np_crit") # Check that analysis has been run up to pedestrian numbers such that # the net effective damping in (at least) one mode has become <0.0 damping_ratios = numpy.array(self.damping_ratios) min_damping = numpy.min(damping_ratios) while min_damping >= 0.0: if verbose: print("Analysis has not identified net damping < 0\n" + "Analysis will be extended to larger pedestrian numbers") # Get key details of Np analysed last time around Np_max = numpy.max(self.N_pedestrians) # max from last run Np_min = numpy.min(self.N_pedestrians) # min from last run Np_step = (Np_max - Np_min) / (len(self.N_pedestrians)) # Define Np_vals for new analysis new_Np_vals = numpy.arange(Np_max+1,Np_max*rerun_factor,Np_step) # Rerun analysis self._run_analysis(new_Np_vals,append_rslts=True) # Recalculate minimum damping damping_ratios = numpy.array(self.damping_ratios) min_damping = numpy.min(damping_ratios) # Determine minumum damping ratio across all modes for each Np value min_damping = numpy.min(damping_ratios,axis=1) # Define interpolation function Np_vals = self.N_pedestrians Np_func = scipy.interpolate.interp1d(x=Np_vals, y=min_damping) # Use root finding function to obtain Np such that damping = 0.0 Np_crit = scipy.optimize.bisect(f=Np_func,a=Np_vals[0],b=Np_vals[-1]) self.Np_crit = Np_crit # Rerun analysis to obtain properties at critical point # 1) take copy of class to pick up same attributes as analysis run crit_rslts = deepcopy(self) # 2) Edit a few key properties for this purpose crit_rslts.store_mtrxs = True # 3) rerun to override results already held in the above object crit_rslts.run(Np_vals=[Np_crit],calc_Np_crit=False,append_rslts=False) # Determine critical mode i.e. mode with damping =0 when Np=Np_crit crit_mode = numpy.argmin(crit_rslts.damping_ratios[0,:]) self.crit_mode = crit_mode # Extract key properties of critical mode at Np=Np_crit self.fd_crit = crit_rslts.damped_freqs[0,crit_mode] self.eta_crit = crit_rslts.damping_ratios[0,crit_mode] self.s_crit = crit_rslts.eigenvalues[0,crit_mode] self.X_crit = crit_rslts.eigenvectors[0,crit_mode,:] self.C_mtrx_crit = crit_rslts.C_mtrx_list[-1][crit_mode] self.M_mtrx_crit = crit_rslts.M_mtrx_list[-1][crit_mode] if verbose: print("Critical no. pedestrians, Np_crit:\t%.0f" % self.Np_crit) print("Index of critical mode:\t\t%d" % self.crit_mode) print("Damped nat freq, critical mode:\t%.3f Hz" % self.fd_crit) print("Damping ratio, critical mode:\t%.2e" % self.eta_crit) return Np_crit
Ancestors (in MRO)
- LatSync_McRobie
- builtins.object
Static methods
def __init__(
self, modalsys, cp_func=None, mp_func=None, store_mtrxs=False, makePlots=False, verbose=True, **kwargs)
Initialise analysis
Required:
modalsys
, instance ofModalSys
class, defining the modal system to be analysed. For situations where TMDs are present, it is assumed that TMDs have already been appended to the parent modal system, e.g. by use of theappend_TMDs()
method indamper.py
.
Optional:
cp_func
, function of form g(fn), where fn is frequency in Hz, to define how negative damping per pedestrian varies with frequency
def __init__(self, modalsys, cp_func=None, mp_func=None, store_mtrxs=False, makePlots=False, verbose=True, **kwargs): """ Initialise analysis *** Required: * `modalsys`, instance of `ModalSys` class, defining the modal system to be analysed. For situations where TMDs are present, it is assumed that TMDs have already been appended to the parent modal system, e.g. by use of the `append_TMDs()` method in `damper.py`. *** Optional: * `cp_func`, function of form g(fn), where fn is frequency in Hz, to define how _negative_ damping per pedestrian varies with frequency """ if verbose: print("Initialising 'LatSync_McRobie' analysis...") self.modalsys = modalsys """ Modal system that form the basis of this analysis. May have subsystems pre-appended (e.g. TMDs) """ if verbose: print("Defining pedestrian added damping function...") if cp_func is None: if verbose: print("Default cp(f) function will be used") cp_func = self.init_cp_func() elif isinstance(cp_func,float) or isinstance(cp_func,int): if verbose: print("cp = %.1f will be used for all f" % cp_func) cp_val = cp_func cp_func = scipy.interpolate.interp1d([0,2.0],[cp_val,cp_val], bounds_error=False, fill_value=cp_val) elif isfunction(cp_func): if verbose: print("User-supplied cp(f) function will be used") else: raise ValueError("Unexpected `cp_func` argument") self.cp_func = cp_func """ Function to define how added damping per pedestrian varies with frequency. Function with single argument 'fn' expected, i.e. cp = g(fn), cp in [Ns/m] """ if verbose: print("Defining pedestrian added mass function...") if mp_func is None: if verbose: print("Default mp(f) function will be used") mp_func = self.init_mp_func() elif isinstance(mp_func,float): if verbose: print("mp = %.1f will be used for all f" % mp_func) mp_val = mp_func mp_func = scipy.interpolate.interp1d([0,2.0],[mp_val,mp_val], bounds_error=False, fill_value=mp_val) elif isfunction(mp_func): if verbose: print("User-supplied mp(f) function will be used") else: raise ValueError("Unexpected `mp_func` argument") self.mp_func = mp_func """ Function to define how added mass per pedestrian varies with frequency. Function with single argument 'fn' expected, i.e. mp = g(fn), cp in [kg] """ self.store_mtrxs = store_mtrxs """ Boolean option, if True then system damping and stiffness matrices will be stored via nested list for all Np values and modes considered (this can result in a large amount of data being stored!) """ # Class variables defined later on self.eigenvalues = None """ Ndarray of complex-valued eigenvalues of shape (N_Np,N_modes), where: * _N_Np_ is number of pedestrian cases considered * _N_modes_ is number of system matrix eigenmodes """ self.damping_ratios = None """ Ndarray of effective damping ratio for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.damped_freqs = None """ Ndarray of damped natural frequencies for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.eigenvectors = None """ Ndarray of complex-valued eigenvalues of shape _(N_Np,N_modes,N_dof)_, where: * _N_Np_ is number of pedestrian cases considered * _N_modes_ is number of system matrix eigenmodes * _N_dof_ is the total number of degrees of freedom / state variables for the system analysed """ self.N_pedestrians = None """ Array listing the set of pedestrian numbers analysed """ self.cp_vals = None """ Ndarray of effective damping rate per pedestrian for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.mp_vals = None """ Ndarray of effective added mass per pedestrian for each mode and pedestrian count considered. _Ndarray shape as per eigenvalues above._ """ self.ped_effect_mtrx = None """ Matrix to represent the smeared effect of one pedestrian on the deck Used to calculate additive system matrices due to pedestrians """ self.Np_crit = None """ Critical number of pedestrians to give zero net damping """ self.crit_mode = None """ Index of (complex) eigenmode that loses stability at Np=Np_crit """ self.fd_crit = None """ Damped natural frequency of mode that loses stability at Np=Np_crit """ self.eta_crit = None """ Effective damping ratio of mode that loses stability at Np=Np_crit (damping ratio = 0.00 expected) """ self.s_crit = None """ Complex-valued pole / eigenvalue at Np=Np_crit """ self.X_crit = None """ Complex-valued eigenvector of mode that loses stability at Np=Np_crit """ self.C_mtrx_crit = None """ Damping matrix (of 'parent' modal system) applicable to mode that loses stability at Np=Np_crit. Will in general include off-diagonal terms due to mode coupling effect of pedestrians, determined per method presented in McRobie's paper """ self.M_mtrx_crit = None """ Mass matrix (of 'parent' modal system) applicable to mode that loses stability at Np=Np_crit. Will in general include off-diagonal terms due to mode coupling effect of pedestrians, determined per method presented in McRobie's paper """ if makePlots: self.plot_cp_func()
def calc_Np_crit(
self, rerun_factor=2.0, verbose=True)
Calculates from results the critical number of pedestrians for the onset of pedestrian-induced lateral instability
def calc_Np_crit(self,rerun_factor=2.0,verbose=True): """ Calculates from results the critical number of pedestrians for the onset of pedestrian-induced lateral instability """ if verbose: print("Determining critical number of pedestrians " + "for onset on instability...") # Check that analysis has been run if not hasattr(self,"N_pedestrians"): raise ValueError("Analysis does not appear to have been run!" + "Cannot calculate Np_crit") # Check that analysis has been run up to pedestrian numbers such that # the net effective damping in (at least) one mode has become <0.0 damping_ratios = numpy.array(self.damping_ratios) min_damping = numpy.min(damping_ratios) while min_damping >= 0.0: if verbose: print("Analysis has not identified net damping < 0\n" + "Analysis will be extended to larger pedestrian numbers") # Get key details of Np analysed last time around Np_max = numpy.max(self.N_pedestrians) # max from last run Np_min = numpy.min(self.N_pedestrians) # min from last run Np_step = (Np_max - Np_min) / (len(self.N_pedestrians)) # Define Np_vals for new analysis new_Np_vals = numpy.arange(Np_max+1,Np_max*rerun_factor,Np_step) # Rerun analysis self._run_analysis(new_Np_vals,append_rslts=True) # Recalculate minimum damping damping_ratios = numpy.array(self.damping_ratios) min_damping = numpy.min(damping_ratios) # Determine minumum damping ratio across all modes for each Np value min_damping = numpy.min(damping_ratios,axis=1) # Define interpolation function Np_vals = self.N_pedestrians Np_func = scipy.interpolate.interp1d(x=Np_vals, y=min_damping) # Use root finding function to obtain Np such that damping = 0.0 Np_crit = scipy.optimize.bisect(f=Np_func,a=Np_vals[0],b=Np_vals[-1]) self.Np_crit = Np_crit # Rerun analysis to obtain properties at critical point # 1) take copy of class to pick up same attributes as analysis run crit_rslts = deepcopy(self) # 2) Edit a few key properties for this purpose crit_rslts.store_mtrxs = True # 3) rerun to override results already held in the above object crit_rslts.run(Np_vals=[Np_crit],calc_Np_crit=False,append_rslts=False) # Determine critical mode i.e. mode with damping =0 when Np=Np_crit crit_mode = numpy.argmin(crit_rslts.damping_ratios[0,:]) self.crit_mode = crit_mode # Extract key properties of critical mode at Np=Np_crit self.fd_crit = crit_rslts.damped_freqs[0,crit_mode] self.eta_crit = crit_rslts.damping_ratios[0,crit_mode] self.s_crit = crit_rslts.eigenvalues[0,crit_mode] self.X_crit = crit_rslts.eigenvectors[0,crit_mode,:] self.C_mtrx_crit = crit_rslts.C_mtrx_list[-1][crit_mode] self.M_mtrx_crit = crit_rslts.M_mtrx_list[-1][crit_mode] if verbose: print("Critical no. pedestrians, Np_crit:\t%.0f" % self.Np_crit) print("Index of critical mode:\t\t%d" % self.crit_mode) print("Damped nat freq, critical mode:\t%.3f Hz" % self.fd_crit) print("Damping ratio, critical mode:\t%.2e" % self.eta_crit) return Np_crit
def calc_pedestrian_effect(
self, f=1.0, num=1000)
Calculates modal damping matrix due to pedestrians, per the method presented in McRobie paper
def calc_pedestrian_effect(self,f=1.0,num=1000): """ Calculates modal damping matrix due to pedestrians, per the method presented in McRobie paper """ # Get system matrix to represent distributed action of pedestrians # along deck modeshapes ped_effect_mtrx = self.ped_effect_mtrx if ped_effect_mtrx is None: modalsys = self.modalsys modeshape_func = modalsys.modeshapeFunc L = modalsys.Ltrack # Evaluate modeshape ordinates at unif dL = L/num x = numpy.linspace(0,L,num,endpoint=False) + dL/2 phi = modeshape_func(x) # Evaluate matrix product of modeshape ordinate matrices phi_product = phi.T @ phi # Return mode-generalised damping matrix ped_effect_mtrx = dL / L * phi_product self.ped_effect_mtrx = ped_effect_mtrx # Evaluate cp and mp given mode natural frequency provided f = numpy.abs(f) # positive frequencies to be used cp = self.cp_func(f) mp = self.mp_func(f) C_pa = cp * ped_effect_mtrx M_pa = mp * ped_effect_mtrx # Return results as dict rslts = {} rslts['C_pa']=C_pa rslts['M_pa']=M_pa rslts['cp']=cp rslts['mp']=mp return rslts
def init_cp_func(
self, conservative=False)
Initialise default function to use to define frequency-dependence of effective added damping per pedestrian (Ns/m per pedestrian).
Note damping taken in conventional sense, such that positive values define added damping i.e. increased attenuation of vibrations. Correspondingly the effect of pedestrians will tend to be represented by a negative damping rate.
This default profile has been defined with regards to the work done by Dallard et al (2006) and John MacDonald (2008). The default curve defined within this function is taken to be a conservative representation of the damping values presented in Figure 9(a) of the latter's paper "Lateral excitation of bridges by balancing pedestrians"
If conservative=True
is used then the negative 300Ns/m value proposed
by Dallard et al will be used. This value is taken to be independent
of frequency.
def init_cp_func(self,conservative=False): """ Initialise default function to use to define frequency-dependence of effective added damping per pedestrian (Ns/m per pedestrian). Note damping taken in conventional sense, such that positive values define added damping i.e. increased attenuation of vibrations. Correspondingly the effect of pedestrians will tend to be represented by a negative damping rate. This default profile has been defined with regards to the work done by Dallard et al (2006) and John MacDonald (2008). The default curve defined within this function is taken to be a conservative representation of the damping values presented in Figure 9(a) of the latter's paper _"Lateral excitation of bridges by balancing pedestrians"_ If `conservative=True` is used then the negative 300Ns/m value proposed by Dallard et al will be used. This value is taken to be independent of frequency. """ if conservative: f_vals = [0.0,2.0] cp_vals = [-300,-300] else: f_vals = [0.0,0.5,1.0,1.5,2.0] cp_vals = [0,-300,-300,0,0] cp_func = scipy.interpolate.interp1d(f_vals,cp_vals,kind='linear', bounds_error=False, fill_value=0.0) return cp_func
def init_mp_func(
self)
Initialise default function to use to define frequency-dependence of effective added mass per pedestrian (kg per pedestrian)
The default curve defined within this function is taken to be a simplified representation of the added mass values presented in Figure 9(b) of John MacDonald's 2008 Proc. of the Royal Society paper "Lateral excitation of bridges by balancing pedestrians"
def init_mp_func(self): """ Initialise default function to use to define frequency-dependence of effective added mass per pedestrian (kg per pedestrian) The default curve defined within this function is taken to be a simplified representation of the added mass values presented in Figure 9(b) of John MacDonald's 2008 Proc. of the Royal Society paper _"Lateral excitation of bridges by balancing pedestrians"_ """ f_vals = [0.0,0.5,0.8,1.25,2.0] deltaM_vals = [+70,-45,-45,0,0] deltaM_func = scipy.interpolate.interp1d(f_vals,deltaM_vals, kind='linear', bounds_error=False, fill_value=0.0) return deltaM_func
def plot_cp_func(
self, ax=None, f_vals=None)
Plots function defining how negative damping per pedestrian varies with frequency
def plot_cp_func(self,ax=None,f_vals=None): """ Plots function defining how negative damping per pedestrian varies with frequency """ cp_func = self.cp_func if f_vals is None: if hasattr(cp_func,'x'): f_vals = cp_func.x else: raise ValueError("`f_vals` must be provided!") cp_vals = cp_func(f_vals) if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() ax.plot(f_vals,cp_vals,label='Defined $c_p(f)$ function') ax.set_xlim([0,ax.get_xlim()[1]]) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("$c_p$ (Ns/m)") ax.set_title("Damping rate per pedestrian, $c_p$") if hasattr(self,'cp_vals'): fd_vals = self.damped_freqs cp_vals = self.cp_vals lines = ax.plot(fd_vals,cp_vals,'r|') for m, line in enumerate(lines): if m==0: line.set_label('$c_p$ values for modes') ax.legend() return fig
def plot_damping_vs_freq(
self, ax=None)
Plot damping ratio against damped natural frequency (per Figure 5 in McRobie's paper)
def plot_damping_vs_freq(self,ax=None): """ Plot damping ratio against damped natural frequency (per Figure 5 in McRobie's paper) """ if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() eta_vals = self.damping_ratios fd_vals = self.damped_freqs ax.plot(eta_vals,fd_vals,'k.',markersize=0.5) ax.plot(eta_vals[::10],fd_vals[::10],'ko',markersize=1.5) ax.plot(eta_vals[0],fd_vals[0],'bo',markersize=3.0) ax.axvline(x=0.0,color='r',alpha=0.3) # denotes stability limit ax.set_ylim([0.0,ax.get_ylim()[1]]) ax.set_xlabel("Damping ratio") ax.set_ylabel("Damped natural frequency (Hz)") ax.set_title("Frequency vs Effective Damping\n") if self.eta_crit is not None: ax.plot(self.eta_crit,self.fd_crit,'r.') return fig
def plot_damping_vs_pedestrians(
self, ax=None)
Plots effective damping ratio of poles against number of pedestrians
def plot_damping_vs_pedestrians(self,ax=None): """ Plots effective damping ratio of poles against number of pedestrians """ Np_vals = self.N_pedestrians eta_vals = self.damping_ratios ax.plot(Np_vals,eta_vals,'k.',markersize=0.5) ax.axhline(y=0.0,color='r',alpha=0.3) ax.set_xlim([0,ax.get_xlim()[1]]) ax.set_xlabel("Number of pedestrians") ax.set_ylabel("Effective damping ratio") ax.set_title("Effect of pedestrians on effective damping") Np_crit = self.Np_crit if Np_crit is not None: ax.axvline(x=Np_crit,color='r',alpha=0.3) ax.plot(Np_crit,self.eta_crit,'r.')
def plot_mp_func(
self, ax=None, f_vals=None)
Plots function defining how negative damping per pedestrian varies with frequency
def plot_mp_func(self,ax=None,f_vals=None): """ Plots function defining how negative damping per pedestrian varies with frequency """ mp_func = self.mp_func if f_vals is None: if hasattr(mp_func,'x'): f_vals = mp_func.x else: raise ValueError("`f_vals` must be provided!") mp_vals = mp_func(f_vals) if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() ax.plot(f_vals,mp_vals,label='Defined $m_p(f)$ function') ax.set_xlim([0,ax.get_xlim()[1]]) ax.set_ylim([-100,+100]) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("$\Delta M$ (kg)") ax.set_title("Added mass per pedestrian, $\Delta M$") if hasattr(self,'mp_vals'): fd_vals = self.damped_freqs mp_vals = self.mp_vals lines = ax.plot(fd_vals,mp_vals,'r|') for m, line in enumerate(lines): if m==0: line.set_label('$\Delta M$ values for modes') ax.legend() return fig
def plot_poles(
self, ax=None)
Plot eigenvalues (poles) on complex plane (per Figure 4 in paper)
def plot_poles(self,ax=None): """ Plot eigenvalues (poles) on complex plane (per Figure 4 in paper) """ if ax is None: fig, ax = plt.subplots() else: fig = ax.get_figure() s_vals = self.eigenvalues ax.plot(real(s_vals),imag(s_vals),'k.',markersize=0.5) ax.plot(real(s_vals[::10]),imag(s_vals[::10]),'ko',markersize=1.5) ax.plot(real(s_vals[0]),imag(s_vals[0]),'bo',markersize=3.0) ax.axvline(x=0.0,color='r',alpha=0.3) # stability limit ax.axhline(y=0.0,color='k',linewidth=0.5) # illustrates Im(z)=0 axis ax.set_xlabel("Real(s)") ax.set_ylabel("Imag(s)") ax.set_title("Eigenvalues of system state matrix") s_crit = self.s_crit if s_crit is not None: ax.plot(real(s_crit),imag(s_crit),'r.') return fig
def plot_results(
self)
Plots results from the above analysis
N.b: Figures intentionally emulate the format of figures in McRobie's paper, as these have been used to validate the above routine
def plot_results(self): """ Plots results from the above analysis N.b: Figures intentionally emulate the format of figures in McRobie's paper, as these have been used to validate the above routine """ fig = plt.figure() fig.set_size_inches((14,9)) gs = gridspec.GridSpec(3, 2, height_ratios=[1.5, 1, 1]) gs.update(hspace=0.7) # Prepare subplots ax = fig.add_subplot(gs[0, 0]) self.plot_damping_vs_freq(ax=ax) ax = fig.add_subplot(gs[0, 1]) self.plot_poles(ax=ax) ax = fig.add_subplot(gs[1:, 0]) self.plot_damping_vs_pedestrians(ax=ax) ax1 = fig.add_subplot(gs[1, 1]) self.plot_cp_func(ax=ax1) ax = fig.add_subplot(gs[2, 1],sharex = ax1) self.plot_mp_func(ax=ax) return fig
def run(
self, Np_vals, verbose=True, calc_Np_crit=True, append_rslts=False, **kwargs)
Run analysis to explore eigenvalues of system state matrix for various assumed pedestrian crowd densities
Required:
Np_vals
, array-like: each value defines total number of pedestrians on bridge; a range of values will typically be provided, to allow exploration of how system eigenproperties (in particular effective damping) varies with pedestrian numbers
def run(self,Np_vals,verbose=True, calc_Np_crit=True,append_rslts=False, **kwargs): """ Run analysis to explore eigenvalues of system state matrix for various assumed pedestrian crowd densities *** Required: * `Np_vals`, _array-like_: each value defines total number of pedestrians on bridge; a range of values will typically be provided, to allow exploration of how system eigenproperties (in particular effective damping) varies with pedestrian numbers """ if verbose: print("Running lat sync eigenvalues analysis...") # Run analysis self._run_analysis(Np_vals,append_rslts,**kwargs) # Calculate critical number of pedestrians for onset on instability if calc_Np_crit: Np_crit = self.calc_Np_crit(verbose=verbose) else: Np_crit = None if verbose: print("Analysis complete!") return Np_crit
Instance variables
var C_mtrx_crit
Damping matrix (of 'parent' modal system) applicable to mode that loses stability at Np=Np_crit.
Will in general include off-diagonal terms due to mode coupling effect of pedestrians, determined per method presented in McRobie's paper
var M_mtrx_crit
Mass matrix (of 'parent' modal system) applicable to mode that loses stability at Np=Np_crit.
Will in general include off-diagonal terms due to mode coupling effect of pedestrians, determined per method presented in McRobie's paper
var N_pedestrians
Array listing the set of pedestrian numbers analysed
var Np_crit
Critical number of pedestrians to give zero net damping
var X_crit
Complex-valued eigenvector of mode that loses stability at Np=Np_crit
var cp_func
Function to define how added damping per pedestrian varies with frequency. Function with single argument 'fn' expected, i.e. cp = g(fn), cp in [Ns/m]
var cp_vals
Ndarray of effective damping rate per pedestrian for each mode and pedestrian count considered.
Ndarray shape as per eigenvalues above.
var crit_mode
Index of (complex) eigenmode that loses stability at Np=Np_crit
var damped_freqs
Ndarray of damped natural frequencies for each mode and pedestrian count considered.
Ndarray shape as per eigenvalues above.
var damping_ratios
Ndarray of effective damping ratio for each mode and pedestrian count considered.
Ndarray shape as per eigenvalues above.
var eigenvalues
Ndarray of complex-valued eigenvalues of shape (N_Np,N_modes), where:
-
N_Np is number of pedestrian cases considered
-
N_modes is number of system matrix eigenmodes
var eigenvectors
Ndarray of complex-valued eigenvalues of shape (N_Np,N_modes,N_dof), where: * N_Np is number of pedestrian cases considered
-
N_modes is number of system matrix eigenmodes
-
N_dof is the total number of degrees of freedom / state variables for the system analysed
var eta_crit
Effective damping ratio of mode that loses stability at Np=Np_crit (damping ratio = 0.00 expected)
var fd_crit
Damped natural frequency of mode that loses stability at Np=Np_crit
var modalsys
Modal system that form the basis of this analysis. May have subsystems pre-appended (e.g. TMDs)
var mp_func
Function to define how added mass per pedestrian varies with frequency. Function with single argument 'fn' expected, i.e. mp = g(fn), cp in [kg]
var mp_vals
Ndarray of effective added mass per pedestrian for each mode and pedestrian count considered.
Ndarray shape as per eigenvalues above.
var ped_effect_mtrx
Matrix to represent the smeared effect of one pedestrian on the deck Used to calculate additive system matrices due to pedestrians
var s_crit
Complex-valued pole / eigenvalue at Np=Np_crit
var store_mtrxs
Boolean option, if True then system damping and stiffness matrices will be stored via nested list for all Np values and modes considered (this can result in a large amount of data being stored!)
class PedestrianDynamics_transientAnalyses
Implements the full set of analyses required to fully-implement the method given in NA.2.44.4 of UK NA to BS EN 1991-2 i.e. transient analyses for both walkers and joggers, for all modes
class PedestrianDynamics_transientAnalyses(dyn_analysis.Multiple): """ Implements the full set of analyses required to fully-implement the method given in NA.2.44.4 of UK NA to BS EN 1991-2 i.e. transient analyses for both walkers and joggers, for all modes """ def __init__(self, modalsys_obj, bridgeClass='A', **kwargs): """ Initialisation function **** Required: * `modalsys_obj`, modal system to which analysis relates *** Optional: * `bridgeClass`, _string_, either 'A', 'B', 'C', or 'D'. Refer Table NA.7 for description of bridge classes Additional keyword arguments may be passed. These should relate to the `__init__()` function of the `Multiple` class. """ # Get mode indexs to loop over # Note modes sorted into ascending frequency but in conjugate pairs # Hence step through modes x2 nModes = modalsys_obj.GetSystemMatrices()["nDOF"] mode_index_list = numpy.arange(0,2*nModes,2).tolist() # Run parent init function super().__init__(classDef=UKNA_BSEN1991_2_walkers_joggers, dynsys_obj=modalsys_obj, bridgeClass=bridgeClass, mode_index=mode_index_list, analysis_type=["walkers","joggers"], verbose=False, **kwargs) # Save key variables as attributes self.bridgeClass = bridgeClass def plot_modal_params(self): fig, axarr = plt.subplots(2) fig.set_size_inches((6,8)) fig.subplots_adjust(hspace=0.4,right=0.8) # Get modal properties as used in analyses nModes = int(len(self.analysis_list)/2) mode_index = numpy.arange(0,nModes,1)+1 fv = [x.fv for x in self.analysis_list][:nModes] eta = [x.eta for x in self.analysis_list][:nModes] # Plot details of system ax = axarr[0] ax.bar(mode_index,fv) ax.set_xlabel("Mode index") ax.set_xticks(mode_index) ax.set_ylabel("$f_{d}$ (Hz)",fontsize=8.0) ax.set_title("Damped natural frequencies") ax = axarr[1] ax.bar(mode_index,eta) ax.set_xlabel("Mode index") ax.set_xticks(mode_index) ax.set_ylim([0,ax.get_ylim()[1]]) ax.set_ylabel("Damping ratio",fontsize=8.0) ax.yaxis.set_major_formatter(FuncFormatter(lambda y, _: '{:.1%}'.format(y))) ax.set_title("Modal damping") # Add equivilent axis for log dec ax2 = ax.twinx() def logdec(eta): return 2*numpy.pi*eta delta_incr = 0.01 new_yticks=numpy.arange(0,ax.get_ylim()[1],delta_incr/(2*numpy.pi)) ax2.set_ylim(ax.get_ylim()) ax2.set_yticks(new_yticks) ax2.set_yticklabels(["%.2f" % d for d in logdec(new_yticks)]) ax2.set_ylabel("Log decrement $\delta$") return fig def plot_response_stats(self, dynsys2plot=None, responses2plot=None, sharey='row', verbose=True): """ Plots statistics of responses computed, across all analyses *** Optional: * `dynsys2plot`, list of subsystem objects to plot responses for. If None (default) responses will be plotted for all subsystems * `responses2plot`, list of indices specifying responses to be plotted. If None (default) all responses will be plotted, for each sub-system * `sharey`, passed to pyplot.subplots() method. Controls whether subplots should share common y scale. """ print("Plotting response statistics...") nModes = int(len(self.analysis_list)/2) mode_index = numpy.arange(0,nModes,1)+1 fig_list = [] if dynsys2plot is None: dynsys2plot = self.dynsys_obj.DynSys_list else: if not (dynsys2plot in self.dynsys_obj.DynSys_list): raise ValueError("'dynsys2plot' does not belong " + "to system analysed!") if not isinstance(dynsys2plot,list): dynsys2plot = [dynsys2plot] # Loop over all systems and subsystems for dynsys_obj in dynsys2plot: # Get stats dict relating to this system stats_dict = self.stats_dict[dynsys_obj] # Get indices of responses to plot response_names = dynsys_obj.output_names if responses2plot is None: nResponses = len(response_names) responses2plot = list(range(nResponses)) else: nResponses = len(responses2plot) # Create new figure fig, axarr = plt.subplots(nResponses,2, sharex=True, sharey=sharey) fig.set_size_inches((14,8)) fig_list.append(fig) # Loop to plot all responses requested for row, r in enumerate(responses2plot): response_name = response_names[r] # Get stats for this response stats_dict_r = {} for stats_name, stats_vals in stats_dict.items(): stats_dict_r[stats_name]=stats_vals[:,:,r] if verbose: print("Plotting stats for response '%s':" % response_name) # Plot response stats for walkers and joggers responses for col, case in zip([0,1],['walkers','joggers']): ax = axarr[row,col] for i, stats_name in enumerate(['absmax','max','min']): stats_vals = stats_dict_r[stats_name][col,:] w=0.1 ax.bar(mode_index+(i-1)*w, stats_vals, width=w, label=stats_name ) ax.legend(loc='lower right',fontsize=6.0) if col==0: ax.set_ylabel(response_name,size=8.0) if row==0: ax.set_title("Responses to %s" % case) if row==nResponses-1: ax.set_xlabel("Mode index") ax.set_xticks(mode_index) fig.suptitle("Statistics of computed responses to pedestrian loading\n" + "UK NA to BS EN 1991-2, NA.2.44.4\n" + "Subsystem: '%s'\n" % dynsys_obj.name, fontsize=10.0) return fig_list
Ancestors (in MRO)
- PedestrianDynamics_transientAnalyses
- dyn_analysis.Multiple
- builtins.object
Static methods
def __init__(
self, modalsys_obj, bridgeClass='A', **kwargs)
Initialisation function
Required:
modalsys_obj
, modal system to which analysis relates
Optional:
bridgeClass
, string, either 'A', 'B', 'C', or 'D'. Refer Table NA.7 for description of bridge classes
Additional keyword arguments may be passed. These should relate to the
__init__()
function of the Multiple
class.
def __init__(self, modalsys_obj, bridgeClass='A', **kwargs): """ Initialisation function **** Required: * `modalsys_obj`, modal system to which analysis relates *** Optional: * `bridgeClass`, _string_, either 'A', 'B', 'C', or 'D'. Refer Table NA.7 for description of bridge classes Additional keyword arguments may be passed. These should relate to the `__init__()` function of the `Multiple` class. """ # Get mode indexs to loop over # Note modes sorted into ascending frequency but in conjugate pairs # Hence step through modes x2 nModes = modalsys_obj.GetSystemMatrices()["nDOF"] mode_index_list = numpy.arange(0,2*nModes,2).tolist() # Run parent init function super().__init__(classDef=UKNA_BSEN1991_2_walkers_joggers, dynsys_obj=modalsys_obj, bridgeClass=bridgeClass, mode_index=mode_index_list, analysis_type=["walkers","joggers"], verbose=False, **kwargs) # Save key variables as attributes self.bridgeClass = bridgeClass
def collate_stats(
self)
Collates computed statistics into a Pandas DataFrame, as follows:
-
Index is a MultiIndex, comprising the following levels:
-
0 : Name of sub-system
-
1 : Name of output / response
-
-
Columns are a MultiIndex comprising the following levels:
-
0 to -2 : Variables of analysis
vals2permute
-
-1 : Statistic name (e.g. 'max', 'min')
-
def collate_stats(self): """ Collates computed statistics into a Pandas DataFrame, as follows: * Index is a MultiIndex, comprising the following levels: * 0 : Name of sub-system * 1 : Name of output / response * Columns are a MultiIndex comprising the following levels: * 0 to -2 : Variables of analysis `vals2permute` * -1 : Statistic name (e.g. 'max', 'min') """ print("Collating statistics...") # Get list of all tstep_results objects associate with Multiple() if self.results_arr is None: raise ValueError("self.results_arr=None! No results to collate!") # Get lists of inputs permuted for analyses kwargs2permute = list(self.vals2permute.keys()) vals2permute = list(self.vals2permute.values()) # Where list contains objects, get their names for i in range(len(vals2permute)): if hasattr(vals2permute[i][0],'name'): vals2permute[i] = [x.name for x in vals2permute[i]] # Get list of systems and subsystems DynSys_list = self.dynsys_obj.DynSys_list stats_df = None for i, (index, results_obj) in enumerate(numpy.ndenumerate(self.results_arr)): # Get combination of vals2permute for results_obj combination = [vals[i] for i, vals in zip(index, vals2permute)] # Get stats from results_obj stats_df_inner = results_obj.get_response_stats_df() # Loop over all sub-systems for df, sys in zip(stats_df_inner,DynSys_list): # Prepend system to index df = pd.concat([df], axis=0, keys=[sys.name], names=['System','Response']) # Prepare combination values to column MultiIndex tuples = [(*combination,col) for col in df.columns] df.columns = pd.MultiIndex.from_tuples(tuples) df.columns.names = [*kwargs2permute,'Statistic'] # Append results to DataFrame stats_df = pd.concat([stats_df, df],axis=1) self.stats_df = stats_df return stats_df
def plot_modal_params(
self)
def plot_modal_params(self): fig, axarr = plt.subplots(2) fig.set_size_inches((6,8)) fig.subplots_adjust(hspace=0.4,right=0.8) # Get modal properties as used in analyses nModes = int(len(self.analysis_list)/2) mode_index = numpy.arange(0,nModes,1)+1 fv = [x.fv for x in self.analysis_list][:nModes] eta = [x.eta for x in self.analysis_list][:nModes] # Plot details of system ax = axarr[0] ax.bar(mode_index,fv) ax.set_xlabel("Mode index") ax.set_xticks(mode_index) ax.set_ylabel("$f_{d}$ (Hz)",fontsize=8.0) ax.set_title("Damped natural frequencies") ax = axarr[1] ax.bar(mode_index,eta) ax.set_xlabel("Mode index") ax.set_xticks(mode_index) ax.set_ylim([0,ax.get_ylim()[1]]) ax.set_ylabel("Damping ratio",fontsize=8.0) ax.yaxis.set_major_formatter(FuncFormatter(lambda y, _: '{:.1%}'.format(y))) ax.set_title("Modal damping") # Add equivilent axis for log dec ax2 = ax.twinx() def logdec(eta): return 2*numpy.pi*eta delta_incr = 0.01 new_yticks=numpy.arange(0,ax.get_ylim()[1],delta_incr/(2*numpy.pi)) ax2.set_ylim(ax.get_ylim()) ax2.set_yticks(new_yticks) ax2.set_yticklabels(["%.2f" % d for d in logdec(new_yticks)]) ax2.set_ylabel("Log decrement $\delta$") return fig
def plot_response_stats(
self, dynsys2plot=None, responses2plot=None, sharey='row', verbose=True)
Plots statistics of responses computed, across all analyses
Optional:
-
dynsys2plot
, list of subsystem objects to plot responses for. If None (default) responses will be plotted for all subsystems -
responses2plot
, list of indices specifying responses to be plotted. If None (default) all responses will be plotted, for each sub-system -
sharey
, passed to pyplot.subplots() method. Controls whether subplots should share common y scale.
def plot_response_stats(self, dynsys2plot=None, responses2plot=None, sharey='row', verbose=True): """ Plots statistics of responses computed, across all analyses *** Optional: * `dynsys2plot`, list of subsystem objects to plot responses for. If None (default) responses will be plotted for all subsystems * `responses2plot`, list of indices specifying responses to be plotted. If None (default) all responses will be plotted, for each sub-system * `sharey`, passed to pyplot.subplots() method. Controls whether subplots should share common y scale. """ print("Plotting response statistics...") nModes = int(len(self.analysis_list)/2) mode_index = numpy.arange(0,nModes,1)+1 fig_list = [] if dynsys2plot is None: dynsys2plot = self.dynsys_obj.DynSys_list else: if not (dynsys2plot in self.dynsys_obj.DynSys_list): raise ValueError("'dynsys2plot' does not belong " + "to system analysed!") if not isinstance(dynsys2plot,list): dynsys2plot = [dynsys2plot] # Loop over all systems and subsystems for dynsys_obj in dynsys2plot: # Get stats dict relating to this system stats_dict = self.stats_dict[dynsys_obj] # Get indices of responses to plot response_names = dynsys_obj.output_names if responses2plot is None: nResponses = len(response_names) responses2plot = list(range(nResponses)) else: nResponses = len(responses2plot) # Create new figure fig, axarr = plt.subplots(nResponses,2, sharex=True, sharey=sharey) fig.set_size_inches((14,8)) fig_list.append(fig) # Loop to plot all responses requested for row, r in enumerate(responses2plot): response_name = response_names[r] # Get stats for this response stats_dict_r = {} for stats_name, stats_vals in stats_dict.items(): stats_dict_r[stats_name]=stats_vals[:,:,r] if verbose: print("Plotting stats for response '%s':" % response_name) # Plot response stats for walkers and joggers responses for col, case in zip([0,1],['walkers','joggers']): ax = axarr[row,col] for i, stats_name in enumerate(['absmax','max','min']): stats_vals = stats_dict_r[stats_name][col,:] w=0.1 ax.bar(mode_index+(i-1)*w, stats_vals, width=w, label=stats_name ) ax.legend(loc='lower right',fontsize=6.0) if col==0: ax.set_ylabel(response_name,size=8.0) if row==0: ax.set_title("Responses to %s" % case) if row==nResponses-1: ax.set_xlabel("Mode index") ax.set_xticks(mode_index) fig.suptitle("Statistics of computed responses to pedestrian loading\n" + "UK NA to BS EN 1991-2, NA.2.44.4\n" + "Subsystem: '%s'\n" % dynsys_obj.name, fontsize=10.0) return fig_list
def plot_stats(
self, stat='absmax', sys=None, **kwargs)
Produces a plot of a given statistic, taken across multiple analyses
Optional:
-
stat
, name of statistic to be plotted (e.g. 'max', 'min') -
sys
, instance ofDynSys
class, or string to denote name of system, used to select system whose outputs are to be plotted. If None then seperate figures will be produced for each subsystem.
See docstring for plot_stats_for_system()
for other keyword
arguments that may be passed.
def plot_stats(self, stat='absmax', sys=None, **kwargs): """ Produces a plot of a given statistic, taken across multiple analyses Optional: * `stat`, name of statistic to be plotted (e.g. 'max', 'min') * `sys`, instance of `DynSys` class, or string to denote name of system, used to select system whose outputs are to be plotted. If None then seperate figures will be produced for each subsystem. _See docstring for `plot_stats_for_system()` for other keyword arguments that may be passed._ """ # Get list of system names to loop over if sys is None: sys_names = [obj.name for obj in self.dynsys_obj.DynSys_list] else: # Get name of system if isinstance(sys,str): sys_name = sys else: sys_name = sys.name # get name from object sys_names = [sys_name] # list of length 1 # Produce a seperate figure for each sub-system responses fig_list = [] fig_list = [] for sys_name in sys_names: _fig_list = self.plot_stats_for_system(sys_name=sys_name, stat=stat, **kwargs) fig_list.append(_fig_list) # Return list of figures, one for each subsystem return fig_list
def plot_stats_for_system(
self, sys_name, stat, max_responses_per_fig=5, subplot_kwargs={})
Produces a plot of a given statistic, taken across multiple analyses for a specified sub-system
Required:
-
sys_name
, string giving name of system -
stat
, string to specify statistic to be plotted. E.g. 'absmax'
Optional:
max_responses_per_fig
, integer to denote maximum number of responses to be plotted in each figure. If None, all responses will be plotted via a single figure. Default value (=4) should give nice plots in most cases.
Users are advised to tweak the appearance
of figures, e.g. using the pyplot.subplots_adjust()
method.
subplot_kwargs
, dict of keyword arguments to be passed topyplot.subplots()
method, to customise subplots (e.g. share axes)
def plot_stats_for_system(self,sys_name,stat, max_responses_per_fig:int=5, subplot_kwargs={} ): """ Produces a plot of a given statistic, taken across multiple analyses for a specified sub-system Required: * `sys_name`, string giving name of system * `stat`, string to specify statistic to be plotted. E.g. 'absmax' Optional: * `max_responses_per_fig`, integer to denote maximum number of responses to be plotted in each figure. If None, all responses will be plotted via a single figure. Default value (=4) should give nice plots in most cases. _Users are advised to tweak the appearance of figures, e.g. using the `pyplot.subplots_adjust()` method._ * `subplot_kwargs`, dict of keyword arguments to be passed to `pyplot.subplots()` method, to customise subplots (e.g. share axes) """ # Re-collate statistics as required if self.stats_df is None: stats_df = self.collate_stats() else: stats_df = self.stats_df # Slice for requested stat try: stats_df = stats_df.xs(stat,level=-1,axis=1) except KeyError: raise KeyError("Invalid statistic selected!") # Get stats for just this system df_thissys = stats_df.xs(sys_name,level=0,axis=0) # Obtain responses names for this subsystem response_names = df_thissys.index.values nResponses = len(response_names) if max_responses_per_fig is None: max_responses_per_fig = nResponses fig_list = [] for _response_names in chunks(response_names,max_responses_per_fig): # Create figure, with one subplot per response fig, axlist = plt.subplots(len(_response_names), sharex=True, **subplot_kwargs) fig_list.append(fig) for i, (r,ax) in enumerate(zip(_response_names,axlist)): # Get series for this response df = df_thissys.loc[r] # Reshape such that index will be x-variable for plot df = df.unstack() # Make plot ax = df.plot(ax=ax,legend=False) ax.set_ylabel(r, rotation=0, fontsize='small', horizontalAlignment='right', verticalAlignment='center') if i==0: # Add legend to figure fig.legend(ax.lines, df.columns,fontsize='x-small') fig.subplots_adjust(left=0.15,right=0.95) fig.align_ylabels() return fig_list
def run(
self, save=True, solveInParallel=False)
Runs multiple dynamic analyses, as defined by __init__
In principle this can be done in parallel, to efficiently use all avaliable cores and give improve runtime.
Parallel processing not yet implemented due to outstanding bug
associated with using multiprocessing
module from within Spyder
def run(self,save=True,solveInParallel=False): """ Runs multiple dynamic analyses, as defined by `__init__` *** In principle this can be done in parallel, to efficiently use all avaliable cores and give improve runtime. _Parallel processing not yet implemented due to outstanding bug associated with using `multiprocessing` module from within Spyder_ """ print("Running multiple `{0}`\n".format(self.analysisType)) tic=timeit.default_timer() # Run analyses using parallel processing (if possible) if solveInParallel: print("Parallel processing not yet implemented due to " + "outstanding bug associated with using "+ "`multiprocessing` module from within Spyder\n"+ "A single-process analysis will be carried out instead.") # Run all pre-defined analyses for i, x in enumerate(self.analysis_list): print("Analysis #%04d of #%04d" % (i+1, len(self.analysis_list))) x.run(saveResults=False) print("")#clear line for emphasis toc=timeit.default_timer() print("Multiple `%s` analysis complete after %.3f seconds!\n" % (self.analysisType,(toc-tic))) # Reshape results objects in ndarray results_obj_list = [x.tstep_obj.results_obj for x in self.analysis_list] reqdShape = self.vals2permute_shape results_arr = numpy.reshape(results_obj_list,reqdShape) self.results_arr = results_arr # Collate statistics self.collate_stats() # Pickle results if save: self.save()
def save(
self, fName=None)
Serialises object to file fName
def save(self,fName=None): """ Serialises object to file `fName` """ fName=self._pickle_fName(fName) with open('{0}'.format(fName), 'wb') as dill_file: print("\nSerialising `{0}` object to `{1}`".format(self.__class__.__name__,fName)) dill.dump(self, dill_file) print("Serialisation complete!\n")
Instance variables
var bridgeClass
class SteadyStateCrowdLoading
Class to implement steady state crowd loading analysis
Will usually be used as base class for analyses to specific codes, e.g. HiVoSS or UK NA to BS EN 1991-2
class SteadyStateCrowdLoading(): """ Class to implement steady state crowd loading analysis _Will usually be used as base class for analyses to specific codes, e.g. HiVoSS or UK NA to BS EN 1991-2_ """ def __init__(self, modalsys_obj, width_func_list=[3.0], modeshapes_fname_arr=["modeshapes_edge1.csv", "modeshapes_edge2.csv"], name=None, verbose=False, makePlot=False): """ Initialisation function *** Required: * `modalsys_obj`, object defining modal system to be analysed * `modeshapes_fname_arr`, array of strings, of shape (nRegions,2) denoting files containing modeshape variation along the edges of each deck region: * Axis 1 to be of shape = 2 (1 file for each edge) * Axis 0 length corresponds to number of deck regions, nRegions * `width_func_list`, list of functions to denote how width (m) of each deck region varies with chainage _along the centreline of the region_: * List length = nRegions required (this is checked) * List of floats may be provided if deck width is constant within each regions. But list length requirement above still applies *** Optional: * calc_Seff`, _boolean_, dictates whether effective span will be computed according to Figure NA.7. If False, overall span will be used (conservative) * `nLcuts`, _integer_, denotes number of transverse cuts to make when evaluating area integrals. Increasing number will increase accuracy, at the cost of longer computation time """ # Create default name for analysis object if name is None: name = modalsys_obj.name + " - Crowd Loading" self.name = name """ Name assigned to class instance """ # Check `width_func_list` input width_func_list = self._check_width_func_list(width_func_list) nRegions = len(width_func_list) self.nRegions = nRegions """ Integer, denotes number deck regions """ # Check `modeshapes_fname_arr` input if modeshapes_fname_arr is None: # Use already-defined modeshope for system mFunc = modalsys_obj.modeshapeFunc mFunc_list = numpy.array([[mFunc,mFunc]]) Ltrack = modalsys_obj.Ltrack Ltrack_list = numpy.array([[Ltrack,Ltrack]]) max_ordinates = numpy.max(numpy.abs(mFunc.y)) else: modeshapes_fname_arr = numpy.array(modeshapes_fname_arr) self._check_modeshapes_fname_arr(modeshapes_fname_arr,nRegions) # Read modeshape data from file, determine interpolation functions mFunc_list, Ltrack_list, max_ordinates = self._read_modeshapes(modalsys_obj, modeshapes_fname_arr) self.max_modeshape = max_ordinates """ Maximum |modeshape| ordinates, evaluated over all deck regions. Array of length equal to number of modes considered """ # Define deck regions deck_regions_list = [] for r in range(nRegions): deckregion_obj = _DeckStrip(name="Region %d" % r, L=Ltrack_list[r], width_func=width_func_list[r], func_edge1=mFunc_list[r,0], func_edge2=mFunc_list[r,1]) deck_regions_list.append(deckregion_obj) self.deck_regions_list = deck_regions_list """ List of DeckStrip classes, defining deck strips """ # Calculate deck area self.deck_area = self.calc_deck_area() """ Total deck area ($m^2$), obtained by integration over all deck strips """ # Get damped natural frequencies for system eig_results = modalsys_obj.CalcEigenproperties() self.f_d = eig_results["f_d"][0::2] """ Damped natural frequencies of system being analysed """ self.eta = eig_results["eta"][0::2] """ Damping ratios for system modes """ # Save other attributes self.modalsys_obj = modalsys_obj """ `ModalSys` object, corresponding to system being analysed """ self.response_names = [] """ Names of responses computed """ self.response_amplitudes = [] """ Amplitudes of defined responses """ if verbose: self.print_details() def run(self, target_mode:int, load_intensity:float=None, forcing_freq:float=None, run_freq_method=True, run_tstep=False,tstep_kwargs:dict={}, verbose=False, makePlot={'calc_modal_areas':False,'calc_load_intensity':True} ): """ Runs calculations to determine steady-state crowd loading responses *** **Required:** * `target_mode`, _integer_, denotes the mode at which resonance is to be targeted *** **Optional:** * `load_intensity`, _float_, denotes the load intensity (in N/m2) to be applied. If _None_, `load_intensity` to be calculated using `calc_load_intensity()` method function. * `forcing_freq`, _float_, defines forcing frequency (in Hz) at which response is to be determined. For ordinary modal systems with no added dampers this value should normally be set to equal the damped natural frequency for the mode given by `target_mode`, as this will ensure the worst-case resonance response is obtained. This is thus the default behaviour. However is some circumstances (e.g. systems with multiple TMDs appended) it is necessary to evaluate responses at various assumed forcing frequencies, hence why this optional input is provided to the user. * `run_freq_method`, _boolean_, if True the amplitude of modal and real-world responses at steady-state will be determine by frequency transfer matrix method. This is the default behaviour, as it allows the steady-state amplitudes of responses to be efficiently calculated. In some cases (e.g. for systems with constraints, for which the `DynSys.CalcFreqResponse()` method is not currently implemented), it may be desirable to carry out a time-stepping analysis instead. This can be done by setting `run_freq_method=False` and `run_tstep=True` (see documentation in next bullet for further details) * `run_tstep`, _boolean_, if True a time-stepping analysis will be run, which can be used to verify the the steady-state resonance solution that is ordinarily obtained. * `tstep_kwargs`, _dict_ of keyword arguments to be passed when initialising time-stepping analysis using the `TStep()` class. Only used if `run_tstep` is True. Refer [documentation](../docs/tstep.html) for `TStep()` class for more details. * `makePlot`, _dict_ or _bool_, used to control which plots are produced [documentation]: """ if verbose: print("\n**** Running steady-state crowd loading analysis ****") # Handle makePlot if isinstance(makePlot,bool): makePlot_dict = {} makePlot_dict['calc_modal_areas']=makePlot makePlot_dict['calc_load_intensity']=makePlot else: makePlot_dict = makePlot self.target_mode = target_mode if verbose: print("Target mode index: %d" % target_mode) # Define forcing frequency if forcing_freq is None: # Set to be damped natural frequency of mode being targeted forcing_freq = numpy.abs(self.f_d[target_mode]) else: # where input value provided via optional arg # Check float type is provided if not isinstance(forcing_freq,float): raise ValueError("Error: float expected for `forcing_freq`") self.forcing_freq = forcing_freq # Obtain area integrals over all deck strips modal_areas = self.calc_modal_areas(target_mode, makePlot=makePlot_dict['calc_modal_areas']) self.modal_areas = modal_areas # Calculate load intensity, if not provided directly if load_intensity is None: load_intensity, d = self.calc_load_intensity(mode_index=target_mode, fv=forcing_freq, verbose=verbose, makePlot=makePlot_dict['calc_load_intensity']) else: d = {} # null dict self.load_intensity = load_intensity """ Uniform load intensity (in N/m2) to be applied to deck regions """ self.load_intensity_rslts_dict = d """ _Dict_ to store key results from load intensity calculation """ if verbose: print("UDL intensity (N/m2): %.2f" % load_intensity) # Multiply by load intensity to obtain modal force amplitudes # Note modal force amplitudes should be regarded as complex variables # (although here they will be real-valued) whereby negative values # imply a force which is 180 degrees out of phase modal_forces = self.load_intensity * modal_areas self.modal_forces = modal_forces # Save settings used to determine which analyses run self.run_freq_method = run_freq_method self.run_tstep = run_tstep if verbose: print("Forcing frequency: %.2f Hz" % forcing_freq) # Run frequency domain method to determine steady-state responses if run_freq_method: # Evaluate frequency transfer functions G(f) relating applied forces # (i.e. modal forces in this application) to DOF accelerations, # displacements and velocities # (i.e. modal disp, vel, accn in this application) # # It is only necessary to evaluate G(f) at resonant frequency for mode # being targeted, although for validation / plotting it may be # advantageous to evaluate G(f) over a wider frequency range # # By definition, G(f1) is a matrix that can be used to relate # sinusoidal modal forces to sinusoidal modal accelerations, # given forcing frequency f1. # # Note: G(f) will in general be complex-valued; the should be # intepreted to denote scaling and phase lag # Calculate matrix mapping applied loads to state DOFs def Calc_Gf(f): rslts = self.modalsys_obj.CalcFreqResponse(fVals=f, C=numpy.zeros((0,)), verbose=verbose) f = rslts["f"] G_f = rslts["G_f"] return f, G_f G_state_vars = Calc_Gf(forcing_freq)[1] G_state_vars = G_state_vars[0,:,:] # Retain only submatrix mapping modal forces to state DOFs nModes = len(modal_forces) G_state_vars = G_state_vars[:,:nModes] self.G_state_vars = G_state_vars # Get submatrix mapping applied forces to DOF accelerations nDOF = int(G_state_vars.shape[0]/3) G_acc = G_state_vars[2*nDOF:,:] # Get submatrix mapping modal forces to m accelerations G_modalacc = G_acc[:nModes,:nModes] # Calculate modal acceleration amplitudes using transfer matrix modal_accn = numpy.ravel(G_modalacc @ modal_forces) self.modal_accn = modal_accn # Adjust phase such that acc for target mode has phase=0 z1 = 1.0 * numpy.exp(1j*numpy.angle(modal_accn[target_mode])) modal_accn = modal_accn / z1 # Get modal acceleration of targeted mode # Note phase=0 due to above adjustment, hence ok to take real part only modal_accn_target = numpy.real(modal_accn[target_mode]) self.modal_accn_target = modal_accn_target # Calculate responses amplitudes using full transfer matrix self.calc_responses() if verbose: self.print_results() print("**** Analysis complete! ****\n") if run_tstep: # Run time-stepping analysis to validate /check the above approach print("***Running time-stepping analysis, as cross-check of "+ "steady-state resonance approach...***") # Define sinusoidal load function consistent with above def load_func(t): return modal_forces * numpy.sin(2*numpy.pi*forcing_freq*t) load_func_dict = {self.modalsys_obj : load_func} # Define and run time-stepping analysis tstep_analysis = tstep.TStep(dynsys_obj=self.modalsys_obj, force_func_dict=load_func_dict, **tstep_kwargs) tstep_results = tstep_analysis.run(verbose=verbose) self.tstep_results = tstep_results """ Results from time-stepping analysis; instance of `TStep_Results()` """ print("***Time-stepping analysis complete!***") def calc_responses(self,output_mtrx=None,output_names=None): """ Calculate steady-steady amplitude response of defined responses """ # Get output matrix and output names, if not passed to this function if output_mtrx is None: output_mtrx, output_names = self.modalsys_obj.GetOutputMtrx(all_systems=True) # Get frequency response matrix for {disp,vel,accn} variables # (note: assumes already calculated!) G_state_vars = self.G_state_vars # Get modal forces modal_forces = self.modal_forces # Calculate responses at steady-state nResponses = output_mtrx.shape[0] if nResponses>0: response_amplitudes = numpy.ravel(output_mtrx @ G_state_vars @ modal_forces) self.response_amplitudes = response_amplitudes self.response_names = output_names else: print("No output matrix provided. " + "Real-world responses have not been evaluated") def print_details(self): """ Prints details about the analysis that has been defined """ print("\nDetails of deck regions:") deck_regions_list = self.deck_regions_list print("nRegions: %d\n" % len(deck_regions_list)) for deckregion_obj in deck_regions_list: deckregion_obj.print_details() print("Total deck area (m2):\n%.2f" % self.deck_area) print("Damped natural frequencies of system (Hz):\n{0}" .format(self.f_d)) def print_results(self): """ Prints key results from analysis """ def printcomplex(vals): """ Function for printing magnitude and phase of complex-valued vals """ print("Magnitudes:\n{0}".format(numpy.abs(vals))) print("Phase angles (rad):\n{0}".format(numpy.angle(vals))) print("") print("Modal forces (N):") printcomplex(self.modal_forces) print("") if self.run_freq_method: print("Modal accelerations (m/s2):") printcomplex(self.modal_accn) print("") print("Responses:") print(self.response_names) printcomplex(self.response_amplitudes) def get_response_amplitudes(self,print_phase=False): """ Returns response amplitudes as pandas dataframe """ data = numpy.abs(self.response_amplitudes) index_list = self.response_names columns_list = ["Amplitude"] if print_phase: phase_angles = numpy.angle(self.response_amplitudes) data = numpy.hstack((data,phase_angles)) columns_list.append("Phase angle") df = pandas.DataFrame(data=data, index=index_list, columns=columns_list) return df def get_response_stats(self): """ Returns response stats for time-stepping analysis carried out """ if not self.run_tstep: raise ValueError("Error: no time-stepping analysis run, " + "cannot retrieve results!") # Get dict of dict of stats stats_dict_dict = self.tstep_results.CalcResponseStats() # Loop over dict and return stats as list of pandas dataframes df_list = [] for sys_obj, stats_dict in stats_dict_dict.items(): sys_name = sys_obj.name response_names = sys_obj.output_names df = pandas.DataFrame(data=stats_dict, index=[sys_name + " - " + x for x in response_names]) df_list.append(df) # Merge pandas dataframes to return a single dataframe # with all stats for all responses for all systems df_merged = pandas.concat(df_list) # Transpose such that stats as indexes, responses as columns df_merged = df_merged.transpose() return df_merged def plot_results(self,overlay_val:float=None): """ Plots results of steady-state response analysis If time-stepping analysis has also been carried out (as a cross-check) results of that analysis will be plotted also *** **Required:** (No arguments required) *** **Optional:** * `overlay_val`, _float_, value to overlay onto response plots. (A common use of this would be to overlay the applicable acceleration limit, against which responses are to be compared) *** **Returns:** _Nested list_ of figure objects """ fig_list = [] if self.run_freq_method: fig = self.plot_results_freq_method(overlay_val=overlay_val) fig_list.append(fig) if self.run_tstep: tstep_fig_list = self.plot_results_tstep_method(overlay_val=overlay_val) fig_list += tstep_fig_list return fig_list def plot_results_freq_method(self,overlay_val=None): """ Plot results of frequency domain analysis method, if avaliable """ if not self.run_freq_method: print("No results avaliable for frequency domain method\n" + "Could not plot results") return None # Produce and configure figure fig, axarr = plt.subplots(2) fig.suptitle("Steady-state pedestrian crowd loading analysis results\n"+ "Target mode: %d" % self.target_mode) fig.set_size_inches((14,8)) fig.subplots_adjust(hspace=0.3) # Produce bar plot of modal accelerations self.plot_results_modal_accn(axarr[0]) # Produce bar plot of outputs self.plot_results_responses(axarr[1],overlay_val=overlay_val) return fig def plot_results_modal_accn(self,ax=None): """ Produce bar plot of steady-state modal accelerations """ # Define axes to plot to, if not done so already if ax is None: fig,ax = plt.subplots() fig.set_size_inches((14,4)) vals = self.modal_accn # Get or define tick labels if hasattr(self.modalsys_obj,'mode_IDs'): ticks = self.modalsys_obj.mode_IDs else: ticks = ["Mode %d" % (i+1) for i in range(self.modalsys_obj.nDOF)] ax.barh(range(len(vals)),numpy.abs(vals),tick_label=ticks) ax.tick_params(labelsize=8) ax.set_title("Steady-state modal acceleration amplitudes (m/$s^2$)") def plot_results_responses(self,ax=None,overlay_val=None): """ Produce bar plot of steady-state responses """ # Define axes to plot to, if not done so already if ax is None: fig,ax = plt.subplots() fig.set_size_inches((14,4)) # Get data to use in plot vals = self.response_amplitudes ticks = self.response_names # Flatten list in case of nested lists (only when multiple subsystems) if vals.shape != ticks.shape: ticks = list(itertools.chain(*self.response_names.tolist())) ax.barh(range(len(vals)),numpy.abs(vals),tick_label=ticks) ax.tick_params(labelsize=8) ax.set_title("Steady-state response amplitudes") if overlay_val is not None: ax.axvline(overlay_val,color='darkorange') def plot_results_tstep_method(self,overlay_val=None): """ Plot results of time-stepping analysis, if avaliable """ if not self.run_tstep: print("No results avaliable for time-stepping method\n" + "Could not plot results") return [] # Define +/- values based on overlay_val, if provided if overlay_val is not None: overlay_val=[+overlay_val,-overlay_val] else: overlay_val=None fig_list = self.tstep_results.PlotResults(useCommonPlot=False, useCommonScale=False, y_overlay=overlay_val) # Overlay expected modal acceleration amplitudes # (from steady-state analysis, if run) if self.run_freq_method: state_results_fig = fig_list[0][0] ax = state_results_fig.get_axes()[3] # plot of modal accelerations accn_val = numpy.abs(self.modal_accn_target) # value to overlay ax.axhline(+accn_val,color='r',linestyle='--',alpha=0.3) ax.axhline(-accn_val,color='r',linestyle='--',alpha=0.3) # Overlay expected response amplitudes # (from steady-state analysis, if run) response_results_fig_list = fig_list[1] axarr = [x.get_axes() for x in response_results_fig_list] axarr = list(itertools.chain(*axarr)) # flatten list val2overlay = numpy.abs(self.response_amplitudes) for ax, val in zip(axarr,val2overlay): ax.axhline(+val,color='r',linestyle='--',alpha=0.3) ax.axhline(-val,color='r',linestyle='--',alpha=0.3) return fig_list def _check_width_func_list(self, width_func_list): """ Performs checks on input `width_func_list` as supplied to `__init__()` """ # Check that list type if not isinstance(width_func_list,list): # Handle case of single float supplied by converting to list if isinstance(width_func_list,float): width_func_list = [width_func_list] else: raise ValueError("`width_func_list`: list expected") # Check all list items are functions or floats for i, width_func in enumerate(width_func_list): if not isfunction(width_func): if not isinstance(width_func,float): raise ValueError("`width_func_list[%d]`" % i + " is not a function, as required") return width_func_list def _check_modeshapes_fname_arr(self, modeshapes_fname_arr, nRegions_expected): # Check shapes nRegions, nEdges = modeshapes_fname_arr.shape if nRegions_expected != nRegions: raise ValueError("Unexpected `nRegions`\n: " + "Conflicts with nRegions inferred " + "from `width_func_list` input") if nEdges != 2: raise ValueError("nEdges = 2 required\n " + "Modeshapes at edges " + "of each regions must be provided via " + "`modeshapes_fname_arr` input") # Check all files exist for index, file_str in numpy.ndenumerate(modeshapes_fname_arr): if not os.path.isfile(file_str): raise ValueError("`modeshapes_fname_list[{0}]`" % index + " does not exist!") def _read_modeshapes(self,modalsys_obj,modeshapes_fname_arr,num=1000): """ Reads modeshapes from file *** Required: * `modeshapes_fname_arr`, array of strings, denoting file paths """ Ltrack_list = [] modeshapeFunc_list = [] absmax_vals_list = [] # Loop through array containing modeshape filenames for fname_list in modeshapes_fname_arr: Ltrack_list_inner = [] modeshapeFunc_list_inner = [] absmax_vals_list_inner = [] for fName in fname_list: mfunc, L = modalsys_obj.DefineModeshapes(fName=fName) # Evaluate to obtain maximum ordinates m_vals = mfunc(numpy.linspace(0,L,num)) absmax_vals = numpy.max(numpy.abs(m_vals),axis=0) absmax_vals_list_inner.append(absmax_vals.tolist()) # Append to inner lists Ltrack_list_inner.append(L) modeshapeFunc_list_inner.append(mfunc) # Append to create list of lists Ltrack_list.append(Ltrack_list_inner) modeshapeFunc_list.append(modeshapeFunc_list_inner) absmax_vals_list.append(absmax_vals_list_inner) # Obtain overall max modeshape ordinate absmax_vals_list = numpy.array(absmax_vals_list) max_ordinate = numpy.max(absmax_vals_list,axis=(0,1)) # Convert lists to arrays Ltrack_list = numpy.array(Ltrack_list) modeshapeFunc_list = numpy.array(modeshapeFunc_list) return modeshapeFunc_list, Ltrack_list, max_ordinate def calc_deck_area(self): """ Returns the total deck area, taken across all deck strips """ total_area = sum([x.calc_area() for x in self.deck_regions_list]) return float(total_area) def calc_modal_areas(self,mode_index,makePlot=False): """ Returns modal area integrals """ return numpy.sum([x.integrate(mode_index,makePlot=makePlot) for x in self.deck_regions_list],axis=0) def calc_load_intensity(self,**kwargs): """ _Method to be overriden by derived classes_ """ raise ValueError("`calc_load_intensity` not implemented!" + "To be overriden by derived class method")
Ancestors (in MRO)
- SteadyStateCrowdLoading
- builtins.object
Static methods
def __init__(
self, modalsys_obj, width_func_list=[3.0], modeshapes_fname_arr=['modeshapes_edge1.csv', 'modeshapes_edge2.csv'], name=None, verbose=False, makePlot=False)
Initialisation function
Required:
-
modalsys_obj
, object defining modal system to be analysed -
modeshapes_fname_arr
, array of strings, of shape (nRegions,2) denoting files containing modeshape variation along the edges of each deck region:-
Axis 1 to be of shape = 2 (1 file for each edge)
-
Axis 0 length corresponds to number of deck regions, nRegions
-
-
width_func_list
, list of functions to denote how width (m) of each deck region varies with chainage along the centreline of the region:-
List length = nRegions required (this is checked)
-
List of floats may be provided if deck width is constant within each regions. But list length requirement above still applies
-
Optional:
-
calc_Seff`, boolean, dictates whether effective span will be computed according to Figure NA.7. If False, overall span will be used (conservative)
-
nLcuts
, integer, denotes number of transverse cuts to make when evaluating area integrals. Increasing number will increase accuracy,
at the cost of longer computation time
def __init__(self, modalsys_obj, width_func_list=[3.0], modeshapes_fname_arr=["modeshapes_edge1.csv", "modeshapes_edge2.csv"], name=None, verbose=False, makePlot=False): """ Initialisation function *** Required: * `modalsys_obj`, object defining modal system to be analysed * `modeshapes_fname_arr`, array of strings, of shape (nRegions,2) denoting files containing modeshape variation along the edges of each deck region: * Axis 1 to be of shape = 2 (1 file for each edge) * Axis 0 length corresponds to number of deck regions, nRegions * `width_func_list`, list of functions to denote how width (m) of each deck region varies with chainage _along the centreline of the region_: * List length = nRegions required (this is checked) * List of floats may be provided if deck width is constant within each regions. But list length requirement above still applies *** Optional: * calc_Seff`, _boolean_, dictates whether effective span will be computed according to Figure NA.7. If False, overall span will be used (conservative) * `nLcuts`, _integer_, denotes number of transverse cuts to make when evaluating area integrals. Increasing number will increase accuracy, at the cost of longer computation time """ # Create default name for analysis object if name is None: name = modalsys_obj.name + " - Crowd Loading" self.name = name """ Name assigned to class instance """ # Check `width_func_list` input width_func_list = self._check_width_func_list(width_func_list) nRegions = len(width_func_list) self.nRegions = nRegions """ Integer, denotes number deck regions """ # Check `modeshapes_fname_arr` input if modeshapes_fname_arr is None: # Use already-defined modeshope for system mFunc = modalsys_obj.modeshapeFunc mFunc_list = numpy.array([[mFunc,mFunc]]) Ltrack = modalsys_obj.Ltrack Ltrack_list = numpy.array([[Ltrack,Ltrack]]) max_ordinates = numpy.max(numpy.abs(mFunc.y)) else: modeshapes_fname_arr = numpy.array(modeshapes_fname_arr) self._check_modeshapes_fname_arr(modeshapes_fname_arr,nRegions) # Read modeshape data from file, determine interpolation functions mFunc_list, Ltrack_list, max_ordinates = self._read_modeshapes(modalsys_obj, modeshapes_fname_arr) self.max_modeshape = max_ordinates """ Maximum |modeshape| ordinates, evaluated over all deck regions. Array of length equal to number of modes considered """ # Define deck regions deck_regions_list = [] for r in range(nRegions): deckregion_obj = _DeckStrip(name="Region %d" % r, L=Ltrack_list[r], width_func=width_func_list[r], func_edge1=mFunc_list[r,0], func_edge2=mFunc_list[r,1]) deck_regions_list.append(deckregion_obj) self.deck_regions_list = deck_regions_list """ List of DeckStrip classes, defining deck strips """ # Calculate deck area self.deck_area = self.calc_deck_area() """ Total deck area ($m^2$), obtained by integration over all deck strips """ # Get damped natural frequencies for system eig_results = modalsys_obj.CalcEigenproperties() self.f_d = eig_results["f_d"][0::2] """ Damped natural frequencies of system being analysed """ self.eta = eig_results["eta"][0::2] """ Damping ratios for system modes """ # Save other attributes self.modalsys_obj = modalsys_obj """ `ModalSys` object, corresponding to system being analysed """ self.response_names = [] """ Names of responses computed """ self.response_amplitudes = [] """ Amplitudes of defined responses """ if verbose: self.print_details()
def calc_deck_area(
self)
Returns the total deck area, taken across all deck strips
def calc_deck_area(self): """ Returns the total deck area, taken across all deck strips """ total_area = sum([x.calc_area() for x in self.deck_regions_list]) return float(total_area)
def calc_load_intensity(
self, **kwargs)
Method to be overriden by derived classes
def calc_load_intensity(self,**kwargs): """ _Method to be overriden by derived classes_ """ raise ValueError("`calc_load_intensity` not implemented!" + "To be overriden by derived class method")
def calc_modal_areas(
self, mode_index, makePlot=False)
Returns modal area integrals
def calc_modal_areas(self,mode_index,makePlot=False): """ Returns modal area integrals """ return numpy.sum([x.integrate(mode_index,makePlot=makePlot) for x in self.deck_regions_list],axis=0)
def calc_responses(
self, output_mtrx=None, output_names=None)
Calculate steady-steady amplitude response of defined responses
def calc_responses(self,output_mtrx=None,output_names=None): """ Calculate steady-steady amplitude response of defined responses """ # Get output matrix and output names, if not passed to this function if output_mtrx is None: output_mtrx, output_names = self.modalsys_obj.GetOutputMtrx(all_systems=True) # Get frequency response matrix for {disp,vel,accn} variables # (note: assumes already calculated!) G_state_vars = self.G_state_vars # Get modal forces modal_forces = self.modal_forces # Calculate responses at steady-state nResponses = output_mtrx.shape[0] if nResponses>0: response_amplitudes = numpy.ravel(output_mtrx @ G_state_vars @ modal_forces) self.response_amplitudes = response_amplitudes self.response_names = output_names else: print("No output matrix provided. " + "Real-world responses have not been evaluated")
def get_response_amplitudes(
self, print_phase=False)
Returns response amplitudes as pandas dataframe
def get_response_amplitudes(self,print_phase=False): """ Returns response amplitudes as pandas dataframe """ data = numpy.abs(self.response_amplitudes) index_list = self.response_names columns_list = ["Amplitude"] if print_phase: phase_angles = numpy.angle(self.response_amplitudes) data = numpy.hstack((data,phase_angles)) columns_list.append("Phase angle") df = pandas.DataFrame(data=data, index=index_list, columns=columns_list) return df
def get_response_stats(
self)
Returns response stats for time-stepping analysis carried out
def get_response_stats(self): """ Returns response stats for time-stepping analysis carried out """ if not self.run_tstep: raise ValueError("Error: no time-stepping analysis run, " + "cannot retrieve results!") # Get dict of dict of stats stats_dict_dict = self.tstep_results.CalcResponseStats() # Loop over dict and return stats as list of pandas dataframes df_list = [] for sys_obj, stats_dict in stats_dict_dict.items(): sys_name = sys_obj.name response_names = sys_obj.output_names df = pandas.DataFrame(data=stats_dict, index=[sys_name + " - " + x for x in response_names]) df_list.append(df) # Merge pandas dataframes to return a single dataframe # with all stats for all responses for all systems df_merged = pandas.concat(df_list) # Transpose such that stats as indexes, responses as columns df_merged = df_merged.transpose() return df_merged
def plot_results(
self, overlay_val=None)
Plots results of steady-state response analysis
If time-stepping analysis has also been carried out (as a cross-check) results of that analysis will be plotted also
Required:
(No arguments required)
Optional:
overlay_val
, float, value to overlay onto response plots. (A common use of this would be to overlay the applicable acceleration limit, against which responses are to be compared)
Returns:
Nested list of figure objects
def plot_results(self,overlay_val:float=None): """ Plots results of steady-state response analysis If time-stepping analysis has also been carried out (as a cross-check) results of that analysis will be plotted also *** **Required:** (No arguments required) *** **Optional:** * `overlay_val`, _float_, value to overlay onto response plots. (A common use of this would be to overlay the applicable acceleration limit, against which responses are to be compared) *** **Returns:** _Nested list_ of figure objects """ fig_list = [] if self.run_freq_method: fig = self.plot_results_freq_method(overlay_val=overlay_val) fig_list.append(fig) if self.run_tstep: tstep_fig_list = self.plot_results_tstep_method(overlay_val=overlay_val) fig_list += tstep_fig_list return fig_list
def plot_results_freq_method(
self, overlay_val=None)
Plot results of frequency domain analysis method, if avaliable
def plot_results_freq_method(self,overlay_val=None): """ Plot results of frequency domain analysis method, if avaliable """ if not self.run_freq_method: print("No results avaliable for frequency domain method\n" + "Could not plot results") return None # Produce and configure figure fig, axarr = plt.subplots(2) fig.suptitle("Steady-state pedestrian crowd loading analysis results\n"+ "Target mode: %d" % self.target_mode) fig.set_size_inches((14,8)) fig.subplots_adjust(hspace=0.3) # Produce bar plot of modal accelerations self.plot_results_modal_accn(axarr[0]) # Produce bar plot of outputs self.plot_results_responses(axarr[1],overlay_val=overlay_val) return fig
def plot_results_modal_accn(
self, ax=None)
Produce bar plot of steady-state modal accelerations
def plot_results_modal_accn(self,ax=None): """ Produce bar plot of steady-state modal accelerations """ # Define axes to plot to, if not done so already if ax is None: fig,ax = plt.subplots() fig.set_size_inches((14,4)) vals = self.modal_accn # Get or define tick labels if hasattr(self.modalsys_obj,'mode_IDs'): ticks = self.modalsys_obj.mode_IDs else: ticks = ["Mode %d" % (i+1) for i in range(self.modalsys_obj.nDOF)] ax.barh(range(len(vals)),numpy.abs(vals),tick_label=ticks) ax.tick_params(labelsize=8) ax.set_title("Steady-state modal acceleration amplitudes (m/$s^2$)")
def plot_results_responses(
self, ax=None, overlay_val=None)
Produce bar plot of steady-state responses
def plot_results_responses(self,ax=None,overlay_val=None): """ Produce bar plot of steady-state responses """ # Define axes to plot to, if not done so already if ax is None: fig,ax = plt.subplots() fig.set_size_inches((14,4)) # Get data to use in plot vals = self.response_amplitudes ticks = self.response_names # Flatten list in case of nested lists (only when multiple subsystems) if vals.shape != ticks.shape: ticks = list(itertools.chain(*self.response_names.tolist())) ax.barh(range(len(vals)),numpy.abs(vals),tick_label=ticks) ax.tick_params(labelsize=8) ax.set_title("Steady-state response amplitudes") if overlay_val is not None: ax.axvline(overlay_val,color='darkorange')
def plot_results_tstep_method(
self, overlay_val=None)
Plot results of time-stepping analysis, if avaliable
def plot_results_tstep_method(self,overlay_val=None): """ Plot results of time-stepping analysis, if avaliable """ if not self.run_tstep: print("No results avaliable for time-stepping method\n" + "Could not plot results") return [] # Define +/- values based on overlay_val, if provided if overlay_val is not None: overlay_val=[+overlay_val,-overlay_val] else: overlay_val=None fig_list = self.tstep_results.PlotResults(useCommonPlot=False, useCommonScale=False, y_overlay=overlay_val) # Overlay expected modal acceleration amplitudes # (from steady-state analysis, if run) if self.run_freq_method: state_results_fig = fig_list[0][0] ax = state_results_fig.get_axes()[3] # plot of modal accelerations accn_val = numpy.abs(self.modal_accn_target) # value to overlay ax.axhline(+accn_val,color='r',linestyle='--',alpha=0.3) ax.axhline(-accn_val,color='r',linestyle='--',alpha=0.3) # Overlay expected response amplitudes # (from steady-state analysis, if run) response_results_fig_list = fig_list[1] axarr = [x.get_axes() for x in response_results_fig_list] axarr = list(itertools.chain(*axarr)) # flatten list val2overlay = numpy.abs(self.response_amplitudes) for ax, val in zip(axarr,val2overlay): ax.axhline(+val,color='r',linestyle='--',alpha=0.3) ax.axhline(-val,color='r',linestyle='--',alpha=0.3) return fig_list
def print_details(
self)
Prints details about the analysis that has been defined
def print_details(self): """ Prints details about the analysis that has been defined """ print("\nDetails of deck regions:") deck_regions_list = self.deck_regions_list print("nRegions: %d\n" % len(deck_regions_list)) for deckregion_obj in deck_regions_list: deckregion_obj.print_details() print("Total deck area (m2):\n%.2f" % self.deck_area) print("Damped natural frequencies of system (Hz):\n{0}" .format(self.f_d))
def print_results(
self)
Prints key results from analysis
def print_results(self): """ Prints key results from analysis """ def printcomplex(vals): """ Function for printing magnitude and phase of complex-valued vals """ print("Magnitudes:\n{0}".format(numpy.abs(vals))) print("Phase angles (rad):\n{0}".format(numpy.angle(vals))) print("") print("Modal forces (N):") printcomplex(self.modal_forces) print("") if self.run_freq_method: print("Modal accelerations (m/s2):") printcomplex(self.modal_accn) print("") print("Responses:") print(self.response_names) printcomplex(self.response_amplitudes)
def run(
self, target_mode, load_intensity=None, forcing_freq=None, run_freq_method=True, run_tstep=False, tstep_kwargs={}, verbose=False, makePlot={'calc_modal_areas': False, 'calc_load_intensity': True})
Runs calculations to determine steady-state crowd loading responses
Required:
target_mode
, integer, denotes the mode at which resonance is to be targeted
Optional:
-
load_intensity
, float, denotes the load intensity (in N/m2) to be applied. If None,load_intensity
to be calculated usingcalc_load_intensity()
method function. -
forcing_freq
, float, defines forcing frequency (in Hz) at which response is to be determined. For ordinary modal systems with no added dampers this value should normally be set to equal the damped natural frequency for the mode given bytarget_mode
, as this will ensure the worst-case resonance response is obtained. This is thus the default behaviour. However is some circumstances (e.g. systems with multiple TMDs appended) it is necessary to evaluate responses at various assumed forcing frequencies, hence why this optional input is provided to the user. -
run_freq_method
, boolean, if True the amplitude of modal and real-world responses at steady-state will be determine by frequency transfer matrix method. This is the default behaviour, as it allows the steady-state amplitudes of responses to be efficiently calculated. In some cases (e.g. for systems with constraints, for which theDynSys.CalcFreqResponse()
method is not currently implemented), it may be desirable to carry out a time-stepping analysis instead. This can be done by settingrun_freq_method=False
andrun_tstep=True
(see documentation in next bullet for further details) -
run_tstep
, boolean, if True a time-stepping analysis will be run, which can be used to verify the the steady-state resonance solution that is ordinarily obtained. -
tstep_kwargs
, dict of keyword arguments to be passed when initialising time-stepping analysis using theTStep()
class. Only used ifrun_tstep
is True. Refer documentation forTStep()
class for more details. -
makePlot
, dict or bool, used to control which plots are produced
def run(self, target_mode:int, load_intensity:float=None, forcing_freq:float=None, run_freq_method=True, run_tstep=False,tstep_kwargs:dict={}, verbose=False, makePlot={'calc_modal_areas':False,'calc_load_intensity':True} ): """ Runs calculations to determine steady-state crowd loading responses *** **Required:** * `target_mode`, _integer_, denotes the mode at which resonance is to be targeted *** **Optional:** * `load_intensity`, _float_, denotes the load intensity (in N/m2) to be applied. If _None_, `load_intensity` to be calculated using `calc_load_intensity()` method function. * `forcing_freq`, _float_, defines forcing frequency (in Hz) at which response is to be determined. For ordinary modal systems with no added dampers this value should normally be set to equal the damped natural frequency for the mode given by `target_mode`, as this will ensure the worst-case resonance response is obtained. This is thus the default behaviour. However is some circumstances (e.g. systems with multiple TMDs appended) it is necessary to evaluate responses at various assumed forcing frequencies, hence why this optional input is provided to the user. * `run_freq_method`, _boolean_, if True the amplitude of modal and real-world responses at steady-state will be determine by frequency transfer matrix method. This is the default behaviour, as it allows the steady-state amplitudes of responses to be efficiently calculated. In some cases (e.g. for systems with constraints, for which the `DynSys.CalcFreqResponse()` method is not currently implemented), it may be desirable to carry out a time-stepping analysis instead. This can be done by setting `run_freq_method=False` and `run_tstep=True` (see documentation in next bullet for further details) * `run_tstep`, _boolean_, if True a time-stepping analysis will be run, which can be used to verify the the steady-state resonance solution that is ordinarily obtained. * `tstep_kwargs`, _dict_ of keyword arguments to be passed when initialising time-stepping analysis using the `TStep()` class. Only used if `run_tstep` is True. Refer [documentation](../docs/tstep.html) for `TStep()` class for more details. * `makePlot`, _dict_ or _bool_, used to control which plots are produced [documentation]: """ if verbose: print("\n**** Running steady-state crowd loading analysis ****") # Handle makePlot if isinstance(makePlot,bool): makePlot_dict = {} makePlot_dict['calc_modal_areas']=makePlot makePlot_dict['calc_load_intensity']=makePlot else: makePlot_dict = makePlot self.target_mode = target_mode if verbose: print("Target mode index: %d" % target_mode) # Define forcing frequency if forcing_freq is None: # Set to be damped natural frequency of mode being targeted forcing_freq = numpy.abs(self.f_d[target_mode]) else: # where input value provided via optional arg # Check float type is provided if not isinstance(forcing_freq,float): raise ValueError("Error: float expected for `forcing_freq`") self.forcing_freq = forcing_freq # Obtain area integrals over all deck strips modal_areas = self.calc_modal_areas(target_mode, makePlot=makePlot_dict['calc_modal_areas']) self.modal_areas = modal_areas # Calculate load intensity, if not provided directly if load_intensity is None: load_intensity, d = self.calc_load_intensity(mode_index=target_mode, fv=forcing_freq, verbose=verbose, makePlot=makePlot_dict['calc_load_intensity']) else: d = {} # null dict self.load_intensity = load_intensity """ Uniform load intensity (in N/m2) to be applied to deck regions """ self.load_intensity_rslts_dict = d """ _Dict_ to store key results from load intensity calculation """ if verbose: print("UDL intensity (N/m2): %.2f" % load_intensity) # Multiply by load intensity to obtain modal force amplitudes # Note modal force amplitudes should be regarded as complex variables # (although here they will be real-valued) whereby negative values # imply a force which is 180 degrees out of phase modal_forces = self.load_intensity * modal_areas self.modal_forces = modal_forces # Save settings used to determine which analyses run self.run_freq_method = run_freq_method self.run_tstep = run_tstep if verbose: print("Forcing frequency: %.2f Hz" % forcing_freq) # Run frequency domain method to determine steady-state responses if run_freq_method: # Evaluate frequency transfer functions G(f) relating applied forces # (i.e. modal forces in this application) to DOF accelerations, # displacements and velocities # (i.e. modal disp, vel, accn in this application) # # It is only necessary to evaluate G(f) at resonant frequency for mode # being targeted, although for validation / plotting it may be # advantageous to evaluate G(f) over a wider frequency range # # By definition, G(f1) is a matrix that can be used to relate # sinusoidal modal forces to sinusoidal modal accelerations, # given forcing frequency f1. # # Note: G(f) will in general be complex-valued; the should be # intepreted to denote scaling and phase lag # Calculate matrix mapping applied loads to state DOFs def Calc_Gf(f): rslts = self.modalsys_obj.CalcFreqResponse(fVals=f, C=numpy.zeros((0,)), verbose=verbose) f = rslts["f"] G_f = rslts["G_f"] return f, G_f G_state_vars = Calc_Gf(forcing_freq)[1] G_state_vars = G_state_vars[0,:,:] # Retain only submatrix mapping modal forces to state DOFs nModes = len(modal_forces) G_state_vars = G_state_vars[:,:nModes] self.G_state_vars = G_state_vars # Get submatrix mapping applied forces to DOF accelerations nDOF = int(G_state_vars.shape[0]/3) G_acc = G_state_vars[2*nDOF:,:] # Get submatrix mapping modal forces to m accelerations G_modalacc = G_acc[:nModes,:nModes] # Calculate modal acceleration amplitudes using transfer matrix modal_accn = numpy.ravel(G_modalacc @ modal_forces) self.modal_accn = modal_accn # Adjust phase such that acc for target mode has phase=0 z1 = 1.0 * numpy.exp(1j*numpy.angle(modal_accn[target_mode])) modal_accn = modal_accn / z1 # Get modal acceleration of targeted mode # Note phase=0 due to above adjustment, hence ok to take real part only modal_accn_target = numpy.real(modal_accn[target_mode]) self.modal_accn_target = modal_accn_target # Calculate responses amplitudes using full transfer matrix self.calc_responses() if verbose: self.print_results() print("**** Analysis complete! ****\n") if run_tstep: # Run time-stepping analysis to validate /check the above approach print("***Running time-stepping analysis, as cross-check of "+ "steady-state resonance approach...***") # Define sinusoidal load function consistent with above def load_func(t): return modal_forces * numpy.sin(2*numpy.pi*forcing_freq*t) load_func_dict = {self.modalsys_obj : load_func} # Define and run time-stepping analysis tstep_analysis = tstep.TStep(dynsys_obj=self.modalsys_obj, force_func_dict=load_func_dict, **tstep_kwargs) tstep_results = tstep_analysis.run(verbose=verbose) self.tstep_results = tstep_results """ Results from time-stepping analysis; instance of `TStep_Results()` """ print("***Time-stepping analysis complete!***")
Instance variables
var deck_area
Total deck area ($m^2$), obtained by integration over all deck strips
var deck_regions_list
List of DeckStrip classes, defining deck strips
var eta
Damping ratios for system modes
var f_d
Damped natural frequencies of system being analysed
var max_modeshape
Maximum |modeshape| ordinates, evaluated over all deck regions. Array of length equal to number of modes considered
var modalsys_obj
ModalSys
object, corresponding to system being analysed
var nRegions
Integer, denotes number deck regions
var name
Name assigned to class instance
var response_amplitudes
Amplitudes of defined responses
var response_names
Names of responses computed
class UKNA_BSEN1991_2_crowd
Class to implement steady state crowd loading analysis to UK NA to BS EN1991-2
class UKNA_BSEN1991_2_crowd(SteadyStateCrowdLoading): """ Class to implement steady state crowd loading analysis to UK NA to BS EN1991-2 """ def __init__(self, bridgeClass:str=None, crowd_density:float=None, load_direction='vertical', **kwargs): """ Initialisation function *** Optional: * `bridgeClass`, string character ('A' to 'D'), defines bridge class per UK NA to BS EN1991-2 * `crowd_density`, float, defines crowd density (P/m2). Note only used if `bridgeClass=None` * `load_direction`, string, either 'vertical' or 'lateral' required """ if bridgeClass is None: if crowd_density is None: raise ValueError("Either `bridgeClass` or `crowd_density` required") else: if not isinstance(bridgeClass,str): raise ValueError("`bridgeClass` to be string character") else: # correct to uppercase bridgeClass = bridgeClass.upper() self.bridgeClass = bridgeClass """ Bridge class per Table NA.7 """ self.crowd_density = crowd_density """ Crowd density (P/m2) """ self.load_direction = load_direction """ Direction of applied crowd loading """ # Run parent init function super().__init__(**kwargs) def get_crowd_density(self, verbose=True, saveAsAttr=True): """ Get crowd density according from Table NA.7, UK NA to BS EN 1991-2, according to `bridgeClass` Returns density expressed as persons/m2 """ # Define crowd density via bridge class, if defined if self.bridgeClass is not None: bridgeClass = self.bridgeClass if bridgeClass == 'A': density=0.0 elif bridgeClass == 'B': density=0.4 elif bridgeClass == 'C': density=0.8 elif bridgeClass == 'D': density=1.5 else: raise ValueError("Invalid 'bridgeClass'!") # Otherwise use pre-defined density else: density = self.crowd_density self.crowd_density = density if verbose: print("Bridge class: '%s'" % bridgeClass) print("Crowd density (persons/m2): %.1f" % density) if density==0: raise ValueError("Crowd density = 0; no analysis required!") return density def calc_load_intensity(self, mode_index:int, fv:float=None, load_direction=None, calc_lambda=True, verbose=True, makePlot=True): """ Calculates required load intensity (N/m2) to UK NA to BS EN 1991-2 ![Cl_2_44_5](../dynsys/img/UKNA_BSEN1991_2_Cl_2_44_5.png) """ if load_direction is None: load_direction = self.load_direction # Retrieve attributes A = self.deck_area max_modeshape = self.max_modeshape[mode_index] # Get crowd density according to bridgeClass # Crowd density, expressed in persons/m2 crowd_density = self.get_crowd_density(verbose=False) # Get total number of pedestrians N = crowd_density * A # Define code inputs if load_direction == 'vertical': F0 = 280.0 # N, refer Table NA.8 # n.b 'walkers' value to be used in crowd loading calcs else: F0 = 70.0 # N, refer Tue 24/07/2018 12:56 Chris Barker email # Derive adjustment factors if fv is None: # If not provided explicitly set to damped frequency of target mode fv = numpy.abs(self.f_d[mode_index]) k = UKNA_BSEN1991_2_Figure_NA_8(fv=fv, analysis_type="walkers", load_direction=load_direction, makePlot=makePlot) log_dec = 2*numpy.pi*self.eta[mode_index] gamma = UKNA_BSEN1991_2_Figure_NA_9(logDec=log_dec, groupType="crowd", makePlot=makePlot) # Effective number of pedestrians parameter, refer NA. 2.44.5(1) # To generalise to deck of variable width, use area ratio as Seff/S if calc_lambda: Aeff = self.modal_areas[mode_index] / (0.634 * max_modeshape) lambda_val = 0.634*Aeff/A else: lambda_val = 1.0 # conservative according to text - but not always! # Calculate load intensity per NA.2.44.5(1) load_intensity = 1.8*(F0/A)*k*((gamma*N/lambda_val)**0.5) # Prepare markdown string of nicely formatted text md_txt = "" md_txt += "Key results from UK NA to BS EN 1991-2 load intensity calculation:" md_txt += "\n" + "F0 = %.1f\t(N)" % F0 md_txt += "\n" + "A = %.2f\t(m2)" % A md_txt += "\n" + "gamma_max = %.3f\t" % max_modeshape md_txt += "\n" + "Aeff = %.2f\t(m2)" % Aeff md_txt += "\n" + "rho = %.2f\t(P/m2)" % crowd_density md_txt += "\n" + "N = %.1f" % N md_txt += "\n" + "fv = %.3f\t(Hz)" % fv md_txt += "\n" + "k = %.3f" % k md_txt += "\n" + "log_dec = %.3f" % log_dec md_txt += "\n" + "gamma = %.3f" % gamma md_txt += "\n" + "lambda = %.3f" % lambda_val md_txt += "\n" + "w = %.3f\t(N/m2)" % load_intensity if verbose: print("\n" + md_txt + "\n") # Prepare dict to return key results results_dict = {} results_dict["F0"]=F0 results_dict["fv"]=fv results_dict["rho"]=crowd_density results_dict["N"]=N results_dict["k"]=k results_dict["lambda"]=lambda_val results_dict["log_dec"]=log_dec results_dict["w"]=load_intensity results_dict["md_txt"]=md_txt return load_intensity, results_dict
Ancestors (in MRO)
- UKNA_BSEN1991_2_crowd
- SteadyStateCrowdLoading
- builtins.object
Static methods
def __init__(
self, bridgeClass=None, crowd_density=None, load_direction='vertical', **kwargs)
Initialisation function
Optional:
-
bridgeClass
, string character ('A' to 'D'), defines bridge class per UK NA to BS EN1991-2 -
crowd_density
, float, defines crowd density (P/m2). Note only used ifbridgeClass=None
-
load_direction
, string, either 'vertical' or 'lateral' required
def __init__(self, bridgeClass:str=None, crowd_density:float=None, load_direction='vertical', **kwargs): """ Initialisation function *** Optional: * `bridgeClass`, string character ('A' to 'D'), defines bridge class per UK NA to BS EN1991-2 * `crowd_density`, float, defines crowd density (P/m2). Note only used if `bridgeClass=None` * `load_direction`, string, either 'vertical' or 'lateral' required """ if bridgeClass is None: if crowd_density is None: raise ValueError("Either `bridgeClass` or `crowd_density` required") else: if not isinstance(bridgeClass,str): raise ValueError("`bridgeClass` to be string character") else: # correct to uppercase bridgeClass = bridgeClass.upper() self.bridgeClass = bridgeClass """ Bridge class per Table NA.7 """ self.crowd_density = crowd_density """ Crowd density (P/m2) """ self.load_direction = load_direction """ Direction of applied crowd loading """ # Run parent init function super().__init__(**kwargs)
def calc_deck_area(
self)
Returns the total deck area, taken across all deck strips
def calc_deck_area(self): """ Returns the total deck area, taken across all deck strips """ total_area = sum([x.calc_area() for x in self.deck_regions_list]) return float(total_area)
def calc_load_intensity(
self, mode_index, fv=None, load_direction=None, calc_lambda=True, verbose=True, makePlot=True)
Calculates required load intensity (N/m2) to UK NA to BS EN 1991-2
def calc_load_intensity(self, mode_index:int, fv:float=None, load_direction=None, calc_lambda=True, verbose=True, makePlot=True): """ Calculates required load intensity (N/m2) to UK NA to BS EN 1991-2 ![Cl_2_44_5](../dynsys/img/UKNA_BSEN1991_2_Cl_2_44_5.png) """ if load_direction is None: load_direction = self.load_direction # Retrieve attributes A = self.deck_area max_modeshape = self.max_modeshape[mode_index] # Get crowd density according to bridgeClass # Crowd density, expressed in persons/m2 crowd_density = self.get_crowd_density(verbose=False) # Get total number of pedestrians N = crowd_density * A # Define code inputs if load_direction == 'vertical': F0 = 280.0 # N, refer Table NA.8 # n.b 'walkers' value to be used in crowd loading calcs else: F0 = 70.0 # N, refer Tue 24/07/2018 12:56 Chris Barker email # Derive adjustment factors if fv is None: # If not provided explicitly set to damped frequency of target mode fv = numpy.abs(self.f_d[mode_index]) k = UKNA_BSEN1991_2_Figure_NA_8(fv=fv, analysis_type="walkers", load_direction=load_direction, makePlot=makePlot) log_dec = 2*numpy.pi*self.eta[mode_index] gamma = UKNA_BSEN1991_2_Figure_NA_9(logDec=log_dec, groupType="crowd", makePlot=makePlot) # Effective number of pedestrians parameter, refer NA. 2.44.5(1) # To generalise to deck of variable width, use area ratio as Seff/S if calc_lambda: Aeff = self.modal_areas[mode_index] / (0.634 * max_modeshape) lambda_val = 0.634*Aeff/A else: lambda_val = 1.0 # conservative according to text - but not always! # Calculate load intensity per NA.2.44.5(1) load_intensity = 1.8*(F0/A)*k*((gamma*N/lambda_val)**0.5) # Prepare markdown string of nicely formatted text md_txt = "" md_txt += "Key results from UK NA to BS EN 1991-2 load intensity calculation:" md_txt += "\n" + "F0 = %.1f\t(N)" % F0 md_txt += "\n" + "A = %.2f\t(m2)" % A md_txt += "\n" + "gamma_max = %.3f\t" % max_modeshape md_txt += "\n" + "Aeff = %.2f\t(m2)" % Aeff md_txt += "\n" + "rho = %.2f\t(P/m2)" % crowd_density md_txt += "\n" + "N = %.1f" % N md_txt += "\n" + "fv = %.3f\t(Hz)" % fv md_txt += "\n" + "k = %.3f" % k md_txt += "\n" + "log_dec = %.3f" % log_dec md_txt += "\n" + "gamma = %.3f" % gamma md_txt += "\n" + "lambda = %.3f" % lambda_val md_txt += "\n" + "w = %.3f\t(N/m2)" % load_intensity if verbose: print("\n" + md_txt + "\n") # Prepare dict to return key results results_dict = {} results_dict["F0"]=F0 results_dict["fv"]=fv results_dict["rho"]=crowd_density results_dict["N"]=N results_dict["k"]=k results_dict["lambda"]=lambda_val results_dict["log_dec"]=log_dec results_dict["w"]=load_intensity results_dict["md_txt"]=md_txt return load_intensity, results_dict
def calc_modal_areas(
self, mode_index, makePlot=False)
Returns modal area integrals
def calc_modal_areas(self,mode_index,makePlot=False): """ Returns modal area integrals """ return numpy.sum([x.integrate(mode_index,makePlot=makePlot) for x in self.deck_regions_list],axis=0)
def calc_responses(
self, output_mtrx=None, output_names=None)
Calculate steady-steady amplitude response of defined responses
def calc_responses(self,output_mtrx=None,output_names=None): """ Calculate steady-steady amplitude response of defined responses """ # Get output matrix and output names, if not passed to this function if output_mtrx is None: output_mtrx, output_names = self.modalsys_obj.GetOutputMtrx(all_systems=True) # Get frequency response matrix for {disp,vel,accn} variables # (note: assumes already calculated!) G_state_vars = self.G_state_vars # Get modal forces modal_forces = self.modal_forces # Calculate responses at steady-state nResponses = output_mtrx.shape[0] if nResponses>0: response_amplitudes = numpy.ravel(output_mtrx @ G_state_vars @ modal_forces) self.response_amplitudes = response_amplitudes self.response_names = output_names else: print("No output matrix provided. " + "Real-world responses have not been evaluated")
def get_crowd_density(
self, verbose=True, saveAsAttr=True)
Get crowd density according from Table NA.7, UK NA to BS EN
1991-2, according to bridgeClass
Returns density expressed as persons/m2
def get_crowd_density(self, verbose=True, saveAsAttr=True): """ Get crowd density according from Table NA.7, UK NA to BS EN 1991-2, according to `bridgeClass` Returns density expressed as persons/m2 """ # Define crowd density via bridge class, if defined if self.bridgeClass is not None: bridgeClass = self.bridgeClass if bridgeClass == 'A': density=0.0 elif bridgeClass == 'B': density=0.4 elif bridgeClass == 'C': density=0.8 elif bridgeClass == 'D': density=1.5 else: raise ValueError("Invalid 'bridgeClass'!") # Otherwise use pre-defined density else: density = self.crowd_density self.crowd_density = density if verbose: print("Bridge class: '%s'" % bridgeClass) print("Crowd density (persons/m2): %.1f" % density) if density==0: raise ValueError("Crowd density = 0; no analysis required!") return density
def get_response_amplitudes(
self, print_phase=False)
Returns response amplitudes as pandas dataframe
def get_response_amplitudes(self,print_phase=False): """ Returns response amplitudes as pandas dataframe """ data = numpy.abs(self.response_amplitudes) index_list = self.response_names columns_list = ["Amplitude"] if print_phase: phase_angles = numpy.angle(self.response_amplitudes) data = numpy.hstack((data,phase_angles)) columns_list.append("Phase angle") df = pandas.DataFrame(data=data, index=index_list, columns=columns_list) return df
def get_response_stats(
self)
Returns response stats for time-stepping analysis carried out
def get_response_stats(self): """ Returns response stats for time-stepping analysis carried out """ if not self.run_tstep: raise ValueError("Error: no time-stepping analysis run, " + "cannot retrieve results!") # Get dict of dict of stats stats_dict_dict = self.tstep_results.CalcResponseStats() # Loop over dict and return stats as list of pandas dataframes df_list = [] for sys_obj, stats_dict in stats_dict_dict.items(): sys_name = sys_obj.name response_names = sys_obj.output_names df = pandas.DataFrame(data=stats_dict, index=[sys_name + " - " + x for x in response_names]) df_list.append(df) # Merge pandas dataframes to return a single dataframe # with all stats for all responses for all systems df_merged = pandas.concat(df_list) # Transpose such that stats as indexes, responses as columns df_merged = df_merged.transpose() return df_merged
def plot_results(
self, overlay_val=None)
Plots results of steady-state response analysis
If time-stepping analysis has also been carried out (as a cross-check) results of that analysis will be plotted also
Required:
(No arguments required)
Optional:
overlay_val
, float, value to overlay onto response plots. (A common use of this would be to overlay the applicable acceleration limit, against which responses are to be compared)
Returns:
Nested list of figure objects
def plot_results(self,overlay_val:float=None): """ Plots results of steady-state response analysis If time-stepping analysis has also been carried out (as a cross-check) results of that analysis will be plotted also *** **Required:** (No arguments required) *** **Optional:** * `overlay_val`, _float_, value to overlay onto response plots. (A common use of this would be to overlay the applicable acceleration limit, against which responses are to be compared) *** **Returns:** _Nested list_ of figure objects """ fig_list = [] if self.run_freq_method: fig = self.plot_results_freq_method(overlay_val=overlay_val) fig_list.append(fig) if self.run_tstep: tstep_fig_list = self.plot_results_tstep_method(overlay_val=overlay_val) fig_list += tstep_fig_list return fig_list
def plot_results_freq_method(
self, overlay_val=None)
Plot results of frequency domain analysis method, if avaliable
def plot_results_freq_method(self,overlay_val=None): """ Plot results of frequency domain analysis method, if avaliable """ if not self.run_freq_method: print("No results avaliable for frequency domain method\n" + "Could not plot results") return None # Produce and configure figure fig, axarr = plt.subplots(2) fig.suptitle("Steady-state pedestrian crowd loading analysis results\n"+ "Target mode: %d" % self.target_mode) fig.set_size_inches((14,8)) fig.subplots_adjust(hspace=0.3) # Produce bar plot of modal accelerations self.plot_results_modal_accn(axarr[0]) # Produce bar plot of outputs self.plot_results_responses(axarr[1],overlay_val=overlay_val) return fig
def plot_results_modal_accn(
self, ax=None)
Produce bar plot of steady-state modal accelerations
def plot_results_modal_accn(self,ax=None): """ Produce bar plot of steady-state modal accelerations """ # Define axes to plot to, if not done so already if ax is None: fig,ax = plt.subplots() fig.set_size_inches((14,4)) vals = self.modal_accn # Get or define tick labels if hasattr(self.modalsys_obj,'mode_IDs'): ticks = self.modalsys_obj.mode_IDs else: ticks = ["Mode %d" % (i+1) for i in range(self.modalsys_obj.nDOF)] ax.barh(range(len(vals)),numpy.abs(vals),tick_label=ticks) ax.tick_params(labelsize=8) ax.set_title("Steady-state modal acceleration amplitudes (m/$s^2$)")
def plot_results_responses(
self, ax=None, overlay_val=None)
Produce bar plot of steady-state responses
def plot_results_responses(self,ax=None,overlay_val=None): """ Produce bar plot of steady-state responses """ # Define axes to plot to, if not done so already if ax is None: fig,ax = plt.subplots() fig.set_size_inches((14,4)) # Get data to use in plot vals = self.response_amplitudes ticks = self.response_names # Flatten list in case of nested lists (only when multiple subsystems) if vals.shape != ticks.shape: ticks = list(itertools.chain(*self.response_names.tolist())) ax.barh(range(len(vals)),numpy.abs(vals),tick_label=ticks) ax.tick_params(labelsize=8) ax.set_title("Steady-state response amplitudes") if overlay_val is not None: ax.axvline(overlay_val,color='darkorange')
def plot_results_tstep_method(
self, overlay_val=None)
Plot results of time-stepping analysis, if avaliable
def plot_results_tstep_method(self,overlay_val=None): """ Plot results of time-stepping analysis, if avaliable """ if not self.run_tstep: print("No results avaliable for time-stepping method\n" + "Could not plot results") return [] # Define +/- values based on overlay_val, if provided if overlay_val is not None: overlay_val=[+overlay_val,-overlay_val] else: overlay_val=None fig_list = self.tstep_results.PlotResults(useCommonPlot=False, useCommonScale=False, y_overlay=overlay_val) # Overlay expected modal acceleration amplitudes # (from steady-state analysis, if run) if self.run_freq_method: state_results_fig = fig_list[0][0] ax = state_results_fig.get_axes()[3] # plot of modal accelerations accn_val = numpy.abs(self.modal_accn_target) # value to overlay ax.axhline(+accn_val,color='r',linestyle='--',alpha=0.3) ax.axhline(-accn_val,color='r',linestyle='--',alpha=0.3) # Overlay expected response amplitudes # (from steady-state analysis, if run) response_results_fig_list = fig_list[1] axarr = [x.get_axes() for x in response_results_fig_list] axarr = list(itertools.chain(*axarr)) # flatten list val2overlay = numpy.abs(self.response_amplitudes) for ax, val in zip(axarr,val2overlay): ax.axhline(+val,color='r',linestyle='--',alpha=0.3) ax.axhline(-val,color='r',linestyle='--',alpha=0.3) return fig_list
def print_details(
self)
Prints details about the analysis that has been defined
def print_details(self): """ Prints details about the analysis that has been defined """ print("\nDetails of deck regions:") deck_regions_list = self.deck_regions_list print("nRegions: %d\n" % len(deck_regions_list)) for deckregion_obj in deck_regions_list: deckregion_obj.print_details() print("Total deck area (m2):\n%.2f" % self.deck_area) print("Damped natural frequencies of system (Hz):\n{0}" .format(self.f_d))
def print_results(
self)
Prints key results from analysis
def print_results(self): """ Prints key results from analysis """ def printcomplex(vals): """ Function for printing magnitude and phase of complex-valued vals """ print("Magnitudes:\n{0}".format(numpy.abs(vals))) print("Phase angles (rad):\n{0}".format(numpy.angle(vals))) print("") print("Modal forces (N):") printcomplex(self.modal_forces) print("") if self.run_freq_method: print("Modal accelerations (m/s2):") printcomplex(self.modal_accn) print("") print("Responses:") print(self.response_names) printcomplex(self.response_amplitudes)
def run(
self, target_mode, load_intensity=None, forcing_freq=None, run_freq_method=True, run_tstep=False, tstep_kwargs={}, verbose=False, makePlot={'calc_modal_areas': False, 'calc_load_intensity': True})
Runs calculations to determine steady-state crowd loading responses
Required:
target_mode
, integer, denotes the mode at which resonance is to be targeted
Optional:
-
load_intensity
, float, denotes the load intensity (in N/m2) to be applied. If None,load_intensity
to be calculated usingcalc_load_intensity()
method function. -
forcing_freq
, float, defines forcing frequency (in Hz) at which response is to be determined. For ordinary modal systems with no added dampers this value should normally be set to equal the damped natural frequency for the mode given bytarget_mode
, as this will ensure the worst-case resonance response is obtained. This is thus the default behaviour. However is some circumstances (e.g. systems with multiple TMDs appended) it is necessary to evaluate responses at various assumed forcing frequencies, hence why this optional input is provided to the user. -
run_freq_method
, boolean, if True the amplitude of modal and real-world responses at steady-state will be determine by frequency transfer matrix method. This is the default behaviour, as it allows the steady-state amplitudes of responses to be efficiently calculated. In some cases (e.g. for systems with constraints, for which theDynSys.CalcFreqResponse()
method is not currently implemented), it may be desirable to carry out a time-stepping analysis instead. This can be done by settingrun_freq_method=False
andrun_tstep=True
(see documentation in next bullet for further details) -
run_tstep
, boolean, if True a time-stepping analysis will be run, which can be used to verify the the steady-state resonance solution that is ordinarily obtained. -
tstep_kwargs
, dict of keyword arguments to be passed when initialising time-stepping analysis using theTStep()
class. Only used ifrun_tstep
is True. Refer documentation forTStep()
class for more details. -
makePlot
, dict or bool, used to control which plots are produced