Opened 4 years ago

Closed 3 years ago

Last modified 3 years ago

#6400 closed defect (fixed)

Issue with algorithm and the NF in PNLib

Reported by: casella Owned by: perost
Priority: high Milestone: 1.18.0
Component: New Instantiation Version: 1.16.2
Keywords: Cc: lochel

Description (last modified by perost)

Please check ​PNlib.Examples.ConTest.Conflict. The simulation fails at runtime. However, it ​runs fine with the OF.

I flattened the model with the OF and NF, and found interesting differences. For example, the OF flattens this algorithm

algorithm
  when P3.enableOut.delayPassed then
    if P3.enableOut.nOut > 0 then
      P3.enableOut.arcWeightSum := 0.0;
      if P3.enableOut.t + (-P3.enableOut.arcWeightSum) - P3.enableOut.minMarks >= -1e-09 or PNlib.Functions.OddsAndEnds.isEqual(P3.enableOut.arcWeightSum, 0.0, 1e-08) then
      else
        if P3.enableOut.enablingType == PNlib.Types.EnablingType.Priority then
          P3.enableOut.arcWeightSum := 0.0;
        else
          P3.enableOut.arcWeightSum := 0.0;
          P3.enableOut.nremTAout := 0;
          P3.enableOut.nTAout := P3.enableOut.nremTAout;
          P3.enableOut.sumEnablingProbTAout := sum({}[{}[1:P3.enableOut.nremTAout]]);
          P3.enableOut.cumEnablingProb[1] := {}[{}[1]] / P3.enableOut.sumEnablingProbTAout;
          for j in 2:P3.enableOut.nremTAout loop
            P3.enableOut.cumEnablingProb[j] := {}[-1 + j] + {}[{}[j]] / P3.enableOut.sumEnablingProbTAout;
          end for;
          for i in 1:P3.enableOut.nTAout loop
            (P3.enableOut.randNum, P3.enableOut.state128) := Modelica.Math.Random.Generators.Xorshift128plus.random({pre(P3.enableOut.state128[1]), pre(P3.enableOut.state128[2]), pre(P3.enableOut.state128[3]), pre(P3.enableOut.state128[4])});
            P3.enableOut.endWhile := false;
            P3.enableOut.k := 1;
            while P3.enableOut.k <= P3.enableOut.nremTAout and not P3.enableOut.endWhile loop
              if P3.enableOut.randNum <= {}[P3.enableOut.k] then
                P3.enableOut.posTE := {}[P3.enableOut.k];
                P3.enableOut.endWhile := true;
              else
                P3.enableOut.k := 1 + P3.enableOut.k;
              end if;
            end while;
            if P3.enableOut.t + (-{}[P3.enableOut.posTE]) - P3.enableOut.minMarks - P3.enableOut.arcWeightSum >= -1e-09 or PNlib.Functions.OddsAndEnds.isEqual({}[i], 0.0, 1e-08) then
              P3.enableOut.arcWeightSum := P3.enableOut.arcWeightSum + {}[P3.enableOut.posTE];
              P3.enableOut.TEout[P3.enableOut.posTE] := true;
            end if;
            P3.enableOut.nremTAout := -1 + P3.enableOut.nremTAout;
            if P3.enableOut.nremTAout > 0 then
              P3.enableOut.sumEnablingProbTAout := sum({}[{}[1:P3.enableOut.nremTAout]]);
              if P3.enableOut.sumEnablingProbTAout > 0.0 then
                P3.enableOut.cumEnablingProb[1] := {}[{}[1]] / P3.enableOut.sumEnablingProbTAout;
                for j in 2:P3.enableOut.nremTAout loop
                  P3.enableOut.cumEnablingProb[j] := {}[-1 + j] + {}[{}[j]] / P3.enableOut.sumEnablingProbTAout;
                end for;
              else
                P3.enableOut.cumEnablingProb[1:P3.enableOut.nremTAout] := fill(1.0 / /*Real*/(P3.enableOut.nremTAout), P3.enableOut.nremTAout);
              end if;
            end if;
          end for;
        end if;
      end if;
    else
      P3.enableOut.arcWeightSum := 0.0;
      P3.enableOut.nremTAout := 0;
      P3.enableOut.nTAout := 0;
      P3.enableOut.k := 0;
      P3.enableOut.posTE := 0;
      P3.enableOut.randNum := 0.0;
      P3.enableOut.state128 := {pre(P3.enableOut.state128[1]), pre(P3.enableOut.state128[2]), pre(P3.enableOut.state128[3]), pre(P3.enableOut.state128[4])};
      P3.enableOut.sumEnablingProbTAout := 0.0;
      P3.enableOut.endWhile := false;
      P3.enableOut.Index := 0;
    end if;
  end when;

