Stochastic Loading Module
vlachos_et_al.cc
1 #define _USE_MATH_DEFINES
2 #include <algorithm>
3 #include <cmath>
4 #include <ctime>
5 #include <memory>
6 #include <numeric>
7 #include <stdexcept>
8 #include <string>
9 #include <vector>
10 // Boost random generator
11 #include <boost/random/mersenne_twister.hpp>
12 #include <boost/random/uniform_real_distribution.hpp>
13 #include <boost/random/variate_generator.hpp>
14 // Eigen dense matrices
15 #include <Eigen/Dense>
16 
17 #include "factory.h"
18 #include "function_dispatcher.h"
19 #include "json_object.h"
20 #include "lognormal_dist.h"
21 #include "normal_dist.h"
22 #include "normal_multivar.h"
23 #include "numeric_utils.h"
24 #include "vlachos_et_al.h"
25 
26 stochastic::VlachosEtAl::VlachosEtAl(double moment_magnitude,
27  double rupture_distance, double vs30,
28  double orientation,
29  unsigned int num_spectra,
30  unsigned int num_sims)
31  : StochasticModel(),
32  moment_magnitude_{moment_magnitude / 6.0},
33  rupture_dist_{(rupture_distance + 5.0) / 30.0},
34  vs30_{vs30 / 450.0},
35  orientation_{orientation},
36  time_step_{0.01},
37  freq_step_{0.2},
38  cutoff_freq_{220.0},
39  num_spectra_{num_spectra},
40  num_sims_{num_sims},
41  seed_value_{std::numeric_limits<int>::infinity()},
42  model_parameters_{18} {
43  model_name_ = "VlachosEtAl";
44  // Factors for site condition based on Vs30
45  double site_soft = 0.0, site_medium = 0.0, site_hard = 0.0;
46  if (vs30 <= 300.0) {
47  site_soft = 1.0;
48  } else if (vs30 <= 450.0) {
49  site_medium = 1.0;
50  } else {
51  site_hard = 1.0;
52  }
53 
54  // Estimated conditional mean values
55  Eigen::VectorXd conditional_means(7);
56  // clang-format off
57  conditional_means <<
58  1.0, moment_magnitude_, std::log(rupture_dist_),
59  moment_magnitude_ * std::log(rupture_dist_), site_soft * std::log(vs30_),
60  site_medium * std::log(vs30_), site_hard * std::log(vs30_);
61  // clang-format on
62 
63  // Restricted Maximum Likelihood method regression coefficients and variance
64  // components of the normal model parameters (Table 3 on page 13)
65  Eigen::MatrixXd beta(18, 7);
66  // clang-format off
67  beta <<
68  -1.1417, 1.0917, 1.9125, -0.9696, 0.0971, 0.3476, -0.6740,
69  1.8052,-1.8381, -3.5874, 3.7895, 0.3236, 0.5497, 0.2876,
70  1.8969,-1.8819, -2.0818, 1.9000, -0.3520, -0.6959, -0.0025,
71  1.6627,-1.6922, -1.2509, 1.1880, -0.5170, -1.0157, -0.1041,
72  3.8703,-3.4745, -0.0816, 0.0166, 0.4904, 0.8697, 0.3179,
73  1.1043,-1.1852, -1.0068, 0.9388, -0.5603, -0.8855, -0.3174,
74  1.1935,-1.2922, -0.7028, 0.6975, -0.6629, -1.1075, -0.4542,
75  1.7895,-1.5014, -0.0300, -0.1306, 0.4526, 0.7132, 0.1522,
76  -3.6404, 3.3189, -0.5316, 0.3874, -0.3757, -0.8334, 0.1006,
77  -2.2742, 2.1454, 0.6315, -0.6620, 0.1093, -0.1028, -0.0479,
78  0.6930, -0.6202, 1.8037, -1.6064, 0.0727, -0.1498, -0.0722,
79  1.3003, -1.2004, -1.2210, 1.0623, -0.0252, 0.1885, 0.0069,
80  0.4604, -0.4087, -0.5057, 0.4486, 0.1073, -0.0219, -0.1352,
81  2.2304, -2.0398, -0.1364, 0.1910, 0.2425, 0.1801, 0.3233,
82  2.3806, -2.2011, -0.3256, 0.2226, -0.0221, 0.0970, 0.0762,
83  0.2057, -0.1714, 0.3385, -0.2229, 0.0802, 0.2649, 0.0396,
84  -7.6011, 6.8507, -2.3609, 0.9201, -0.7508, -0.7903, -0.6204,
85  -6.3472, 5.8241, 3.2994, -2.8774, -0.1411, -0.5298, -0.0203;
86  // clang-format on
87 
88  // Variance of model parameters (Table 3 on page 13)
89  Eigen::VectorXd variance(18);
90  // clang-format off
91  variance <<
92  0.90, 0.80, 0.78, 0.74, 0.66, 0.73, 0.72, 0.70, 0.69,
93  0.78, 0.90, 0.90, 0.90, 0.90, 0.80, 0.90, 0.35, 0.80;
94  // clang-format on
95 
96  // Estimated correlation matrix (Table A1 on page 24)
97  Eigen::MatrixXd correlation_matrix(18, 18);
98  // clang-format off
99  correlation_matrix <<
100  1.0000, 0.0382, -0.0912, -0.0701, -0.0214, -0.0849, -0.0545, -0.0185, 0.0270, -0.0122, 0.0059, -0.0344, -0.0342, 0.0409, -0.0137, -0.0168, -0.0990, -0.6701,
101  0.0382, 1.0000, -0.1159, -0.1856, 0.0681, -0.2018, -0.2765, -0.0304, -0.1719, -0.1157, -0.0347, -0.0277, -0.0189, 0.0357, 0.0657, -0.0070, 0.3690, -0.0510,
102  -0.0912, -0.1159, 1.0000, 0.9467, 0.4123, 0.4815, 0.4240, 0.2120, 0.1070, -0.1898, 0.0506, -0.0661, -0.0380, 0.0260, 0.0506, -0.0317, -0.0278, 0.0245,
103  -0.0701, -0.1856, 0.9467, 1.0000, 0.4075, 0.4891, 0.4940, 0.2285, 0.2009, -0.1709, 0.0365, -0.0579, -0.0999, 0.0467, 0.0410, 0.0027, -0.0966, 0.0631,
104  -0.0214, 0.0681, 0.4123, 0.4075, 1.0000, 0.1772, 0.1337, 0.7315, -0.0066, -0.2787, 0.0703, -0.0541, -0.0453, 0.1597, 0.0792, 0.0220, 0.0606, -0.0844,
105  -0.0849, -0.2018, 0.4815, 0.4891, 0.1772, 1.0000, 0.9448, 0.3749, 0.1682, -0.0831, 0.0124, -0.1236, -0.0346, -0.0054, 0.0877, -0.0197, -0.0867, 0.0281,
106  -0.0545, -0.2765, 0.4240, 0.4940, 0.1337, 0.9448, 1.0000, 0.3530, 0.2305, -0.0546, -0.0223, -0.0782, -0.0872, 0.0074, 0.0999, 0.0066, -0.1358, 0.0626,
107  -0.0185, -0.0304, 0.2120, 0.2285, 0.7315, 0.3749, 0.3530, 1.0000, 0.1939, -0.0617, -0.0017, -0.0942, -0.0332, 0.0813, 0.0810, -0.0032, -0.0870, -0.0599,
108  0.0270, -0.1719, 0.1070, 0.2009, -0.0066, 0.1682, 0.2305, 0.1939, 1.0000, -0.1851, -0.2073, -0.0756, -0.1637, -0.0865, 0.0699, -0.0485, -0.2153, 0.0320,
109  -0.0122, -0.1157, -0.1898, -0.1709, -0.2787, -0.0831, -0.0546, -0.0617, -0.1851, 1.0000, 0.2139, 0.0769, 0.1391, 0.0769, -0.1838, 0.0377, -0.1615, 0.1000,
110  0.0059, -0.0347, 0.0506, 0.0365, 0.0703, 0.0124, -0.0223, -0.0017, -0.2073, 0.2139, 1.0000, -0.1102, -0.0530, 0.0791, 0.0012, 0.0090, -0.0236, 0.0037,
111  -0.0344, -0.0277, -0.0661, -0.0579, -0.0541, -0.1236, -0.0782, -0.0942, -0.0756, 0.0769, -0.1102, 1.0000, -0.2562, -0.0406, 0.3154, 0.0065, -0.0093, -0.0354,
112  -0.0342, -0.0189, -0.0380, -0.0999, -0.0453, -0.0346, -0.0872, -0.0332, -0.1637, 0.1391, -0.0530, -0.2562, 1.0000, -0.1836, -0.1624, -0.5646, 0.0216, 0.0243,
113  0.0409, 0.0357, 0.0260, 0.0467, 0.1597, -0.0054, 0.0074, 0.0813, -0.0865, 0.0769, 0.0791, -0.0406, -0.1836, 1.0000, 0.1624, 0.1989, 0.0549, -0.0411,
114  -0.0137, 0.0657, 0.0506, 0.0410, 0.0792, 0.0877, 0.0999, 0.0810, 0.0699, -0.1838, 0.0012, 0.3154, -0.1624, 0.1624, 1.0000, 0.1552, 0.0844, -0.0637,
115  -0.0168, -0.0070, -0.0317, 0.0027, 0.0220, -0.0197, 0.0066, -0.0032, -0.0485, 0.0377, 0.0090, 0.0065, -0.5646, 0.1989, 0.1552, 1.0000, 0.0058, 0.0503,
116  -0.0990, 0.3690, -0.0278, -0.0966, 0.0606, -0.0867, -0.1358, -0.0870, -0.2153, -0.1615, -0.0236, -0.0093, 0.0216, 0.0549, 0.0844, 0.0058, 1.0000, -0.0930,
117  -0.6701, -0.0510, 0.0245, 0.0631, -0.0844, 0.0281, 0.0626, -0.0599, 0.0320, 0.1000, 0.0037, -0.0354, 0.0243, -0.0411, -0.0637, 0.0503, -0.0930, 1.0000;
118  // clang-format on
119 
120  // Mean of transformed normal model parameters (described by Eq. 25 on page 12)
121  means_ = beta * conditional_means;
122 
123  // Convert the standard deviation and correlation to covariance
124  covariance_ = numeric_utils::corr_to_cov(correlation_matrix,
125  (variance.array().sqrt()).matrix());
126 
127  // Generate realizations of model parameters
128  sample_generator_ =
129  Factory<numeric_utils::RandomGenerator>::instance()->create(
130  "MultivariateNormal");
131  sample_generator_->generate(parameter_realizations_, means_, covariance_,
132  num_spectra_);
133  parameter_realizations_.transposeInPlace();
134 
135  // Create distributions for model parameters
136  model_parameters_[0] =
137  Factory<stochastic::Distribution, double, double>::instance()->create(
138  "LognormalDist", std::move(-1.735), std::move(0.523));
139  model_parameters_[1] =
140  Factory<stochastic::Distribution, double, double>::instance()->create(
141  "LognormalDist", std::move(1.009), std::move(0.422));
142  model_parameters_[2] =
143  Factory<stochastic::Distribution, double, double>::instance()->create(
144  "NormalDist", std::move(0.249), std::move(1.759));
145  model_parameters_[3] =
146  Factory<stochastic::Distribution, double, double>::instance()->create(
147  "NormalDist", std::move(0.768), std::move(1.958));
148  model_parameters_[4] =
149  Factory<stochastic::Distribution, double, double>::instance()->create(
150  "LognormalDist", std::move(2.568), std::move(0.557));
151  model_parameters_[5] =
152  Factory<stochastic::Distribution, double, double>::instance()->create(
153  "NormalDist", std::move(0.034), std::move(1.471));
154  model_parameters_[6] =
155  Factory<stochastic::Distribution, double, double>::instance()->create(
156  "NormalDist", std::move(0.441), std::move(1.733));
157  model_parameters_[7] =
158  Factory<stochastic::Distribution, double, double>::instance()->create(
159  "LognormalDist", std::move(3.356), std::move(0.473));
160  model_parameters_[8] =
161  Factory<stochastic::Distribution, double, double>::instance()->create(
162  "BetaDist", std::move(2.516), std::move(9.174));
163  model_parameters_[9] =
164  Factory<stochastic::Distribution, double, double>::instance()->create(
165  "BetaDist", std::move(3.582), std::move(15.209));
166  model_parameters_[10] =
167  Factory<stochastic::Distribution, double, double>::instance()->create(
168  "LognormalDist", std::move(0.746), std::move(0.404));
169  model_parameters_[11] =
170  Factory<stochastic::Distribution, double, double, double>::instance()
171  ->create("StudentstDist", std::move(0.205), std::move(0.232),
172  std::move(7.250));
173  model_parameters_[12] =
174  Factory<stochastic::Distribution, double, double>::instance()->create(
175  "InverseGaussianDist", std::move(0.499), std::move(0.213));
176  model_parameters_[13] =
177  Factory<stochastic::Distribution, double, double>::instance()->create(
178  "LognormalDist", std::move(0.702), std::move(0.435));
179  model_parameters_[14] =
180  Factory<stochastic::Distribution, double, double, double>::instance()
181  ->create("StudentstDist", std::move(0.792), std::move(0.157),
182  std::move(4.223));
183  model_parameters_[15] =
184  Factory<stochastic::Distribution, double, double>::instance()->create(
185  "InverseGaussianDist", std::move(0.350), std::move(0.170));
186  model_parameters_[16] =
187  Factory<stochastic::Distribution, double, double>::instance()->create(
188  "LognormalDist", std::move(9.470), std::move(1.317));
189  model_parameters_[17] =
190  Factory<stochastic::Distribution, double, double>::instance()->create(
191  "LognormalDist", std::move(3.658), std::move(0.375));
192 
193  // Standard normal distribution with mean at 0.0 and standard deviation of 1.0
194  auto std_normal_dist =
195  Factory<stochastic::Distribution, double, double>::instance()->create(
196  "NormalDist", std::move(0.0), std::move(1.0));
197 
198  physical_parameters_.resize(parameter_realizations_.rows(),
199  parameter_realizations_.cols());
200 
201  // Transform sample normal model parameters to physical space
202  for (unsigned int i = 0; i < parameter_realizations_.rows(); ++i) {
203  for (unsigned int j = 0; j < model_parameters_.size(); ++j) {
204  physical_parameters_(i, j) =
205  (model_parameters_[j]->inv_cumulative_dist_func(
206  std_normal_dist->cumulative_dist_func(
207  std::vector<double>{parameter_realizations_(i, j)})))[0];
208  }
209  }
210 }
211 
212 
213 stochastic::VlachosEtAl::VlachosEtAl(double moment_magnitude,
214  double rupture_distance, double vs30,
215  double orientation,
216  unsigned int num_spectra,
217  unsigned int num_sims,
218  int seed_value)
219  : StochasticModel(),
220  moment_magnitude_{moment_magnitude / 6.0},
221  rupture_dist_{(rupture_distance + 5.0) / 30.0},
222  vs30_{vs30 / 450.0},
223  orientation_{orientation},
224  time_step_{0.01},
225  freq_step_{0.2},
226  cutoff_freq_{220.0},
227  num_spectra_{num_spectra},
228  num_sims_{num_sims},
229  seed_value_{seed_value},
230  model_parameters_{18} {
231  model_name_ = "VlachosEtAl";
232  // Factors for site condition based on Vs30
233  double site_soft = 0.0, site_medium = 0.0, site_hard = 0.0;
234  if (vs30 <= 300.0) {
235  site_soft = 1.0;
236  } else if (vs30 <= 450.0) {
237  site_medium = 1.0;
238  } else {
239  site_hard = 1.0;
240  }
241 
242  // Estimated conditional mean values
243  Eigen::VectorXd conditional_means(7);
244  // clang-format off
245  conditional_means <<
246  1.0, moment_magnitude_, std::log(rupture_dist_),
247  moment_magnitude_ * std::log(rupture_dist_), site_soft * std::log(vs30_),
248  site_medium * std::log(vs30_), site_hard * std::log(vs30_);
249  // clang-format on
250 
251  // Restricted Maximum Likelihood method regression coefficients and variance
252  // components of the normal model parameters (Table 3 on page 13)
253  Eigen::MatrixXd beta(18, 7);
254  // clang-format off
255  beta <<
256  -1.1417, 1.0917, 1.9125, -0.9696, 0.0971, 0.3476, -0.6740,
257  1.8052,-1.8381, -3.5874, 3.7895, 0.3236, 0.5497, 0.2876,
258  1.8969,-1.8819, -2.0818, 1.9000, -0.3520, -0.6959, -0.0025,
259  1.6627,-1.6922, -1.2509, 1.1880, -0.5170, -1.0157, -0.1041,
260  3.8703,-3.4745, -0.0816, 0.0166, 0.4904, 0.8697, 0.3179,
261  1.1043,-1.1852, -1.0068, 0.9388, -0.5603, -0.8855, -0.3174,
262  1.1935,-1.2922, -0.7028, 0.6975, -0.6629, -1.1075, -0.4542,
263  1.7895,-1.5014, -0.0300, -0.1306, 0.4526, 0.7132, 0.1522,
264  -3.6404, 3.3189, -0.5316, 0.3874, -0.3757, -0.8334, 0.1006,
265  -2.2742, 2.1454, 0.6315, -0.6620, 0.1093, -0.1028, -0.0479,
266  0.6930, -0.6202, 1.8037, -1.6064, 0.0727, -0.1498, -0.0722,
267  1.3003, -1.2004, -1.2210, 1.0623, -0.0252, 0.1885, 0.0069,
268  0.4604, -0.4087, -0.5057, 0.4486, 0.1073, -0.0219, -0.1352,
269  2.2304, -2.0398, -0.1364, 0.1910, 0.2425, 0.1801, 0.3233,
270  2.3806, -2.2011, -0.3256, 0.2226, -0.0221, 0.0970, 0.0762,
271  0.2057, -0.1714, 0.3385, -0.2229, 0.0802, 0.2649, 0.0396,
272  -7.6011, 6.8507, -2.3609, 0.9201, -0.7508, -0.7903, -0.6204,
273  -6.3472, 5.8241, 3.2994, -2.8774, -0.1411, -0.5298, -0.0203;
274  // clang-format on
275 
276  // Variance of model parameters (Table 3 on page 13)
277  Eigen::VectorXd variance(18);
278  // clang-format off
279  variance <<
280  0.90, 0.80, 0.78, 0.74, 0.66, 0.73, 0.72, 0.70, 0.69,
281  0.78, 0.90, 0.90, 0.90, 0.90, 0.80, 0.90, 0.35, 0.80;
282  // clang-format on
283 
284  // Estimated correlation matrix (Table A1 on page 24)
285  Eigen::MatrixXd correlation_matrix(18, 18);
286  // clang-format off
287  correlation_matrix <<
288  1.0000, 0.0382, -0.0912, -0.0701, -0.0214, -0.0849, -0.0545, -0.0185, 0.0270, -0.0122, 0.0059, -0.0344, -0.0342, 0.0409, -0.0137, -0.0168, -0.0990, -0.6701,
289  0.0382, 1.0000, -0.1159, -0.1856, 0.0681, -0.2018, -0.2765, -0.0304, -0.1719, -0.1157, -0.0347, -0.0277, -0.0189, 0.0357, 0.0657, -0.0070, 0.3690, -0.0510,
290  -0.0912, -0.1159, 1.0000, 0.9467, 0.4123, 0.4815, 0.4240, 0.2120, 0.1070, -0.1898, 0.0506, -0.0661, -0.0380, 0.0260, 0.0506, -0.0317, -0.0278, 0.0245,
291  -0.0701, -0.1856, 0.9467, 1.0000, 0.4075, 0.4891, 0.4940, 0.2285, 0.2009, -0.1709, 0.0365, -0.0579, -0.0999, 0.0467, 0.0410, 0.0027, -0.0966, 0.0631,
292  -0.0214, 0.0681, 0.4123, 0.4075, 1.0000, 0.1772, 0.1337, 0.7315, -0.0066, -0.2787, 0.0703, -0.0541, -0.0453, 0.1597, 0.0792, 0.0220, 0.0606, -0.0844,
293  -0.0849, -0.2018, 0.4815, 0.4891, 0.1772, 1.0000, 0.9448, 0.3749, 0.1682, -0.0831, 0.0124, -0.1236, -0.0346, -0.0054, 0.0877, -0.0197, -0.0867, 0.0281,
294  -0.0545, -0.2765, 0.4240, 0.4940, 0.1337, 0.9448, 1.0000, 0.3530, 0.2305, -0.0546, -0.0223, -0.0782, -0.0872, 0.0074, 0.0999, 0.0066, -0.1358, 0.0626,
295  -0.0185, -0.0304, 0.2120, 0.2285, 0.7315, 0.3749, 0.3530, 1.0000, 0.1939, -0.0617, -0.0017, -0.0942, -0.0332, 0.0813, 0.0810, -0.0032, -0.0870, -0.0599,
296  0.0270, -0.1719, 0.1070, 0.2009, -0.0066, 0.1682, 0.2305, 0.1939, 1.0000, -0.1851, -0.2073, -0.0756, -0.1637, -0.0865, 0.0699, -0.0485, -0.2153, 0.0320,
297  -0.0122, -0.1157, -0.1898, -0.1709, -0.2787, -0.0831, -0.0546, -0.0617, -0.1851, 1.0000, 0.2139, 0.0769, 0.1391, 0.0769, -0.1838, 0.0377, -0.1615, 0.1000,
298  0.0059, -0.0347, 0.0506, 0.0365, 0.0703, 0.0124, -0.0223, -0.0017, -0.2073, 0.2139, 1.0000, -0.1102, -0.0530, 0.0791, 0.0012, 0.0090, -0.0236, 0.0037,
299  -0.0344, -0.0277, -0.0661, -0.0579, -0.0541, -0.1236, -0.0782, -0.0942, -0.0756, 0.0769, -0.1102, 1.0000, -0.2562, -0.0406, 0.3154, 0.0065, -0.0093, -0.0354,
300  -0.0342, -0.0189, -0.0380, -0.0999, -0.0453, -0.0346, -0.0872, -0.0332, -0.1637, 0.1391, -0.0530, -0.2562, 1.0000, -0.1836, -0.1624, -0.5646, 0.0216, 0.0243,
301  0.0409, 0.0357, 0.0260, 0.0467, 0.1597, -0.0054, 0.0074, 0.0813, -0.0865, 0.0769, 0.0791, -0.0406, -0.1836, 1.0000, 0.1624, 0.1989, 0.0549, -0.0411,
302  -0.0137, 0.0657, 0.0506, 0.0410, 0.0792, 0.0877, 0.0999, 0.0810, 0.0699, -0.1838, 0.0012, 0.3154, -0.1624, 0.1624, 1.0000, 0.1552, 0.0844, -0.0637,
303  -0.0168, -0.0070, -0.0317, 0.0027, 0.0220, -0.0197, 0.0066, -0.0032, -0.0485, 0.0377, 0.0090, 0.0065, -0.5646, 0.1989, 0.1552, 1.0000, 0.0058, 0.0503,
304  -0.0990, 0.3690, -0.0278, -0.0966, 0.0606, -0.0867, -0.1358, -0.0870, -0.2153, -0.1615, -0.0236, -0.0093, 0.0216, 0.0549, 0.0844, 0.0058, 1.0000, -0.0930,
305  -0.6701, -0.0510, 0.0245, 0.0631, -0.0844, 0.0281, 0.0626, -0.0599, 0.0320, 0.1000, 0.0037, -0.0354, 0.0243, -0.0411, -0.0637, 0.0503, -0.0930, 1.0000;
306  // clang-format on
307 
308  // Mean of transformed normal model parameters (described by Eq. 25 on page 12)
309  means_ = beta * conditional_means;
310 
311  // Convert the standard deviation and correlation to covariance
312  covariance_ = numeric_utils::corr_to_cov(correlation_matrix,
313  (variance.array().sqrt()).matrix());
314 
315  // Generate realizations of model parameters
316  sample_generator_ =
317  Factory<numeric_utils::RandomGenerator, int>::instance()->create(
318  "MultivariateNormal", std::move(seed_value_));
319  sample_generator_->generate(parameter_realizations_, means_, covariance_,
320  num_spectra_);
321  parameter_realizations_.transposeInPlace();
322 
323  // Create distributions for model parameters
324  model_parameters_[0] =
325  Factory<stochastic::Distribution, double, double>::instance()->create(
326  "LognormalDist", std::move(-1.735), std::move(0.523));
327  model_parameters_[1] =
328  Factory<stochastic::Distribution, double, double>::instance()->create(
329  "LognormalDist", std::move(1.009), std::move(0.422));
330  model_parameters_[2] =
331  Factory<stochastic::Distribution, double, double>::instance()->create(
332  "NormalDist", std::move(0.249), std::move(1.759));
333  model_parameters_[3] =
334  Factory<stochastic::Distribution, double, double>::instance()->create(
335  "NormalDist", std::move(0.768), std::move(1.958));
336  model_parameters_[4] =
337  Factory<stochastic::Distribution, double, double>::instance()->create(
338  "LognormalDist", std::move(2.568), std::move(0.557));
339  model_parameters_[5] =
340  Factory<stochastic::Distribution, double, double>::instance()->create(
341  "NormalDist", std::move(0.034), std::move(1.471));
342  model_parameters_[6] =
343  Factory<stochastic::Distribution, double, double>::instance()->create(
344  "NormalDist", std::move(0.441), std::move(1.733));
345  model_parameters_[7] =
346  Factory<stochastic::Distribution, double, double>::instance()->create(
347  "LognormalDist", std::move(3.356), std::move(0.473));
348  model_parameters_[8] =
349  Factory<stochastic::Distribution, double, double>::instance()->create(
350  "BetaDist", std::move(2.516), std::move(9.714));
351  model_parameters_[9] =
352  Factory<stochastic::Distribution, double, double>::instance()->create(
353  "BetaDist", std::move(3.582), std::move(15.209));
354  model_parameters_[10] =
355  Factory<stochastic::Distribution, double, double>::instance()->create(
356  "LognormalDist", std::move(0.746), std::move(0.404));
357  model_parameters_[11] =
358  Factory<stochastic::Distribution, double, double, double>::instance()
359  ->create("StudentstDist", std::move(0.205), std::move(0.232),
360  std::move(7.250));
361  model_parameters_[12] =
362  Factory<stochastic::Distribution, double, double>::instance()->create(
363  "InverseGaussianDist", std::move(0.499), std::move(0.213));
364  model_parameters_[13] =
365  Factory<stochastic::Distribution, double, double>::instance()->create(
366  "LognormalDist", std::move(0.702), std::move(0.435));
367  model_parameters_[14] =
368  Factory<stochastic::Distribution, double, double, double>::instance()
369  ->create("StudentstDist", std::move(0.792), std::move(0.157),
370  std::move(4.223));
371  model_parameters_[15] =
372  Factory<stochastic::Distribution, double, double>::instance()->create(
373  "InverseGaussianDist", std::move(0.350), std::move(0.170));
374  model_parameters_[16] =
375  Factory<stochastic::Distribution, double, double>::instance()->create(
376  "LognormalDist", std::move(9.470), std::move(1.317));
377  model_parameters_[17] =
378  Factory<stochastic::Distribution, double, double>::instance()->create(
379  "LognormalDist", std::move(3.658), std::move(0.375));
380 
381  // Standard normal distribution with mean at 0.0 and standard deviation of 1.0
382  auto std_normal_dist =
383  Factory<stochastic::Distribution, double, double>::instance()->create(
384  "NormalDist", std::move(0.0), std::move(1.0));
385 
386  physical_parameters_.resize(parameter_realizations_.rows(),
387  parameter_realizations_.cols());
388 
389  // Transform sample normal model parameters to physical space
390  for (unsigned int i = 0; i < parameter_realizations_.rows(); ++i) {
391  for (unsigned int j = 0; j < model_parameters_.size(); ++j) {
392  physical_parameters_(i, j) =
393  (model_parameters_[j]->inv_cumulative_dist_func(
394  std_normal_dist->cumulative_dist_func(
395  std::vector<double>{parameter_realizations_(i, j)})))[0];
396  }
397  }
398 }
399 
400 utilities::JsonObject stochastic::VlachosEtAl::generate(
401  const std::string& event_name, bool g_units) {
402 
403  // Pool of acceleration time histories based on number of spectra and
404  // simulations requested
405  std::vector<std::vector<std::vector<double>>> acceleration_pool(
406  num_spectra_,
407  std::vector<std::vector<double>>(num_sims_, std::vector<double>()));
408 
409  // Generate family of time histories for each spectrum. Family size is
410  // specified by requested number of simulations per spectra.
411  try {
412  for (unsigned int i = 0; i < num_spectra_; ++i) {
413  time_history_family(acceleration_pool[i], physical_parameters_.row(i));
414  }
415  } catch (const std::exception& e) {
416  std::cerr << e.what();
417  throw;
418  }
419 
420  // Create JsonObject for events
421  auto events = utilities::JsonObject();
422  std::vector<utilities::JsonObject> events_array(num_spectra_ * num_sims_);
423 
424  // Add pattern information for JSON
425  auto pattern_x = utilities::JsonObject();
426  auto pattern_y = utilities::JsonObject();
427  pattern_x.add_value("type", "UniformAcceleration");
428  pattern_x.add_value("timeSeries", "accel_x");
429  pattern_x.add_value("dof", 1);
430  pattern_y.add_value("type", "UniformAcceleration");
431  pattern_y.add_value("timeSeries", "accel_y");
432  pattern_y.add_value("dof", 2);
433 
434  // Create JSON for specific event
435  auto event_data = utilities::JsonObject();
436  // Loop over spectra
437  for (unsigned int i = 0; i < num_spectra_; ++i) {
438  // Loop over different simulations for current spectra
439  for (unsigned int j = 0; j < num_sims_; ++j) {
440  event_data.add_value("name", event_name + "_Spectra" + std::to_string(i) +
441  "_Sim" + std::to_string(j));
442  event_data.add_value("type", "Seismic");
443  event_data.add_value("dT", time_step_);
444  event_data.add_value("numSteps", acceleration_pool[i][j].size());
445  event_data.add_value(
446  "pattern", std::vector<utilities::JsonObject>{pattern_x, pattern_y});
447 
448  // Rotate accelerations, if necessary
449  std::vector<double> x_accels(acceleration_pool[i][j].size());
450  std::vector<double> y_accels(acceleration_pool[i][j].size());
451  rotate_acceleration(acceleration_pool[i][j], x_accels, y_accels, g_units);
452 
453  // Add time histories for x and y directions to event
454  auto time_history_x = utilities::JsonObject();
455  auto time_history_y = utilities::JsonObject();
456  time_history_x.add_value("name", "accel_x");
457  time_history_x.add_value("type", "Value");
458  time_history_x.add_value("dT", time_step_);
459  time_history_x.add_value("data", x_accels);
460  time_history_y.add_value("name", "accel_y");
461  time_history_y.add_value("type", "Value");
462  time_history_y.add_value("dT", time_step_);
463  time_history_y.add_value("data", y_accels);
464  event_data.add_value("timeSeries", std::vector<utilities::JsonObject>{
465  time_history_x, time_history_y});
466  events_array[i * num_sims_ + j] = event_data;
467  event_data.clear();
468  }
469  }
470 
471  events.add_value("Events", events_array);
472 
473  return events;
474 }
475 
476 bool stochastic::VlachosEtAl::generate(const std::string& event_name,
477  const std::string& output_location,
478  bool g_units) {
479  bool status = true;
480 
481  // Generate pool of acceleration time histories
482  try{
483  auto json_output = generate(event_name, g_units);
484  json_output.write_to_file(output_location);
485  } catch (const std::exception& e) {
486  std::cerr << e.what();
487  status = false;
488  throw;
489  }
490 
491  return status;
492 }
493 
494 bool stochastic::VlachosEtAl::time_history_family(
495  std::vector<std::vector<double>>& time_histories,
496  const Eigen::VectorXd& parameters) const {
497  bool status = true;
498  auto identified_parameters = identify_parameters(parameters);
499 
500  unsigned int num_times =
501  static_cast<unsigned int>(std::ceil(identified_parameters[17] / time_step_)) + 1;
502  unsigned int num_freqs =
503  static_cast<unsigned int>(std::ceil(cutoff_freq_ / freq_step_)) + 1;
504 
505  std::vector<double> times(num_times);
506  std::vector<double> frequencies(num_freqs);
507  double total_time = (num_times - 1) * time_step_;
508 
509  for (unsigned int i = 0; i < times.size(); ++i) {
510  times[i] = i * time_step_ / total_time;
511  }
512  times[0] = 1E-6;
513 
514  for (unsigned int i = 0; i < frequencies.size(); ++i) {
515  frequencies[i] = i * freq_step_;
516  }
517 
518  // Calculate non-dimensional energy accumulation
519  auto energy = energy_accumulation(
520  std::vector<double>{identified_parameters[0], identified_parameters[1]},
521  times);
522 
523  // Calculate mode 1 and 2 dominant frequencies
524  auto mode_1_freqs = modal_frequencies(
525  std::vector<double>{identified_parameters[2], identified_parameters[3],
526  identified_parameters[4]},
527  energy);
528  auto mode_2_freqs = modal_frequencies(
529  std::vector<double>{identified_parameters[5], identified_parameters[6],
530  identified_parameters[7]},
531  energy);
532 
533  // Calculate 2nd modal participation factor
534  auto mode_2_participation = modal_participation_factor(
535  std::vector<double>{identified_parameters[10], identified_parameters[11],
536  identified_parameters[12], identified_parameters[13],
537  identified_parameters[14], identified_parameters[15]},
538  energy);
539 
540  // Calculate amplitude modulating function
541  auto amplitude_modulation = amplitude_modulating_function(
542  identified_parameters[17], identified_parameters[16],
543  std::vector<double>{identified_parameters[0], identified_parameters[1]},
544  times);
545 
546  // Parameters for high-pass Butterworth filter
547  int filter_order = 4;
548  double norm_cutoff_freq = 0.20;
549 
550  // Calculate energy content of the Butterworth filter transfer function
551  std::vector<double> highpass_butter_energy(frequencies.size());
552  double freq_ratio_sq;
553  for (unsigned int i = 0; i < frequencies.size(); ++i) {
554  freq_ratio_sq = std::pow(frequencies[i] / (2.0 * M_PI * norm_cutoff_freq),
555  2 * filter_order);
556  highpass_butter_energy[i] = freq_ratio_sq / (1.0 + freq_ratio_sq);
557  }
558 
559  // Calculate the evolutionary power spectrum with unit variance at
560  // each time step
561  Eigen::MatrixXd power_spectrum(times.size(), frequencies.size());
562 
563  for (unsigned int i = 0; i < times.size(); ++i) {
564  power_spectrum.row(i) =
565  kt_2(std::vector<double>{mode_1_freqs[i], identified_parameters[8], 1.0,
566  mode_2_freqs[i], identified_parameters[9],
567  mode_2_participation[i]},
568  frequencies, highpass_butter_energy);
569 
570  double freq_domain_integral =
571  2.0 * numeric_utils::trapazoid_rule(power_spectrum.row(i), freq_step_);
572 
573  power_spectrum.row(i) =
574  power_spectrum.row(i) * amplitude_modulation[i] / freq_domain_integral;
575 
576  }
577 
578  // Get coefficients for highpass Butterworth filter
579  int num_samples =
580  static_cast<int>(std::round(1.5 * static_cast<double>(filter_order) /
581  (2.0 * norm_cutoff_freq)) /
582  time_step_ +
583  1);
584 
585  auto hp_butter =
586  Dispatcher<std::vector<std::vector<double>>, int, double>::instance()
587  ->dispatch("HighPassButter", filter_order,
588  norm_cutoff_freq / (1.0 / time_step_ / 2.0));
589 
590  // Calculate filter impulse response for calculated number of samples
591  auto impulse_response =
592  Dispatcher<std::vector<double>, std::vector<double>, std::vector<double>,
593  int, int>::instance()
594  ->dispatch("ImpulseResponse", hp_butter[0], hp_butter[1],
595  filter_order, num_samples);
596 
597  try {
598  // Generate family of time histories
599  for (unsigned int i = 0; i < num_sims_; ++i) {
600  simulate_time_history(time_histories[i], power_spectrum);
601  post_process(time_histories[i], impulse_response);
602  }
603  } catch (const std::exception& e) {
604  std::cerr << e.what();
605  status = false;
606  throw;
607  }
608 
609  return status;
610 }
611 
612 void stochastic::VlachosEtAl::simulate_time_history(
613  std::vector<double>& time_history,
614  const Eigen::MatrixXd& power_spectrum) const {
615  unsigned int num_times = power_spectrum.rows(),
616  num_freqs = power_spectrum.cols();
617 
618  time_history.resize(num_times, 0.0);
619 
620  std::vector<double> times(num_times);
621  std::vector<double> frequencies(num_freqs);
622 
623  for (unsigned int i = 0; i < times.size(); ++i) {
624  times[i] = i * time_step_;
625  }
626 
627  for (unsigned int i = 0; i < frequencies.size(); ++i) {
628  frequencies[i] = i * freq_step_;
629  }
630 
631  static unsigned int history_seed = static_cast<unsigned int>(std::time(nullptr));
632  history_seed = history_seed + 10;
633 
634  auto generator =
635  seed_value_ != std::numeric_limits<int>::infinity()
636  ? boost::random::mt19937(static_cast<unsigned int>(seed_value_ + 10))
637  : boost::random::mt19937(history_seed);
638 
639  boost::random::uniform_real_distribution<> distribution(0.0, 2.0 * M_PI);
640  boost::random::variate_generator<boost::random::mt19937&,
641  boost::random::uniform_real_distribution<>>
642  angle_gen(generator, distribution);
643 
644  std::vector<double> phase_angle(num_freqs, 0.0);
645 
646  for (auto & angle : phase_angle) {
647  angle = angle_gen();
648  }
649 
650  // Loop over all frequencies and times to calculate time history
651  for (unsigned int i = 0; i < num_times; ++i) {
652  for (unsigned int j = 0; j < num_freqs; ++j) {
653  time_history[i] =
654  time_history[i] +
655  std::sqrt(power_spectrum(i, j)) *
656  std::cos(frequencies[j] * times[i] + phase_angle[j]);
657  }
658  time_history[i] = 2.0 * std::sqrt(freq_step_) * time_history[i];
659  }
660 }
661 
662 bool stochastic::VlachosEtAl::post_process(
663  std::vector<double>& time_history,
664  const std::vector<double>& filter_imp_resp) const {
665 
666  bool status = true;
667  double time_hann_2 = 1.0;
668 
669  Eigen::VectorXd window = Eigen::VectorXd::Ones(time_history.size());
670  unsigned int window1_size =
671  static_cast<unsigned int>(time_hann_2 / time_step_ + 1);
672  unsigned int window2_size = static_cast<unsigned int>((window1_size - 1) / 2);
673 
674  // Check if input time history length is sufficient
675  if (time_history.size() < window1_size) {
676  throw std::runtime_error(
677  "\nERROR: in stochastic::VlachosEtAl::post_process: Input time history "
678  "too short for Hanning Window size\n");
679  }
680 
681  // Get Hanning window of length window1_size
682  auto hann_window =
683  Dispatcher<Eigen::VectorXd, unsigned int>::instance()->dispatch(
684  "HannWindow", window1_size);
685 
686  window.head(window2_size) = hann_window.head(window2_size);
687  window.tail(window2_size) = hann_window.head(window2_size).reverse();
688 
689  // Calculate mean of time history
690  double mean = std::accumulate(time_history.begin(), time_history.end(), 0.0) /
691  static_cast<double>(time_history.size());
692 
693  // Apply window
694  for (unsigned int i = 0; i < time_history.size(); ++i) {
695  time_history[i] = window[i] * (time_history[i] - mean);
696  }
697 
698  // Apply 4th order Butterworth filter
699  std::vector<double> filtered_history(filter_imp_resp.size() +
700  time_history.size() - 1);
701  try {
702  numeric_utils::convolve_1d(filter_imp_resp, time_history, filtered_history);
703  } catch (const std::exception& e) {
704  std::cerr << e.what();
705  status = false;
706  throw;
707  }
708 
709  // Copy filtered results to time_history
710  time_history = filtered_history;
711 
712  return status;
713 }
714 
715 Eigen::VectorXd stochastic::VlachosEtAl::identify_parameters(
716  const Eigen::VectorXd& initial_params) const {
717  // Initialize non-dimensional cumulative energy
718  std::vector<double> energy(static_cast<unsigned int>(1.0 / 0.05) + 1, 0.0);
719 
720  for (unsigned int i = 1; i < energy.size(); ++i) {
721  energy[i] = energy[i - 1] + 0.05;
722  }
723 
724  // Initialze mode 1 parameters and frequencies
725  std::vector<double> mode_1_params = {initial_params(2),
726  initial_params(3),
727  initial_params(4)};
728  auto mode_1_freqs = modal_frequencies(mode_1_params, energy);
729 
730  // Initialize mode 2 parameters and frequencies
731  std::vector<double> mode_2_params = {initial_params(5),
732  initial_params(6),
733  initial_params(7)};
734  auto mode_2_freqs = modal_frequencies(mode_2_params, energy);
735 
736  double mode_1_mean = initial_params(11), mode_2_mean = initial_params(14);
737 
738  // Standard normal distribution with mean at 0.0 and standard deviation of 1.0
739  auto std_normal_dist =
740  Factory<stochastic::Distribution, double, double>::instance()->create(
741  "NormalDist", std::move(0.0), std::move(1.0));
742 
743  Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> realizations(
744  initial_params.size(), 1);
745  Eigen::VectorXd transformed_realizations = initial_params;
746 
747  // Check if any mode 1 dominant frequencies across all non-dimensional energy
748  // values are greater than the corresponding values for mode 2
749  bool freq_comparison = false;
750  for (unsigned int i = 0; i < mode_1_freqs.size(); ++i) {
751  if (mode_1_freqs[i] > mode_2_freqs[i]) {
752  freq_comparison = true;
753  break;
754  }
755  }
756 
757  // Iterate until suitable parameter values have been identified
758  while (freq_comparison || (mode_1_mean > mode_2_mean)) {
759 
760  // Generate realizations of parameters
761  sample_generator_->generate(realizations, means_, covariance_, 1);
762 
763  // Transform parameter realizations to physical space
764  for (unsigned int i = 0; i < initial_params.size(); ++i) {
765  transformed_realizations(i) =
766  (model_parameters_[i]->inv_cumulative_dist_func(
767  std_normal_dist->cumulative_dist_func(
768  std::vector<double>{realizations(i, 0)})))[0];
769  }
770 
771  // Calculate dominant modal frequencies
772  mode_1_freqs =
773  modal_frequencies(std::vector<double>{transformed_realizations(2),
774  transformed_realizations(3),
775  transformed_realizations(4)},
776  energy);
777 
778  mode_2_params = {transformed_realizations(5), transformed_realizations(6),
779  transformed_realizations(7)};
780  mode_2_freqs =
781  modal_frequencies(std::vector<double>{transformed_realizations(5),
782  transformed_realizations(6),
783  transformed_realizations(7)},
784  energy);
785 
786  mode_1_mean = transformed_realizations(11);
787  mode_2_mean = transformed_realizations(14);
788 
789  // Check if any mode 1 dominant frequencies across all non-dimensional energy values
790  // are greater than the corresponding values for mode 2
791  freq_comparison = false;
792  for (unsigned int i = 0; i < mode_1_freqs.size(); ++i) {
793  if (mode_1_freqs[i] > mode_2_freqs[i]) {
794  freq_comparison = true;
795  break;
796  }
797  }
798  }
799 
800  return transformed_realizations;
801 }
802 
803 std::vector<double> stochastic::VlachosEtAl::modal_frequencies(
804  const std::vector<double>& parameters,
805  const std::vector<double>& energy) const {
806  std::vector<double> frequencies(energy.size());
807 
808  for (unsigned int i = 0; i < energy.size(); ++i) {
809  frequencies[i] = parameters[2] * std::pow(0.5 + energy[i], parameters[0]) *
810  std::pow(1.5 - energy[i], parameters[1]);
811  }
812 
813  return frequencies;
814 }
815 
816 std::vector<double> stochastic::VlachosEtAl::energy_accumulation(
817  const std::vector<double>& parameters,
818  const std::vector<double>& times) const {
819  std::vector<double> accumulated_energy(times.size());
820 
821  for (unsigned int i = 0; i < times.size(); ++i) {
822  accumulated_energy[i] =
823  std::exp(-std::pow(times[i] / parameters[0], -parameters[1])) /
824  std::exp(-std::pow(1.0 / parameters[0], -parameters[1]));
825  }
826 
827  return accumulated_energy;
828 }
829 
830 std::vector<double> stochastic::VlachosEtAl::modal_participation_factor(
831  const std::vector<double>& parameters,
832  const std::vector<double>& energy) const {
833  std::vector<double> participation_factor(energy.size());
834 
835  for (unsigned int i = 0; i < energy.size(); ++i) {
836  participation_factor[i] = std::pow(
837  10.0,
838  parameters[0] * std::exp(-std::pow(
839  (energy[i] - parameters[1]) / parameters[2], 2)) +
840  parameters[3] *
841  std::exp(
842  -std::pow((energy[i] - parameters[4]) / parameters[5], 2)) -
843  2.0);
844  }
845 
846  return participation_factor;
847 }
848 
849 std::vector<double> stochastic::VlachosEtAl::amplitude_modulating_function(
850  double duration, double total_energy, const std::vector<double>& parameters,
851  const std::vector<double>& times) const {
852  std::vector<double> func_vals(times.size());
853 
854  double mult_term = total_energy * parameters[1] / (parameters[0] * duration);
855  double exponent_1 = std::pow(1.0 / parameters[0], -parameters[1]);
856 
857  for (unsigned int i = 0; i < times.size(); ++i) {
858  func_vals[i] = mult_term *
859  std::exp(exponent_1 - std::pow(times[i] / parameters[0], -parameters[1])) *
860  std::pow(times[i] / parameters[0], -1.0 - parameters[1]);
861  }
862 
863  return func_vals;
864 }
865 
866 Eigen::VectorXd stochastic::VlachosEtAl::kt_2(
867  const std::vector<double>& parameters,
868  const std::vector<double>& frequencies,
869  const std::vector<double>& highpass_butter) const {
870  Eigen::VectorXd power_spectrum(frequencies.size());
871  double mode1 = 0, mode2 = 0;
872 
873  for (unsigned int i = 0; i < frequencies.size(); ++i) {
874  mode1 = parameters[2] *
875  (1.0 + 4.0 * parameters[1] * parameters[1] *
876  std::pow(frequencies[i] / parameters[0], 2)) /
877  (std::pow(1 - std::pow(frequencies[i] / parameters[0], 2), 2) +
878  4.0 * parameters[1] * parameters[1] *
879  std::pow(frequencies[i] / parameters[0], 2));
880  mode2 = parameters[5] *
881  (1.0 + 4.0 * parameters[4] * parameters[4] *
882  std::pow(frequencies[i] / parameters[3], 2)) /
883  (std::pow(1 - std::pow(frequencies[i] / parameters[3], 2), 2) +
884  4.0 * parameters[4] * parameters[4] *
885  std::pow(frequencies[i] / parameters[3], 2));
886 
887  power_spectrum[i] = highpass_butter[i] * (mode1 + mode2);
888  }
889 
890  return power_spectrum;
891 }
892 
893 void stochastic::VlachosEtAl::rotate_acceleration(
894  const std::vector<double>& acceleration, std::vector<double>& x_accels,
895  std::vector<double>& y_accels, bool g_units) const {
896 
897  x_accels.resize(acceleration.size());
898  y_accels.resize(acceleration.size());
899 
900  double conversion_factor = g_units ? 100.0 * 9.81 : 100.0;
901 
902  // No orientation specified to acceleration oriented along x-axis
903  if (std::abs(orientation_) < 1E-6) {
904  for (unsigned int i = 0; i < acceleration.size(); ++i) {
905  // Division by conversion_factor to convert either to m/s^2 or g
906  x_accels[i] = acceleration[i] / conversion_factor;
907  }
908  y_accels.assign(acceleration.size(), 0.0);
909  // Rotate accelerations to match orientation
910  } else {
911  for (unsigned int i = 0; i < acceleration.size(); ++i) {
912  // Division by conversion_factor to convert either to m/s^2 or g
913  x_accels[i] =
914  acceleration[i] * std::cos(orientation_ * M_PI / 180.0) / conversion_factor;
915  y_accels[i] =
916  acceleration[i] * std::sin(orientation_ * M_PI / 180.0) / conversion_factor;
917  }
918  }
919 }
numeric_utils::convolve_1d
bool convolve_1d(const std::vector< double > &input_x, const std::vector< double > &input_y, std::vector< double > &response)
Definition: numeric_utils.cc:21
stochastic::VlachosEtAl::post_process
bool post_process(std::vector< double > &time_history, const std::vector< double > &filter_imp_resp) const
Definition: vlachos_et_al.cc:662
stochastic::VlachosEtAl::energy_accumulation
std::vector< double > energy_accumulation(const std::vector< double > &parameters, const std::vector< double > &times) const
Definition: vlachos_et_al.cc:816
stochastic::VlachosEtAl::kt_2
Eigen::VectorXd kt_2(const std::vector< double > &parameters, const std::vector< double > &frequencies, const std::vector< double > &highpass_butter) const
Definition: vlachos_et_al.cc:866
Factory::create
std::shared_ptr< Tbaseclass > create(const std::string &key, Targs &&...args)
Definition: factory.h:70
Dispatcher::dispatch
Treturntype dispatch(const std::string &key, Targs...args)
Definition: function_dispatcher.h:47
Factory::instance
static Factory * instance()
Definition: factory.h:49
Dispatcher::instance
static Dispatcher * instance()
Definition: function_dispatcher.h:21
stochastic::VlachosEtAl::modal_frequencies
std::vector< double > modal_frequencies(const std::vector< double > &parameters, const std::vector< double > &energy) const
Definition: vlachos_et_al.cc:803
stochastic::VlachosEtAl::amplitude_modulating_function
std::vector< double > amplitude_modulating_function(double duration, double total_energy, const std::vector< double > &parameters, const std::vector< double > &times) const
Definition: vlachos_et_al.cc:849
stochastic::VlachosEtAl::time_history_family
bool time_history_family(std::vector< std::vector< double >> &time_histories, const Eigen::VectorXd &parameters) const
Definition: vlachos_et_al.cc:494
stochastic::VlachosEtAl::rotate_acceleration
void rotate_acceleration(const std::vector< double > &acceleration, std::vector< double > &x_accels, std::vector< double > &y_accels, bool g_units) const
Definition: vlachos_et_al.cc:893
stochastic::VlachosEtAl::identify_parameters
Eigen::VectorXd identify_parameters(const Eigen::VectorXd &initial_params) const
Definition: vlachos_et_al.cc:715
stochastic::VlachosEtAl::generate
utilities::JsonObject generate(const std::string &event_name, bool g_units=false) override
Definition: vlachos_et_al.cc:400
numeric_utils::trapazoid_rule
double trapazoid_rule(const std::vector< double > &input_vector, double spacing)
Definition: numeric_utils.cc:64
numeric_utils::corr_to_cov
Eigen::MatrixXd corr_to_cov(const Eigen::MatrixXd &corr, const Eigen::VectorXd &std_dev)
Definition: numeric_utils.cc:8
stochastic::VlachosEtAl::simulate_time_history
void simulate_time_history(std::vector< double > &time_history, const Eigen::MatrixXd &power_spectrum) const
Definition: vlachos_et_al.cc:612
stochastic::VlachosEtAl::modal_participation_factor
std::vector< double > modal_participation_factor(const std::vector< double > &parameters, const std::vector< double > &energy) const
Definition: vlachos_et_al.cc:830
Generated by   doxygen 1.8.11