Iinitialization doesn't converge even when near solution

[crupp@…]

Here is a simple model of a mass-spring system that doesn't converge in initialization even though the start values are very near the solution. Eventually the correct solution is identified, but LOG_NLS_V indicates non-convergence and gives no indication how it got the final solution. Setting the start position of the masses to r_0[1]=0.4 and 0.7 even terminates the solve because it evaluates when the spring lengths zero, resulting in a divide-by-zero error. Tested on r24700.

loadString("
model chain
  parameter Real k = 1.0e3;
  parameter Real L = 1.0;
  parameter Real Lstretch = 1.001;
  inner Modelica.Mechanics.MultiBody.World world(g = 0.0);
  Modelica.Mechanics.MultiBody.Forces.Spring spring1(c = k, s_unstretched = L / 3, fixedRotationAtFrame_b = true);
  Modelica.Mechanics.MultiBody.Forces.Spring spring2(c = k, s_unstretched = L / 3, fixedRotationAtFrame_b = true);
  Modelica.Mechanics.MultiBody.Forces.Spring spring3(c = k, s_unstretched = L / 3, fixedRotationAtFrame_b = false);
  Modelica.Mechanics.MultiBody.Parts.PointMass mass1(m = 0.5, r_0(start = {L / 3, 0.0, 0.0}, each fixed = false), v_0(start = {0.0, 0.0, 0.0}, each fixed = true), a_0(start = {0.0, 0.0, 0.0}, each fixed = true));
  Modelica.Mechanics.MultiBody.Parts.PointMass mass2(m = 0.5, r_0(start = {2 * L / 3, 0.0, 0.0}, each fixed = false), v_0(start = {0.0, 0.0, 0.0}, each fixed = true), a_0(start = {0.0, 0.0, 0.0}, each fixed = true));
  Modelica.Mechanics.MultiBody.Parts.Fixed fixed(r = {Lstretch, 0, 0}, animation = false);
equation
  connect(world.frame_b, spring1.frame_a);
  connect(spring1.frame_b, spring2.frame_a);
  connect(spring1.frame_b, mass1.frame_a);
  connect(spring2.frame_b, spring3.frame_a);
  connect(spring2.frame_b, mass2.frame_a);
  connect(spring3.frame_b, fixed.frame_b);
end chain;
");

simulate(chain, simflags="-lv=LOG_NLS_V");

I can see several reasons why the optimizer is not converging based on the information from LOG_NLS_V:
1) In the first Newton iteration, the position of mass2 is moved from 2/3 to 1.001. This is a huge step relative to both its original value and its nominal value.
2) nls variable spring1.spring.s_rel starts at 0.0 because it is doesn't have a defined start value. Ideally it would have a start value of 1/3. It doesn't show up with +d=initialization, so a user would have no reason to believe it needs to be set in the first place.
3) Simiarly, nls variable spring2.spring.f doesn't have a manual start value and hence starts at 0.0, which is very close to the correct value. Also, even though it is reported with +d=initialization, it is unreasonable to have a user define start values for cut forces since these are almost always unknown. Being required to set start values for such obscure (and random) variables makes it very difficult to use OpenModelica and will make it impossible for industry-size models.

A few potential solutions that may help this and other models:
1) Setting the initial damping value for the damped Newton scheme (command line option)
2) Setting a relative max step size for the iteration variables based on their nominal value (command line option)
3) Using a different optimizer to solve the system (e.g., IPOPT), Newton is often overly sensitive to start values
4) Pre-calculate the values of non-state variables (or variables without start values) to get a better initial guess to the Newton algorithm (e.g., spring1.spring.s_rel and spring2.spring.f can be directly calculated via the regular equations and the state variables)

In any case, allowing more control over the initialization root finding algorithm would really help for some classes of problems.

[crupp@…]

LOG_NLS_V output

[Francesco Casella]

Initializes correctly in OpenModelica-v1.13.0-dev-63-g453e1c7-64bit