Performance issue and/or endless loop in periodic extrapolation

[rruusu]

Description

The while loops used to wrap the input u for periodic extrapolation in several functions like ModelicaTableAdditions_CombiTable1D_getValue() can cause severe performance issues and, in extreme cases, infinite loops.

Code Reference

In ModelicaTableAdditions_CombiTable1D_getValue():

/* Periodic extrapolation */
if (tableID->extrapolation == PERIODIC) {
    const double T = uMax - uMin;

    if (u < uMin) {
        do {
            u += T;
        } while (u < uMin);
    }
    else if (u > uMax) {
        do {
            u -= T;
        } while (u > uMax);
    }
    last = findRowIndex(table, nRow, nCol, tableID->last, u);
    tableID->last = last;
}

Impact & Root Cause

1. Performance Bottleneck: With a large value for u, the loop can take a significant amount of time. This can happen during step-size adjustments even in well-behaved models.
2. Infinite Loop: As u and T are double-precision floats, they are limited by a 53-bit mantissa. If the absolute value of u is larger than T * 2^53, adding or subtracting T has no effect on the value of u, resulting in an endless loop that will freeze the whole simulation.

Proposed Solution

This O(n) logic can be replaced with an O(1) modulo operation using fmod. This eliminates both the performance penalty and the floating-point infinite loop risk.

For example:

/* Periodic extrapolation */
if (tableID->extrapolation == PERIODIC) {
    const double T = uMax - uMin;
    
    if (u < uMin || u > uMax) {
        double rem = fmod(u - uMin, T);
        /* fmod can return negative values, so adjust to strictly positive modulo */
        if (rem < 0.0) {
            rem += T;
        }
        u = uMin + rem;
    }
    
    last = findRowIndex(table, nRow, nCol, tableID->last, u);
    tableID->last = last;
}

Mitigation

This loop can be avoided in Modelica by adding the modulus operation to the assignment of the table input u directly.

From:

table.u = <equation for u>;

To:;

table.u = table.u_min + mod(<equation_for_u>, table.u_max - table.u_min);

If the table is continuously periodic, a noEvent(...) can be added around the modulus to avoid triggering of unnecessary events.

[rruusu]

It seems I just missed the opportunity to report this to https://github.com/tbeu/ModelicaStandardTables, which is the actual code base in which I was bitten by this, using Modelica Standard Library 4.0.0.

The corresponding location there is in function ModelicaStandardTables_CombiTable1D_getValue() (and others).

[tbeu]

In the spirit of TDD can you give me a concrete example table where the issue occurs. Thanks.

[rruusu]

Here is an artificial case that forces the infinite loop issue to occur at 5 seconds of simulation time.

model PeriodicTableTest
  Modelica.Blocks.Tables.CombiTable1Ds table(table=[0, 0; 1, 1; 2, 0], extrapolation=Modelica.Blocks.Types.Extrapolation.Periodic);
  Real u(start=0.0);
  Real y;
equation
  der(u) = 1.0;
  table.u = u;
  y = table.y[1];
  
  // Trigger an infinite loop in table output calculation
  when time >= 5.0 then
    reinit(u, 1e16);
  end when;
annotation(
    experiment(StartTime = 0, StopTime = 10, Tolerance = 1e-06, Interval = 0.02));
end PeriodicTableTest;

If you drop the value to reinit(u, 1e8), it will just massively slow down at the 5 seconds of simulation time.

After changing table.u = u to table.u = table.u_min + mod(u, table.u_max - table.u_min) to mitigate, neither an infinite loop nor a slowdown occurs. With reinit(u, 1e16) the result from 5 s on will be all zeros, though, as the integrator cannot add anything to the u variable at that magnitude either.