Source code for openquake.hazardlib.gsim.chiou_youngs_2008
# -*- coding: utf-8 -*-# vim: tabstop=4 shiftwidth=4 softtabstop=4## Copyright (C) 2012-2023 GEM Foundation## OpenQuake is free software: you can redistribute it and/or modify it# under the terms of the GNU Affero General Public License as published# by the Free Software Foundation, either version 3 of the License, or# (at your option) any later version.## OpenQuake is distributed in the hope that it will be useful,# but WITHOUT ANY WARRANTY; without even the implied warranty of# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the# GNU Affero General Public License for more details.## You should have received a copy of the GNU Affero General Public License# along with OpenQuake. If not, see <http://www.gnu.org/licenses/>."""Module exports :class:`ChiouYoungs2008`."""importnumpyasnpfromopenquake.hazardlib.gsim.baseimportGMPE,CoeffsTablefromopenquake.hazardlibimportconstfromopenquake.hazardlib.imtimportPGA,PGV,SAdef_get_ln_y_ref(ctx,C):""" Get an intensity on a reference soil. Implements eq. 13a. """# reverse faulting flagFrv=np.zeros_like(ctx.mag)Frv[(30<=ctx.rake)&(ctx.rake<=150)]=1# normal faulting flagFnm=np.zeros_like(ctx.mag)Fnm[(-120<=ctx.rake)&(ctx.rake<=-60)]=1# hanging wall flagFhw=np.where(ctx.rx>=0,1,0)# aftershock flag. always zero since we only consider main shockAS=0ln_y_ref=(# first line of eq. 13aC['c1']+(C['c1a']*Frv+C['c1b']*Fnm+C['c7']*(ctx.ztor-4))*(1-AS)+(C['c10']+C['c7a']*(ctx.ztor-4))*AS# second line+C['c2']*(ctx.mag-6)+((C['c2']-C['c3'])/C['cn'])*np.log(1+np.exp(C['cn']*(C['cm']-ctx.mag)))# third line+C['c4']*np.log(ctx.rrup+C['c5']*np.cosh(C['c6']*np.fmax(ctx.mag-C['chm'],0)))# fourth line+(C['c4a']-C['c4'])*np.log(np.sqrt(ctx.rrup**2+C['crb']**2))# fifth line+(C['cg1']+C['cg2']/(np.cosh(np.fmax(ctx.mag-C['cg3'],0))))*ctx.rrup# sixth line+C['c9']*Fhw*np.tanh(ctx.rx*(np.cos(np.radians(ctx.dip))**2)/C['c9a'])*(1-np.sqrt(ctx.rjb**2+ctx.ztor**2)/(ctx.rrup+0.001)))returnln_y_refdef_get_mean(ctx,C,ln_y_ref,exp1,exp2):""" Add site effects to an intensity. Implements eq. 13b. """# we do not support estimating of basin depth and instead# rely on it being available (since we require it).z1pt0=ctx.z1pt0# we consider random variables being zero since we want# to find the exact mean value.eta=epsilon=0ln_y=(# first line of eq. 13bln_y_ref+C['phi1']*np.log(ctx.vs30/1130).clip(-np.inf,0)# second line+C['phi2']*(exp1-exp2)*np.log((np.exp(ln_y_ref)+C['phi4'])/C['phi4'])# third line+C['phi5']*(1.0-1.0/np.cosh(C['phi6']*(z1pt0-C['phi7']).clip(0,np.inf)))+C['phi8']/np.cosh(0.15*(z1pt0-15).clip(0,np.inf))# fourth line+eta+epsilon)returnln_ydef_set_stddevs(sig,tau,phi,ctx,C,ln_y_ref,exp1,exp2):""" Get standard deviation for a given intensity on reference soil. Implements equations 19, 20 and 21 for inter-event, intra-event and total standard deviations respectively. """# aftershock flag is zero, we consider only main shock.AS=0Fmeasured=np.where(ctx.vs30measured,True,False)Finferred=~Fmeasured# eq. 19 to calculate inter-event standard errormag_test=np.clip(ctx.mag,5.0,7.0)-5.0t=C['tau1']+(C['tau2']-C['tau1'])/2*mag_test# b and c coeffs from eq. 10b=C['phi2']*(exp1-exp2)c=C['phi4']y_ref=np.exp(ln_y_ref)# eq. 20NL=b*y_ref/(y_ref+c)sigma=(# first line of eq. 20(C['sig1']+0.5*(C['sig2']-C['sig1'])*mag_test+C['sig4']*AS)# second line*np.sqrt((C['sig3']*Finferred+0.7*Fmeasured)+(1+NL)**2))# eq. 21sig[:]=np.sqrt(((1+NL)**2)*(t**2)+sigma**2)tau[:]=np.abs((1+NL)*t)phi[:]=sigma
[docs]classChiouYoungs2008(GMPE):""" Implements GMPE developed by Brian S.-J. Chiou and Robert R. Youngs and published as "An NGA Model for the Average Horizontal Component of Peak Ground Motion and Response Spectra" (2008, Earthquake Spectra, Volume 24, No. 1, pages 173-215). """#: Supported tectonic region type is active shallow crust, see page 174.DEFINED_FOR_TECTONIC_REGION_TYPE=const.TRT.ACTIVE_SHALLOW_CRUST#: Supported intensity measure types are spectral acceleration,#: peak ground velocity and peak ground acceleration, see tables#: at pages 198 and 199.DEFINED_FOR_INTENSITY_MEASURE_TYPES={PGA,PGV,SA}#: Supported intensity measure component is orientation-independent#: measure :attr:`~openquake.hazardlib.const.IMC.GMRotI50`, see page 174.DEFINED_FOR_INTENSITY_MEASURE_COMPONENT=const.IMC.GMRotI50#: Supported standard deviation types are inter-event, intra-event#: and total, see chapter "Variance model".DEFINED_FOR_STANDARD_DEVIATION_TYPES={const.StdDev.TOTAL,const.StdDev.INTER_EVENT,const.StdDev.INTRA_EVENT}DEFINED_FOR_REFERENCE_VELOCITY=1130.#: Required site parameters are Vs30 (eq. 13b), Vs30 measured flag (eq. 20)#: and Z1.0 (eq. 13b).REQUIRES_SITES_PARAMETERS={'vs30','vs30measured','z1pt0'}#: Required rupture parameters are magnitude, rake (eq. 13a and 13b),#: dip (eq. 13a) and ztor (eq. 13a).REQUIRES_RUPTURE_PARAMETERS={'dip','rake','mag','ztor'}#: Required distance measures are RRup, Rjb and Rx (all are in eq. 13a).REQUIRES_DISTANCES={'rrup','rjb','rx'}
[docs]defcompute(self,ctx:np.recarray,imts,mean,sig,tau,phi):""" See :meth:`superclass method <.base.GroundShakingIntensityModel.compute>` for spec of input and result values. """form,imtinenumerate(imts):# extracting dictionary of coefficients specific to required# intensity measure type.C=self.COEFFS[imt]# intensity on a reference soil is used for both mean# and stddev calculations.ln_y_ref=_get_ln_y_ref(ctx,C)# exp1 and exp2 are parts of eq. 10 and eq. 13b,# calculate it once for both.exp1=np.exp(C['phi3']*(ctx.vs30.clip(-np.inf,1130)-360))exp2=np.exp(C['phi3']*(1130-360))mean[m]=_get_mean(ctx,C,ln_y_ref,exp1,exp2)_set_stddevs(sig[m],tau[m],phi[m],ctx,C,ln_y_ref,exp1,exp2)
