Initialization of BranchingDynamicPipes with homotopy methods failing

[AnHeuermann]

While changing the way OpenModelica initializes when the homotopy-operator is used inside a model we noticed that Modelica.Fluid.Examples.BranchingDynamicPipes behaves not as expected and is triggering a bunch of asserts.

When using homotopy method for the initialization problem (Advanced.OnlyUseHomotopyMethod = true in Dymola or -homotopyOnFirstTry in OpenModelica) we get an initialization failure.

@casella stated:

I got an initialization failure, with some error messages similar to the ones reported by omc, plus this additional diagnostic message:

Jacobian inverse norm estimate: 1.74013e+012
Condition number estimate: 1.93572e+014
1-norm of the residual = 21864.4
The estimates indicate that the Jacobian is close to singular, suggesting that there is no solution.

I guess this indicates something's wrong in the parametrization of the system, which results in a singular or nearly singular Jacobian at lambda = 0. I'll investigate this further, but I'd say this is a problem with the model, not with the tools.

Changing the parametrization of used homotopy operators inside the model, so that the used homotopy solvers converge to the right solution, could be a solution. At the same time it should not be dependent of the chosen tool to solve the model.

See also OpenModelica ticket 5139.

[HansOlsson]

I believe the homotopy used in the flowModel should revised - or removed.

The reason is that the model has the settable parameters dp_nominal, m_flow_nominal and the computed parameter dp_fric_nominal (which doesn't seem important).

They are (almost) only used for homotopy as far as I can see, and people not using homotopy (as in this case) will not set them to good values (because they don't matter). Thus when homotopy is activated it starts from bad values, and fails - which is counter to the idea of using homotopy to start from known good values.

[HansOlsson]

To clarify I see three possibilities for correcting this general issue for the flowModel (assuming I have understood it):

• Just remove the homotopy-part.
• Replace the homotopy-specific parameters with existing parameters in some way.
• Add a "parameter Boolean activateHomotopyParameters=false;" and make the previously mentioned parameters conditional enabled on this one, and only use the homotopy if this parameter is set. Thus as default there is no homotopy, but users can activate it - with homotopy-specific parameters, so no loss of functionality.

Note 1: There might also be some error in the homotopy-equations; if so correct that. Note 2: This indicates that adding "homotopy" in models is no panacea; a different idea would be to revise ideas of probability-1 homotopy - but that will require more (especially for fluid systems).

[casella]

What I can say is that I have been using homotopy consistently in the context of thermal power plant modelling for more than 10 years, and I found it extremely successful in that area. As soon as I started using it, initializing complex plant models in steady-state became almost trivial (as long as you get the model right, with the proper parameter values, so that a solution actually existst), while previously it was a true pain in the neck.

I wish I something rigorous as probability-1 homotopy was available in this field, but I'm not aware of anything like, and I have no idea whether this is theoretically possible at all.

The rationale of using homotopy for momentum equations is to make them linear, based on nominal parameters, so that they can be solved independently (or almost independently) from the energy balance equations. The coupling of momentum and energy equations in compressible fluids creates horrendous nonlinearity, and homotopy allows to get rid of it to get a first rough-cut estimation of all the iteration variables without the need of specifying them manually one by one, which is terrible.

The same argument goes for using homotopy to initially fix the heat transfer coefficients to some nominal value, rather than having them depend on the flow rate and fluid properties, which causes more horrendously nonlinear behaviour and more pain getting convergence on large (i.e., realistic) models. In this case, homotopy makes the energy equations nearly linear, once the simplified momentum equations have been solved, that is. I guess this was planned for Modelica.Fluid, but I never had time to do that, because I was doing all my modelling using other libraries (e.g. ThermoPower).

That said, I understand that Modelica.Fluid covers a lot more ground than thermal power plants, and may be used in many different contexts, some of which actually don't need this feature at all. So, I agree that we should make this feature optional and not built-in without any way to escape except disabling homotopy in the tool, which is not a good idea because you may need it in other parts of your model.

So, I also agree with @HansOlsson that we should definitely make homotopy optional in these models, as we already did for the opAmp and Pump models in #3291. Then, one is free to turn it on if needed, or to ignore it completely if the solver is good enough or the modeller lucky enough not to get convergence issues.

In order to make this backward compatible, I think we just need to proceed with the third option proposed by @HansOlsson. Given the highly structured nature of the Fluid library, I guess this requires a minimal amount of code changes in some base class. Figuring out how to do it right may not be trivial, because sometimes you have unanticipated consequences in child classes and you need to take care of them.

In principle, the change could not be really 100% backwards compatible, since there could be models around that actually depend on homotopy to initialize successfully, that would be broken when ported to MSL 4.0.0. I guess this is OK as long as we explain this issue in the release notes. After all, adding some homotopy parameters is not big deal, and I don't think there will be so many models affected by this change.

I think we should do this for 4.0.0. You should allow me a couple of weeks to get this done, because one needs to be careful not to break everything up - I need to do some analysis and testing to get thing done properly. I would consider this a bug fix and not a new feature, because ultimately that's what it is: this thing should fix existing models that don't work in some tools.

@beutlich, @dietmarw is this OK with you?

[dietmarw]

I agree with your analysis but it is definitely too late to be included now for 4.0.0. There are just too many moving parts that can get broken by such a big change. In addition rushing exactly this now will definitely end in tears. So rather aim for the 4.0.1 bugfix release.

[dietmarw]

Just to clarify, I'm not against fixing this particular model as reported in this ticket for 4.0.0 if it can be fixed "locally" but would wait with the more general overhaul of the Fluid library until 4.0.1

[casella]

I'm happy moving this to 4.0.1. @dietmarw, do we already have a tag for this milestone? Can you please create one, so I can move it there?

[beutlich]

In principle, the change could not be really 100% backwards compatible, since there could be models around that actually depend on homotopy to initialize successfully, that would be broken when ported to MSL 4.0.0.

That would also disqualify it from being fixed by MSL v4.0.1. Reschudling for MSL v5.0.0.