OpenModelica computes incorrect results if equations are slightly reformulated

[Michael Wetter]

For the code section below, OpenModelica computes different values for the variables wrong and hMed_out:

  XiVap = XiActual[i_x];
  wrong    = (hActual -(1-XiActual[i_x]) * hCon);
  hMed_out = (hActual -(1-XiVap)         * hCon);

This code is in ​https://github.com/iea-annex60/modelica-annex60/blob/issue17_enthalpySensor/Annex60/Fluid/Sensors/LatentEnthalpyFlowRate.mo#L69

The problem can be reproduced by running testOpenModelica.sh as on ​https://github.com/iea-annex60/modelica-annex60/tree/issue17_enthalpySensor

A detailed issue report is accessible at ​https://github.com/iea-annex60/modelica-annex60/issues/17

[Michael Wetter]

Result plot

[Willi Braun]

The behaviour of that model has changed to the one, when the bug was reported. So now we break the simulation, since assertion are triggered now after the initialization:

"Simulation execution failed for model: Annex60.Fluid.Sensors.Examples.MoistAirEnthalpyFlowRate
Error: Inputs differ by more than threShold
  time       = 0.12
  u1         = 45119.7
  u2         = 47313.7
  abs(u1-u2) = 2193.96
  threShold  = 0.01
stdout            | warning | Integrator attempt to handle a problem with a called assert.
Error: Inputs differ by more than threShold
  time       = 0.06
  u1         = 2630.29
  u2         = 2618.61
  abs(u1-u2) = 11.6822
  threShold  = 0.01
stdout            | info    | model terminate | Simulation terminated by an assert at time: 0.06

Also it simulates correctly if we initialize with the result file of dymola.

As far as I see the issue is related to the following model and there the initial algorithm section:

model MassFractionsTwoPort "Ideal two port sensor for mass fraction"
  extends Modelica.Fluid.Sensors.BaseClasses.PartialFlowSensor;
  extends Modelica.Icons.RotationalSensor;
  Modelica.Blocks.Interfaces.RealOutput Xi "Mass fraction in port medium";
  parameter String substanceName = "water" "Name of mass fraction";

protected 
  parameter Integer ind(fixed=false) 
    "Index of species in vector of independent mass fractions";
initial algorithm 
  ind:= -1;
  for i in 1:Medium.nC loop
    if ( Modelica.Utilities.Strings.isEqual(Medium.substanceNames[i], substanceName)) then
      ind := i;
    end if;
  end for;
  assert(ind > 0, "Mass fraction '" + substanceName + "' is not present in medium '"
         + Medium.mediumName + "'.\n"
         + "Check sensor parameter and medium model.");
equation 
  if allowFlowReversal then
     Xi = Modelica.Fluid.Utilities.regStep(port_a.m_flow, port_b.Xi_outflow[ind], port_a.Xi_outflow[ind], m_flow_small);
  else
     Xi = port_b.Xi_outflow[ind];
  end if;
end MassFractionsTwoPort;

[Lennart Ochel]

Replying to wbraun:

As far as I see the issue is related to the following model and there the initial algorithm section:

model MassFractionsTwoPort "Ideal two port sensor for mass fraction"
  extends Modelica.Fluid.Sensors.BaseClasses.PartialFlowSensor;
  extends Modelica.Icons.RotationalSensor;
  Modelica.Blocks.Interfaces.RealOutput Xi "Mass fraction in port medium";
  parameter String substanceName = "water" "Name of mass fraction";

protected 
  parameter Integer ind(fixed=false) 
    "Index of species in vector of independent mass fractions";
initial algorithm 
  ind:= -1;
  for i in 1:Medium.nC loop
    if ( Modelica.Utilities.Strings.isEqual(Medium.substanceNames[i], substanceName)) then
      ind := i;
    end if;
  end for;
  assert(ind > 0, "Mass fraction '" + substanceName + "' is not present in medium '"
         + Medium.mediumName + "'.\n"
         + "Check sensor parameter and medium model.");
equation 
  if allowFlowReversal then
     Xi = Modelica.Fluid.Utilities.regStep(port_a.m_flow, port_b.Xi_outflow[ind], port_a.Xi_outflow[ind], m_flow_small);
  else
     Xi = port_b.Xi_outflow[ind];
  end if;
end MassFractionsTwoPort;

I cannot see any problem with this initial algorithm. It is handled correctly.

[Willi Braun]

Replying to lochel:

Replying to wbraun:

As far as I see the issue is related to the following model and there the initial algorithm section:

model MassFractionsTwoPort "Ideal two port sensor for mass fraction"
  extends Modelica.Fluid.Sensors.BaseClasses.PartialFlowSensor;
  extends Modelica.Icons.RotationalSensor;
  Modelica.Blocks.Interfaces.RealOutput Xi "Mass fraction in port medium";
  parameter String substanceName = "water" "Name of mass fraction";

protected 
  parameter Integer ind(fixed=false) 
    "Index of species in vector of independent mass fractions";
initial algorithm 
  ind:= -1;
  for i in 1:Medium.nC loop
    if ( Modelica.Utilities.Strings.isEqual(Medium.substanceNames[i], substanceName)) then
      ind := i;
    end if;
  end for;
  assert(ind > 0, "Mass fraction '" + substanceName + "' is not present in medium '"
         + Medium.mediumName + "'.\n"
         + "Check sensor parameter and medium model.");
equation 
  if allowFlowReversal then
     Xi = Modelica.Fluid.Utilities.regStep(port_a.m_flow, port_b.Xi_outflow[ind], port_a.Xi_outflow[ind], m_flow_small);
  else
     Xi = port_b.Xi_outflow[ind];
  end if;
end MassFractionsTwoPort;

I cannot see any problem with this initial algorithm. It is handled correctly.

Then please download the library and check that.

[Willi Braun]

Okay, perhaps I was wrong, it's not the initialization.
It's like stated in the beginning, the issue is:

  XiVap = XiActual[i_x];
  wrong    = (hActual -(1-XiActual[i_x]) * hCon);
  hMed_out = (hActual -(1-XiVap)         * hCon);

[Willi Braun]

fixed in r19578.