Unit issue in Modelica.Mechanics.Translational.Components.ElastoGap

[HansOlsson]

The commit 0315615 changed ElastoGap to be unit-unsafe. In particular the previous s_ref parameter was removed, and the non-linear exponent n=0 was replaced with the general exponent n=1 (roughly it replaced (s_rel-s_rel0)*abs((s_rel-s_rel0)/s_ref)^n by (s_rel-s_rel0)^n).

The solution seems to be to restore s_ref (e.g., parameter SI.Position s_ref=1 annotation(enable=n>1 or n<1);) and replace c*abs(s_rel - s_rel0)^n by c*abs((s_rel - s_rel0)/s_ref)^n*s_ref. The extra factor s_ref ensures that c keeps its previous unit (and especially doesn't have a unit depending on the real parameter n).

My assumption is that if someone is investigating non-linear elasto gaps the additional freedom of setting s_ref is actually helpful.

[HansOlsson]

The unit-errors are:
The current model works if n=1, but for a general n the factor (s_rel-s_rel0)^n doesn't make sense in terms of units. For integer values it would make sense, but would require a change of unit for c (except for n=1).

[henrikt-ma]

parameter SI.Position s_ref=1 annotation(enable=n>1 or n<1);) 

Is conditioning on the typically non-evaluated, Real-valued, n such a good idea? Isn't it better to use the description of s_ref to say that it only matters for n != 1?

The solution seems to be to restore s_ref (e.g., parameter SI.Position s_ref=1 annotation(enable=n>1 or n<1);) and replace c*abs(s_rel - s_rel0)^n by c*abs((s_rel - s_rel0)/s_ref)^n*s_ref. The extra factor s_ref ensures that c keeps its previous unit (and especially doesn't have a unit depending on the real parameter n).

While multiplying with s_ref would "restore" the unit, wouldn't it result in an unexpected scaling effect for other values than 1?

[tobolar]

What about to introduce a unit conversion like:

constant SI.Length convUnitLenght = 1;
Real e_delta;

e_delta = (s_rel - s_rel0)/convUnitLenght;
f_c = smooth(1, noEvent(if contact then -c*convUnitLenght*abs(e_delta)^n else 0));

[HansOlsson]

What about to introduce a unit conversion like:

constant SI.Length convUnitLenght = 1;
Real e_delta;

e_delta = (s_rel - s_rel0)/convUnitLenght;
f_c = smooth(1, noEvent(if contact then -c*convUnitLenght*abs(e_delta)^n else 0));

That would also work, but I thought exposing it as a parameter (as it was a long time ago) would be more useful; for the non-linear case.

[tobolar]

@HansOlsson You mean s_ref?

I agree with @henrikt-ma

While multiplying with s_ref would "restore" the unit, wouldn't it result in an unexpected scaling effect for other values than 1?

and referring the model documentation, this wouldn't be the intended behaviour, right?

[HansOlsson]

@HansOlsson You mean s_ref?

I agree with @henrikt-ma

While multiplying with s_ref would "restore" the unit, wouldn't it result in an unexpected scaling effect for other values than 1?

and referring the model documentation, this wouldn't not the intended behaviour, right?

It's more complicated.

For me the ideal solution is something like:

  parameter SI.Force cNew;
  parameter SI.Length s_ref = 1;
  Real e_delta(unit="1");
equation
  e_delta = (s_rel - s_rel0)/s_ref;
  f_c = smooth(1, noEvent(if contact then -cNew*abs(e_delta)^n else 0));

This is a clean solution (and likely the only clean solution):

• No unit-trickery.
• All parameters are physical with an appropriate unit:
  • The s_ref is a physical parameter representing the appropriate scale for the compression.
  • The cNew represent the force at the reference scale (when s_rel=s_rel0+s_ref).
  • (For a linear spring only cNew/s_ref matters).
• It is backwards compatible, except that cNew has a new unit compared to c (but the values work the same if you don't change s_ref; and as default s_ref is of course not changed, since it is a new parameter).

Using a conversion factor so that users give a value in the unit in "N/m" and we then internally convert that value to "N/m2" (for n=2) is just subverting the unit-system, and gives two problems:

• It gives the impression that unit-checking is something that shall be circumvented, not something that help users. We have previously seen that this often hides additional issues in models - so this is clearly not something we should encourage users to do (and not do ourselves).
• With displayUnit one can provide values in other units (and that is one of the intended uses), but if a user switches to inches and pound-force the conversion isn't scaling invariant when n>1; so it's very unclear what will happen.

Thus in a future major version we should give users access to the clean solution, and to me the question is only how to handle the intermediate step as cleanly as possible in a minor release. One possibility is to have:

  parameter SI.Force cNew = c*s_ref annotation(Dialog(group="Non linear"));
  parameter SI.Length s_ref = 1 annotation(Dialog(group="Non linear"));
  parameter Real(unit="N/m") c;

[tobolar]

@HansOlsson I like your approach. It simply changes f = c * delta_x^n to f = f_ref * k^n. Consequently,

• cNew shall better be called like f_ref since it is indeed a force. Adding some examples in the documentation would be helpful - I can care about.

Regarding current temporary clean solution:

• shouldn't be cNew = c*s_ref^n better?
• why not SI.TranslationalSpringConstant c?

We can optionally put cNew and s_ref into "Advanced" tab for now.

[HansOlsson]

@HansOlsson I like your approach. It simply changes f = c * delta_x^n to f = f_ref * k^n. Consequently,

• cNew shall better be called like f_ref since it is indeed a force. Adding some examples in the documentation would be helpful - I can care about.

Good.

Regarding current temporary clean solution:

• shouldn't be cNew = c*s_ref^n better?

Only f_ref = c*s_ref is unit-correct with current c. I agree that the other formula would be better, except for the compatibility.

• why not SI.TranslationalSpringConstant c?

Yes, I just didn't remember that name.

We can optionally put cNew and s_ref into "Advanced" tab for now.

Good enough.