Opened 10 years ago

Closed 10 years ago

#2510 closed defect (fixed)

OpenModelica computes incorrect results if equations are slightly reformulated

Reported by: mwetter Owned by: wbraun
Priority: high Milestone: 1.9.1
Component: Backend Version: trunk
Keywords: Cc: lochel, wbraun

Description

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

Attachments (1)

results.png (25.4 KB) - added by mwetter 10 years ago.
Result plot

Download all attachments as: .zip

Change History (7)

Changed 10 years ago by mwetter

Result plot

comment:1 Changed 10 years ago by sjoelund.se

  • Cc lochel wbraun added

comment:2 follow-up: Changed 10 years ago by wbraun

  • Owner changed from probably noone to lochel
  • Status changed from new to assigned

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;

comment:3 in reply to: ↑ 2 ; follow-up: Changed 10 years ago by 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.

Last edited 10 years ago by lochel (previous) (diff)

comment:4 in reply to: ↑ 3 Changed 10 years ago by wbraun

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.

comment:5 Changed 10 years ago by wbraun

  • Owner changed from lochel to wbraun

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);

comment:6 Changed 10 years ago by wbraun

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

fixed in r19578.

Note: See TracTickets for help on using tickets.