Numerical data of results shown in *Source time functions of earthquakes based on a stochastic differential equation*:
https://doi.org/10.1038/s41598-022-07873-2
See the above paper for details of case A and case B.
For each cases, the Bessel processes ($X_t^{(1)}$ and $X_t^{(2)}$) and their convolution (STF) are given in ascii format.
$X_t^{(1)}$ has shorter duration than $X_t^{(2)}$.
$n$-th STF (`STF/${n}.txt`) is a convolution of $n$-th $X_t^{(1)}$ (`Xt1/${n}.txt`) and $n$-th $X_t^{(2)}$ (`Xt2/${n}.txt`), where $n=1, ..., 1000$.
STFs are not normalized by their duration or total moment, so the normalization is required to obtain plots in the paper.
The Bessel process as a 1-D array is obtained by the following Julia-1.6.1 function using DifferentialEquations.jl package,
function BesselProcess(Tmin,Tmax)
Xₜ = zeros(Tmax)
(d, σ, u₀, dt) = (0, 1, 1e-3, 1e-6)
while true
f(u,p,t) = ifelse(u≤0, 0, 0.5*(d-1)/u)
g(u,p,t) = ifelse(u≤0, 0, σ)
prob = SDEProblem(f,g,u₀,(0, Tmax*dt))
sol = solve(prob,SRIW1(),dt=dt,adaptive=false)
Xₜ = max.(0, sol.u)
if Xₜ[Tmin-1] > 0 && Xₜ[Tmax] == 0
break
end
end
return Xₜ
end
where `(Tmin, Tmax) = (1000, 2000)` and `(Tmin, Tmax) = (200, 2000)` for case A and B, respectively, and `d=0` is equivalent to $b=1$ for GR law (see eq.(5) in the paper). Note that $d$ must be close to $0$, otherwise, `break` may not work because the duration $T$ does not satisfy $T_\textrm{min} < T < T_\textrm{max}$ probabilistically (see eq.(3) in the paper).
