Source code for openquake.hmtk.faults.mfd.characteristic

# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4

#
# LICENSE
#
# Copyright (C) 2010-2023 GEM Foundation, G. Weatherill, M. Pagani,
# D. Monelli.
#
# The Hazard Modeller's Toolkit 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.
#
# You should have received a copy of the GNU Affero General Public License
# along with OpenQuake. If not, see <http://www.gnu.org/licenses/>
#
# DISCLAIMER
#
# The software Hazard Modeller's Toolkit (openquake.hmtk) provided herein
# is released as a prototype implementation on behalf of
# scientists and engineers working within the GEM Foundation (Global
# Earthquake Model).
#
# It is distributed for the purpose of open collaboration and in the
# hope that it will be useful to the scientific, engineering, disaster
# risk and software design communities.
#
# The software is NOT distributed as part of GEM's OpenQuake suite
# (https://www.globalquakemodel.org/tools-products) and must be considered as a
# separate entity. The software provided herein is designed and implemented
# by scientific staff. It is not developed to the design standards, nor
# subject to same level of critical review by professional software
# developers, as GEM's OpenQuake software suite.
#
# Feedback and contribution to the software is welcome, and can be
# directed to the hazard scientific staff of the GEM Model Facility
# (hazard@globalquakemodel.org).
#
# The Hazard Modeller's Toolkit (openquake.hmtk) is therefore distributed
# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
# for more details.
#
# The GEM Foundation, and the authors of the software, assume no
# liability for use of the software.

"""
Module :mod: openquake.hmtk.faults.mfd.characteristic implements
:class:Characteristic the simple characteristic earthquake calculator
of recurrence.
"""

import numpy as np
from scipy.stats import truncnorm
from math import fabs
from openquake.hmtk.faults.mfd.base import _scale_moment, BaseMFDfromSlip


[docs]class Characteristic(BaseMFDfromSlip): """ Class to implement the characteristic earthquake model assuming a truncated Gaussian distribution :param str mfd_model: Type of magnitude frequency distribution :param float mfd_weight: Weight of the mfd distribution (for subsequent logic tree processing) :param float bin_width: Width of the magnitude bin (rates are given for the centre point) :param float mmin: Minimum magnitude :param float mmax: Maximum magnitude :param float mmax_sigma: Uncertainty on maximum magnitude :param float lower_bound: Lower bound of Gaussian distribution (as number of standard deviations) :param float upper_bound: Upper bound of Gaussian distribution (as number of standard deviations) :param float sigma: Standard deviation (in magnitude units) of the Gaussian distribution :param numpy.ndarray occurrence_rate: Activity rates for magnitude in the range mmin to mmax in steps of bin_width """
[docs] def setUp(self, mfd_conf): """ Input core configuration parameters as specified in the configuration file :param dict mfd_conf: Configuration file containing the following attributes: * 'Model_Weight' - Logic tree weight of model type (float) * 'MFD_spacing' - Width of MFD bin (float) * 'Minimum_Magnitude' - Minimum magnitude of activity rates (float) * 'Maximum_Magnitude' - Characteristic magnituded (float) (if not defined will use scaling relation) * 'Maximum_Magnitude_Uncertainty' - Uncertainty on maximum magnitude (If not defined and the MSR has a sigma term then this will be taken from sigma) * 'Lower_Bound' - Lower bound in terms of number of sigma (float) * 'Upper_Bound' - Upper bound in terms of number of sigma (float) * 'Sigma' - Standard deviation (in magnitude units) of distribution """ self.mfd_model = "Characteristic" self.mfd_weight = mfd_conf["Model_Weight"] self.bin_width = mfd_conf["MFD_spacing"] self.mmin = None self.mmax = None self.mmax_sigma = None self.lower_bound = mfd_conf["Lower_Bound"] self.upper_bound = mfd_conf["Upper_Bound"] self.sigma = mfd_conf["Sigma"] self.occurrence_rate = None
[docs] def get_mmax(self, mfd_conf, msr, rake, area): """ Gets the mmax for the fault - reading directly from the config file or using the msr otherwise :param dict mfd_config: Configuration file (see setUp for paramters) :param msr: Instance of :class: nhlib.scalerel :param float rake: Rake of the fault (in range -180 to 180) :param float area: Area of the fault surface (km^2) """ if mfd_conf["Maximum_Magnitude"]: self.mmax = mfd_conf["Maximum_Magnitude"] else: self.mmax = msr.get_median_mag(area, rake) self.mmax_sigma = mfd_conf.get( "Maximum_Magnitude_Uncertainty", None ) or msr.get_std_dev_mag(None, rake)
[docs] def get_mfd(self, slip, area, shear_modulus=30.0): """ Calculates activity rate on the fault :param float slip: Slip rate in mm/yr :param fault_width: Width of the fault (km) :param float disp_length_ratio: Displacement to length ratio (dimensionless) :param float shear_modulus: Shear modulus of the fault (GPa) :returns: * Minimum Magnitude (float) * Bin width (float) * Occurrence Rates (numpy.ndarray) """ # Working in Nm so convert: shear_modulus - GPa -> Nm # area - km ** 2. -> m ** 2. # slip - mm/yr -> m/yr moment_rate = ( (shear_modulus * 1.0e9) * (area * 1.0e6) * (slip / 1000.0) ) moment_mag = _scale_moment(self.mmax, in_nm=True) characteristic_rate = moment_rate / moment_mag if self.sigma and (fabs(self.sigma) > 1e-5): self.mmin = self.mmax + (self.lower_bound * self.sigma) mag_upper = self.mmax + (self.upper_bound * self.sigma) mag_range = np.arange( self.mmin, mag_upper + self.bin_width, self.bin_width ) self.occurrence_rate = characteristic_rate * ( truncnorm.cdf( mag_range + (self.bin_width / 2.0), self.lower_bound, self.upper_bound, loc=self.mmax, scale=self.sigma, ) - truncnorm.cdf( mag_range - (self.bin_width / 2.0), self.lower_bound, self.upper_bound, loc=self.mmax, scale=self.sigma, ) ) else: # Returns only a single rate self.mmin = self.mmax self.occurrence_rate = np.array([characteristic_rate], dtype=float) return self.mmin, self.bin_width, self.occurrence_rate