Skip to content

Commit e5fed1e

Browse files
committed
FIX: Avoid divide by zero in KinSrcRamp slip time function (zero slip or impulse duration).
1 parent 4dd5347 commit e5fed1e

File tree

2 files changed

+9
-4
lines changed

2 files changed

+9
-4
lines changed

CHANGES.md

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -26,6 +26,7 @@ The version numbers are in the form `MAJOR.MINOR.PATCH`, where major releases in
2626
* Add compatibility for reading PETSc HDF5 format
2727
* **Fixed**
2828
* Account for processes without cells in initial condition patch when verifying configuration.
29+
* Avoid divide by zero for `KinSrcRamp` when final slip is zero.
2930
* Fix pythia import in `pylith_eqinfo`.
3031

3132
## Version 4.2.0 (2025-01-15)

libsrc/pylith/faults/KinSrcRamp.cc

Lines changed: 8 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -193,7 +193,9 @@ pylith::faults::KinSrcRamp::slipRateFn(const PylithInt dim,
193193
const double finalSlipMag = dim == 2 ?
194194
sqrt(pow(finalSlip[0],2) + pow(finalSlip[1],2)) :
195195
sqrt(pow(finalSlip[0],2) + pow(finalSlip[1],2) + pow(finalSlip[2],2));
196-
const double slipCoef = 2.0 * _KinSrcRamp::maxAcc(finalSlipMag, riseTime, impulseDuration) / impulseDuration;
196+
const double slipCoef = impulseDuration > 0.0 ?
197+
2.0 * _KinSrcRamp::maxAcc(finalSlipMag, riseTime, impulseDuration) / impulseDuration :
198+
0.0;
197199

198200
double slipRateTimeFn = 0.0;
199201
if (t-t0 < 0.0) {
@@ -217,7 +219,7 @@ pylith::faults::KinSrcRamp::slipRateFn(const PylithInt dim,
217219
} // if/else
218220

219221
for (PylithInt i = 0; i < dim; ++i) {
220-
slipRate[i] = slipCoef * slipRateTimeFn * finalSlip[i] / finalSlipMag;
222+
slipRate[i] = finalSlipMag > 0.0 ? slipCoef * slipRateTimeFn * finalSlip[i] / finalSlipMag : 0.0
221223
} // for
222224
} // slipRateFn
223225

@@ -266,7 +268,9 @@ pylith::faults::KinSrcRamp::slipAccFn(const PylithInt dim,
266268
const double finalSlipMag = dim == 2 ?
267269
sqrt(pow(finalSlip[0],2) + pow(finalSlip[1],2)) :
268270
sqrt(pow(finalSlip[0],2) + pow(finalSlip[1],2) + pow(finalSlip[2],2));
269-
const double slipCoef = 2.0 * _KinSrcRamp::maxAcc(finalSlipMag, riseTime, impulseDuration) / impulseDuration;
271+
const double slipCoef = impulseDuration > 0.0 ?
272+
2.0 * _KinSrcRamp::maxAcc(finalSlipMag, riseTime, impulseDuration) / impulseDuration :
273+
0.0;
270274

271275
double slipAccTimeFn = 0.0;
272276
if (t-t0 < 0.0) {
@@ -286,7 +290,7 @@ pylith::faults::KinSrcRamp::slipAccFn(const PylithInt dim,
286290
} // if/else
287291

288292
for (PylithInt i = 0; i < dim; ++i) {
289-
slipAcc[i] = slipCoef * slipAccTimeFn * finalSlip[i] / finalSlipMag;
293+
slipAcc[i] = finalSlipMag > 0.0 ? slipCoef * slipAccTimeFn * finalSlip[i] / finalSlipMag : 0.0;
290294
} // for
291295
} // slipAccFn
292296

0 commit comments

Comments
 (0)