while the NF flattens it as:

algorithm
  P3.enableIn.TEin := {false};
  when P3.enableIn.delayPassed then
    P3.enableIn.disTAin := PNlib.Functions.OddsAndEnds.boolAnd(P3.enableIn.TAein, P3.enableIn.disTransition);
    P3.enableIn.arcWeightSum := PNlib.Functions.OddsAndEnds.conditionalSum(P3.enableIn.arcWeight, P3.enableIn.disTAin);
    if P3.enableIn.t + P3.enableIn.arcWeightSum - P3.enableIn.maxMarks <= 1e-09 or PNlib.Functions.OddsAndEnds.isEqual(P3.enableIn.arcWeightSum, 0.0, 1e-08) then
      P3.enableIn.TEin := P3.enableIn.TAein;
    else
      P3.enableIn.TEin := PNlib.Functions.OddsAndEnds.boolAnd(P3.enableIn.TAein, not P3.enableIn.disTransition);
      if P3.enableIn.enablingType == PNlib.Types.EnablingType.Priority then
        P3.enableIn.arcWeightSum := 0.0;
        for i in 1:1 loop
          P3.enableIn.Index := Modelica.Math.Vectors.find(/*Real*/(i), /*Real[1]*/(P3.enableIn.enablingPrio), 0.0);
          if P3.enableIn.Index > 0 and P3.enableIn.disTAin[P3.enableIn.Index] and (P3.enableIn.t + P3.enableIn.arcWeightSum + P3.enableIn.arcWeight[P3.enableIn.Index] - P3.enableIn.maxMarks <= 1e-09 or PNlib.Functions.OddsAndEnds.isEqual(P3.enableIn.arcWeight[P3.enableIn.Index], 0.0, 1e-08)) then
            P3.enableIn.TEin[P3.enableIn.Index] := true;
            P3.enableIn.arcWeightSum := P3.enableIn.arcWeightSum + P3.enableIn.arcWeight[P3.enableIn.Index];
          end if;
        end for;
      else
        P3.enableIn.arcWeightSum := 0.0;
        P3.enableIn.remTAin := {0};
        P3.enableIn.nremTAin := 0;
        for i in 1:1 loop
          if P3.enableIn.disTAin[i] then
            P3.enableIn.nremTAin := P3.enableIn.nremTAin + 1;
            P3.enableIn.remTAin[P3.enableIn.nremTAin] := i;
          end if;
        end for;
        P3.enableIn.nTAin := P3.enableIn.nremTAin;
        P3.enableIn.sumEnablingProbTAin := sum(P3.enableIn.enablingProb[P3.enableIn.remTAin[1:P3.enableIn.nremTAin]]);
        P3.enableIn.cumEnablingProb := {0.0};
        P3.enableIn.cumEnablingProb[1] := P3.enableIn.enablingProb[P3.enableIn.remTAin[1]] / P3.enableIn.sumEnablingProbTAin;
        for j in 2:P3.enableIn.nremTAin loop
          P3.enableIn.cumEnablingProb[j] := P3.enableIn.cumEnablingProb[j - 1] + P3.enableIn.enablingProb[P3.enableIn.remTAin[j]] / P3.enableIn.sumEnablingProbTAin;
        end for;
        for i in 1:P3.enableIn.nTAin loop
          (P3.enableIn.randNum, P3.enableIn.state128) := Modelica.Math.Random.Generators.Xorshift128plus.random(pre(P3.enableIn.state128)) "uniform distributed random number";
          P3.enableIn.endWhile := false;
          P3.enableIn.k := 1;
          while P3.enableIn.k <= P3.enableIn.nremTAin and not P3.enableIn.endWhile loop
            if P3.enableIn.randNum <= P3.enableIn.cumEnablingProb[P3.enableIn.k] then
              P3.enableIn.posTE := P3.enableIn.remTAin[P3.enableIn.k];
              P3.enableIn.endWhile := true;
            else
              P3.enableIn.k := P3.enableIn.k + 1;
            end if;
          end while;
          if P3.enableIn.t + P3.enableIn.arcWeightSum + P3.enableIn.arcWeight[P3.enableIn.posTE] - P3.enableIn.maxMarks <= 1e-09 or PNlib.Functions.OddsAndEnds.isEqual(P3.enableIn.arcWeight[i], 0.0, 1e-08) then
            P3.enableIn.arcWeightSum := P3.enableIn.arcWeightSum + P3.enableIn.arcWeight[P3.enableIn.posTE];
            P3.enableIn.TEin[P3.enableIn.posTE] := true;
          end if;
          P3.enableIn.nremTAin := P3.enableIn.nremTAin - 1;
          if P3.enableIn.nremTAin > 0 then
            P3.enableIn.remTAin := PNlib.Functions.OddsAndEnds.deleteElementInt(P3.enableIn.remTAin, P3.enableIn.k);
            P3.enableIn.cumEnablingProb := {0.0};
            P3.enableIn.sumEnablingProbTAin := sum(P3.enableIn.enablingProb[P3.enableIn.remTAin[1:P3.enableIn.nremTAin]]);
            if P3.enableIn.sumEnablingProbTAin > 0.0 then
              P3.enableIn.cumEnablingProb[1] := P3.enableIn.enablingProb[P3.enableIn.remTAin[1]] / P3.enableIn.sumEnablingProbTAin;
              for j in 2:P3.enableIn.nremTAin loop
                P3.enableIn.cumEnablingProb[j] := P3.enableIn.cumEnablingProb[j - 1] + P3.enableIn.enablingProb[P3.enableIn.remTAin[j]] / P3.enableIn.sumEnablingProbTAin;
              end for;
            else
              P3.enableIn.cumEnablingProb[1:P3.enableIn.nremTAin] := fill(1.0 / /*Real*/(P3.enableIn.nremTAin), P3.enableIn.nremTAin);
            end if;
          end if;
        end for;
      end if;
    end if;
  end when;
  for i in 1:1 loop
    P3.enableIn.TEin_[i] := P3.enableIn.TEin[i] and P3.enableIn.active[i];
  end for;

