Summary
MoorDyn_Init() raises FE_OVERFLOW inside the stationary initial-condition solver. The overflow is an intermediate value that is discarded, so the result is correct, but a host program that runs with floating-point traps enabled (as several CFD codes do to catch NaNs early) is killed inside MoorDyn_Init (SIGILL on macOS arm64, SIGFPE on Linux with feenableexcept).
Version: MoorDyn-C v2.7.1 (tag v2.7.1, commit 319ff99), built with CMake, default options (bundled Eigen, PYTHON_WRAPPER=OFF), AppleClang on macOS 26 arm64. The same expression is in master.
Where
moordyn::CFL::cfl2dt(const real cfl, const real v) in source/Util/CFL.hpp returns cfl * length() / v. The CFL base class initialises _l to std::numeric_limits<real>::max() (CFL.hpp line 62) and only lines set it (Line.cpp line 251), so for points, rods and bodies the product cfl * max() overflows to +inf.
StationaryScheme::Step() in source/Time.cpp (around line 206) calls it for every object to limit the next step:
real v = 0.5 * dt * _error;
for (auto obj : lines) new_dt = (std::min)(new_dt, obj->cfl2dt(cfl, v));
for (auto obj : points) new_dt = (std::min)(new_dt, obj->cfl2dt(cfl, v));
...
The +inf from a point is discarded by std::min, which is why the solver still works, but the overflow flag is raised on every stationary step. A backtrace with the trap enabled (RelWithDebInfo build):
stop reason = EXC_BAD_INSTRUCTION
frame #0: libmoordyn.2.dylib`moordyn::CFL::cfl2dt(this=..., cfl=0.0453, v=0.00408) const at CFL.hpp:90 [inlined]
frame #1: moordyn::time::StationaryScheme::Step(double&) (Time.cpp)
frame #2: moordyn::MoorDyn::icStationary()
frame #3: moordyn::MoorDyn::Init(double const*, double const*, bool)
frame #4: MoorDyn_Init
Reproducer
A single line between two fixed points in a light fluid (an overhead conductor in air), external wave kinematics on. No other dependency than libmoordyn.
span.txt:
MoorDyn-C input for a single fixed-fixed conductor span
----------------------- LINE TYPES ------------------------------------------
TypeName Diam Mass/m EA BA/-zeta EI Cd Ca CdAx CaAx
(name) (m) (kg/m) (N) (N-s/-) (N-m^2) (-) (-) (-) (-)
drake 0.0281 1.628 3e+07 -0.5 0 1 1.0 0.0 0.0
---------------------- POINT PROPERTIES --------------------------------
ID Type X Y Z Mass Volume CdA Ca
(#) (-) (m) (m) (m) (kg) (m^3) (m^2) (-)
1 Fixed 0.0 0.0 -100 0 0 0 0
2 Fixed 300 0.0 -100 0 0 0 0
---------------------- LINES ----------------------------------------
ID LineType AttachA AttachB UnstrLen NumSegs LineOutputs
(#) (name) (#) (#) (m) (-) (-)
1 drake 1 2 301.5 20 -
---------------------- OPTIONS -----------------------------------------
0 writeLog Write a log file
0.001 dtM time step to use in the line integration (s)
9.81 g gravity (m/s^2)
1.2 WtrDnsty fluid density (kg/m^3)
1000 WtrDpth depth of the flat bottom (m)
1 WaveKin the fluid kinematics are provided through the API (-)
0 ICgenDynamic stationary initial-condition solver
1 disableOutput
1 disableOutTime
------------------------- need this line --------------------------------------
trap_repro.c:
#include <moordyn/MoorDyn2.h>
#include <fenv.h>
#include <stdio.h>
#include <string.h>
#pragma STDC FENV_ACCESS ON
static void enable_overflow_trap (void)
{
#if defined(__aarch64__) && defined(__APPLE__)
fenv_t env; fegetenv(&env); env.__fpcr |= __fpcr_trap_overflow; fesetenv(&env);
#elif defined(__GLIBC__)
feenableexcept(FE_OVERFLOW);
#endif
}
int main (int argc, char** argv)
{
if (argc > 2 && strcmp(argv[2], "trap") == 0) { enable_overflow_trap(); }
feclearexcept(FE_ALL_EXCEPT);
MoorDyn s = MoorDyn_Create(argv[1]);
MoorDyn_SetVerbosity(s, MOORDYN_ERR_LEVEL);
const int rc = MoorDyn_Init(s, NULL, NULL);
printf("MoorDyn_Init rc = %d; FE_OVERFLOW raised: %s; FE_DIVBYZERO %s; FE_INVALID %s\n", rc,
fetestexcept(FE_OVERFLOW) ? "yes" : "no", fetestexcept(FE_DIVBYZERO) ? "yes" : "no",
fetestexcept(FE_INVALID) ? "yes" : "no");
MoorDyn_Close(s);
return 0;
}$ cc -I$PREFIX/include -L$PREFIX/lib -lmoordyn -Wl,-rpath,$PREFIX/lib trap_repro.c -o trap_repro
$ ./trap_repro span.txt
MoorDyn_Init rc = 0; FE_OVERFLOW raised: yes; FE_DIVBYZERO no; FE_INVALID no
$ ./trap_repro span.txt trap
Illegal instruction: 4
Only the overflow flag is raised; the invalid and divide-by-zero traps pass.
Suggested fix
Keep the sentinel out of the arithmetic, for example in CFL.hpp:
virtual inline real cfl2dt(const real cfl, const real v) const
{
if (_l == (std::numeric_limits<real>::max)()) return _l;
return cfl * length() / v;
}or evaluate length() / v * cfl, which does not overflow for the sentinel when v >= cfl. Either way StationaryScheme::Step gets the same new_dt as today without raising the flag.
Happy to open a pull request with the one-line guard if you prefer.
Summary
MoorDyn_Init()raisesFE_OVERFLOWinside the stationary initial-condition solver. The overflow is an intermediate value that is discarded, so the result is correct, but a host program that runs with floating-point traps enabled (as several CFD codes do to catch NaNs early) is killed insideMoorDyn_Init(SIGILL on macOS arm64, SIGFPE on Linux withfeenableexcept).Version: MoorDyn-C v2.7.1 (tag
v2.7.1, commit 319ff99), built with CMake, default options (bundled Eigen,PYTHON_WRAPPER=OFF), AppleClang on macOS 26 arm64. The same expression is inmaster.Where
moordyn::CFL::cfl2dt(const real cfl, const real v)insource/Util/CFL.hppreturnscfl * length() / v. TheCFLbase class initialises_ltostd::numeric_limits<real>::max()(CFL.hpp line 62) and only lines set it (Line.cppline 251), so for points, rods and bodies the productcfl * max()overflows to+inf.StationaryScheme::Step()insource/Time.cpp(around line 206) calls it for every object to limit the next step:The
+inffrom a point is discarded bystd::min, which is why the solver still works, but the overflow flag is raised on every stationary step. A backtrace with the trap enabled (RelWithDebInfo build):Reproducer
A single line between two fixed points in a light fluid (an overhead conductor in air), external wave kinematics on. No other dependency than
libmoordyn.span.txt:trap_repro.c:Only the overflow flag is raised; the invalid and divide-by-zero traps pass.
Suggested fix
Keep the sentinel out of the arithmetic, for example in
CFL.hpp:or evaluate
length() / v * cfl, which does not overflow for the sentinel whenv >= cfl. Either wayStationaryScheme::Stepgets the samenew_dtas today without raising the flag.Happy to open a pull request with the one-line guard if you prefer.