ModelicaByExample.Components.LotkaVolterra.Examples.ThreeSpecies_Quiescent doesn't simulate

[Adrian Pop]

There is only this model in ModelicaByExample that does not simulate:
​https://test.openmodelica.org/libraries/ModelicaByExample/files/ModelicaByExample.Components.LotkaVolterra.Examples.ThreeSpecies_Quiescent.sim

The model is here:
​https://github.com/xogeny/ModelicaBook
or in the OpenModelica installation

[Lennart Ochel]

Actually, I do not see any issue. I think the model is not that well-constructed. The initial solution (0, 0, 0) for all three populations seem to be valid.

[1] Real foxes.population(start=10, nominal=1) = -4.70198e-038 (pre: 10)
[2] Real rabbits.population(start=10, nominal=1) = 0 (pre: 10)
[3] Real wolves.population(start=10, nominal=1) = 2.35099e-038 (pre: 10)

Due to numerical reasons, an assertion is triggered.
BTW: It would be good if the assertion would be a bit more helpful (see #3074). That means in this particular case, it would be good to show the assertion itself as well (foxes.population >= 0.0).

[Lennart Ochel]

I just noticed that the Simulation output is not in the description:

assert            | debug   | Population must be greater than zero
assert            | info    | simulation terminated by an assertion at initialization

[Adrian Pop]

Michael, can you fix this model?

[Lennart Ochel]

There is a non-linear system during initialization that need to be solved. Actually, I did not check if there is another suitable solution for this model. But if that is the case, it would probably help to modify the following start attributes:

Notification: List of all iteration variables (DAE kind: initialization)
Iteration variables of torn nonlinear equation system:
  wolves.population:VARIABLE(start = 10.0 protected = true )  "Population in this region".ModelicaByExample.Components.LotkaVolterra.Examples.ThreeSpecies_Quiescent, .ModelicaByExample.Components.LotkaVolterra.Components.RegionalPopulation$wolves, .Real type: Real 
  foxes.population:VARIABLE(start = 10.0 protected = true )  "Population in this region".ModelicaByExample.Components.LotkaVolterra.Examples.ThreeSpecies_Quiescent, .ModelicaByExample.Components.LotkaVolterra.Components.RegionalPopulation$foxes, .Real type: Real 
  rabbits.population:VARIABLE(start = 10.0 protected = true )  "Population in this region".ModelicaByExample.Components.LotkaVolterra.Examples.ThreeSpecies_Quiescent, .ModelicaByExample.Components.LotkaVolterra.Components.RegionalPopulation$rabbits, .Real type: Real

You get the above output using +d=iterationVars as omc debug flag.

[Lennart Ochel]

A suitable solution seem to be:

[1] Real foxes.population = 1.66667
[2] Real rabbits.population = 26.6667
[3] Real wolves.population = 3.33333

This one can be achieved e.g. if either start value of foxes.population is set to 2 or start value of rabbits.population is set to 100.

[Adrian Pop]

Lennart, maybe you can fork ​https://github.com/xogeny/ModelicaBook on github and then do a pull request :)

[Lennart Ochel]

Replying to adrpo:

Lennart, maybe you can fork ​https://github.com/xogeny/ModelicaBook on github and then do a pull request :)

​https://github.com/xogeny/ModelicaBook/pull/228

[Martin Sjölund]

Why not add the assertions to the initial NLS so the NLS solver knows the solution is out of bounds? :)

[Lennart Ochel]

Replying to sjoelund.se:

Why not add the assertions to the initial NLS so the NLS solver knows the solution is out of bounds? :)

In that case the NLS solver would still try to find the trivial solution (0, 0, 0) as it is valid.

[Martin Sjölund]

Not if it knows to look for constraints or you feed it back high residuals to show that the region is invalid :)

[Lennart Ochel]

The zero-solution is valid. I just think it is not the one of interest. Hence, one has to either specify additional constrains or set proper start values.

[Martin Sjölund]

It's not valid if it triggers an assertion...

[Lennart Ochel]

The assertion is not part of the alg. loop. Therefore, nothing is wrong with the calculated solution. If the NLS solver should be aware of the additional constraint on population, I propose following changes:

ModelicaByExample/Components/LotkaVolterra/Components/RegionalPopulation.mo

@@ -13 +13 @@
 model RegionalPopulation "Population of animals in a specific region"
     annotation (Placement(transformation(extent={{-10,90},{10,110}}),
         iconTransformation(extent={{-10,90},{10,110}})));
 protected
-  Real population(start=10) = species.population "Population in this region";
+  Real population(start=10, min=0) = species.population "Population in this region";
 initial equation
   if init==InitializationOptions.FixedPopulation then

ModelicaByExample/Components/LotkaVolterra/Components/RegionalPopulation.mo

@@ -23 +23 @@
 initial equation
   end if;
 equation
   der(population) = species.rate;
-  assert(population>=0, "Population must be greater than zero");
   annotation (Diagram(graphics={Rectangle(
           extent={{-100,100},{100,-100}},
           lineColor={0,127,0},

[Martin Sjölund]

Replying to lochel:

The assertion is not part of the alg. loop.

You are wrong here. It is part of the alg. loop. It is just not part of the alg. loop generated by OpenModelica.

[Martin Sjölund]

Note that in Modelica, a min-attribute on a Real should expand to a corresponding assertion.
If you would argue that a min-attribute is part of the algebraic loop, so should assertions that depend only on variables in that same algebraic loop.
Related would be setting an initial equation assert(phi > -pi and phi < pi, "..."); because you want all tools to initialize to the same angle.

[Lennart Ochel]

Replying to sjoelund.se:

Replying to lochel:

The assertion is not part of the alg. loop.

You are wrong here. It is part of the alg. loop. It is just not part of the alg. loop generated by OpenModelica.

Hmm ... I might agree. In that case we should open a ticket about it :-)

[Lennart Ochel]

Merged ​#228

[Dietmar Winkler]

It doesn't make sense to keep closed ticket in the "Future" milestone that were simply forgotten to assign to the correct milestone in the past.

[Martin Sjölund]

Removing the pre1.9.4 milestone in favor of 1.9.4.