It seems to me that what follows when P3.enableIn.delayPassed is quite different in the two cases. Do I miss something?

Change History (8)

comment:1 Changed 4 years ago by perost

  • Description modified (diff)

You're looking at two different when-algorithms, the OF one is for P3.enableOut.delayPassed while the NF one is for P3.enableIn.delayPassed.

The when-algorithms for P3.enableIn.delayPassed are actually quite similar in the OF and the NF, while the one for P3.enableOut.delayPassed is simplified by the NF to just:

  when P3.enableOut.delayPassed then
    P3.enableOut.arcWeightSum := 0.0;
    P3.enableOut.nremTAout := 0;
    P3.enableOut.nTAout := 0;
    P3.enableOut.k := 0;
    P3.enableOut.posTE := 0;
    P3.enableOut.randNum := 0.0;
    P3.enableOut.state128 := pre(P3.enableOut.state128);
    P3.enableOut.sumEnablingProbTAout := 0.0;
    P3.enableOut.endWhile := false;
    P3.enableOut.Index := 0;
  end when;

comment:2 Changed 4 years ago by casella

Ops, sorry...

The question remains open, what is the difference that is making the runtime fail?

comment:3 Changed 4 years ago by lochel

The example works fine when selecting newton as non-linear solver. The example also works fine if Runge–Kutta is selected instead of dassl.

The problem in the example is, that the solution of the non-linear system changes in a discrete way based on the value of P1.t.

solution with P1.t == 0.0

[ 1]      T3.decreasingFactorOut[1]  =                1
[ 2]      T2.decreasingFactorOut[1]  =                1
[ 3]          T2.instantaneousSpeed  =       0.66666667
[ 4]          T3.instantaneousSpeed  =       0.33333333
[ 5]              T1.speedSumOut[1]  =                1
[ 6]               T3.speedSumIn[1]  =                1

solution with P1.t > 0.0

[ 1]      T3.decreasingFactorOut[1]  =                1
[ 2]      T2.decreasingFactorOut[1]  =                1
[ 3]          T2.instantaneousSpeed  =                2
[ 4]          T3.instantaneousSpeed  =                1
[ 5]              T1.speedSumOut[1]  =                3
[ 6]               T3.speedSumIn[1]  =                1

comment:4 Changed 4 years ago by lochel

In other words, it works fine with __OpenModelica_simulationFlags(nls = "newton", s = "rungekutta").

comment:5 Changed 4 years ago by casella

What is strange is that this is not an issue with the OF. Any idea why that is the case?

comment:6 Changed 3 years ago by lochel

My changes got merged into the PNlib and we should now reach 100% coverage.

comment:7 Changed 3 years ago by lochel

  • Resolution set to fixed
  • Status changed from new to closed

comment:8 Changed 3 years ago by casella

Very good, we have another library with 100% support in OpenModelica :)

Note: See TracTickets for help on using tickets.