Opened 4 years ago

Closed 4 years ago

Last modified 4 years ago

#5927 closed defect (fixed)

Viewing State Variables of Linearised Model

Reported by: kiko.w2g@… Owned by: Karim.Abdelhak
Priority: blocker Milestone: 1.16.0
Component: Run-time Version: v1.16.0-dev
Keywords: State Variables, Linearisation Cc:

Description

Raised ticket to add state variable information to the output file of the linearisation of a model. Currently the linearisation file [from linearize(ModelName) in CLI] outputs the x0 vector around the linearisation point but does not inform what the variable name of these states are corresponding to the model. Without adding the flag -d=stateselection in Tools|Options|Simulation|Additional_Translation_Flags, the user cannot know the states of the system. Also note that the output of the flag does not match the order of the x0 vector from the linearisation file.

Attachments (1)

linearize_test.mos (502 bytes) - added by casella 4 years ago.

Download all attachments as: .zip

Change History (45)

comment:1 Changed 4 years ago by casella

  • Component changed from OMEdit to Backend
  • Owner changed from adeas31 to Karim.Abdelhak
  • Status changed from new to assigned

We had a discussion today in the developers' meeting about this. @Karim.Abdelhak will look into that, the easiest solution is to add some comment lines to the Matlab function with the list of input, state, and output variable names.

Would that be ok for you?

comment:2 Changed 4 years ago by Karim.Abdelhak

I added some changes in PR829, just need to check if it passes the testsuite. Probably i need to adapt some linearization tests.

Unfortunately i could not return it as an array of strings as i originally wanted, since string arrays have some general problems right now. They will be dumped as a plain comment in the model description, but can be easily adapted once we fix the string arrays.

I tested it with multiple models, here is an example.
Original Model:

model simple_test
 Real x1(start=1);
 Real x2(start=2);
 parameter Real a=6,b=2,c=4;
 input Real u = sin(0);
 output Real y;
equation
 der(x1) = x1*(a-b*x1-x2);
 der(x2) = x2*(c-x1-x2);
 y = x1 * u + x2 * u;
end simple_test; 

Linearization:

model linear_simple__test
  "states[2] = [x1, x2]
   inputs[1] = [u]
   outputs[1] = [y]"
  parameter Integer n = 2 "number of states";
  parameter Integer p = 1 "number of inputs";
  parameter Integer q = 1 "number of outputs";
  parameter Integer nz = 2 "data recovery variables";

  parameter Real x0[2] = {1, 2};
  parameter Real u0[1] = {0};
  parameter Real z0[2] = {0, 0};

  parameter Real A[n, n] = [-2.665600749850023e-07, -1.000000001354202; -2.000000002708404, -1.00000010131423];
  parameter Real B[n, p] = [0; 0];
  parameter Real C[q, n] = [0, 0];
  parameter Real D[q, p] = [3];
  parameter Real Cz[nz, n] = [0, 0; 0, 0];
  parameter Real Dz[nz, p] = [1; 3];

  Real x[n](start=x0);
  input Real u[p](start=u0);
  output Real y[q];
  output Real z[nz];

  Real 'x_x1' = x[1];
  Real 'x_x2' = x[2];
  Real 'u_u' = u[1];
  Real 'y_y' = y[1];
  Real 'z_u' = z[1];
  Real 'z_y' = z[2];
equation
  der(x) = A * x + B * u;
  y = C * x + D * u;
  z = Cz * x + Dz * u;
end linear_simple__test;

Is this what you would expect?
(If you want to pull the branch from the PR to test your own stuff just tell me!)

comment:3 Changed 4 years ago by Karim.Abdelhak

Only linearization models were affected and it broke nothing so i pushed the changes.

Feel free to close the ticket if this is what you wanted!

comment:4 follow-up: Changed 4 years ago by adrpo

I'm not sure if we test OMPython and OMJulia usage of linearization API, hopefully those still work.

comment:5 follow-up: Changed 4 years ago by casella

Looks good, using aliases for inputs, states and outputs with the variable names is a good idea. I assume in case of hierarchical models you will substitute the "." with a "_".

However, I'm not sure what you mean by "recovery variables". Are we computing all the internal variables besides the outputs? For models of real-life size this may actually blow up the linearized model. Our average power plant models has 4-8 inputs, 6-12 outputs, 100-300 states, and maybe over 20,000 variables. Matrices Cz and Dz will be huge.

At least I would make the output of such variables optional and off by default.

comment:6 in reply to: ↑ 4 Changed 4 years ago by Karim.Abdelhak

Replying to adrpo:

I'm not sure if we test OMPython and OMJulia usage of linearization API, hopefully those still work.

I only added a header comments to the generated files containing the information about states, inputs outputs etc. It should be alright...

But shouldn't we add test cases for OMPython and OMJulia in general here?

comment:7 in reply to: ↑ 5 ; follow-up: Changed 4 years ago by Karim.Abdelhak

Replying to casella:

Looks good, using aliases for inputs, states and outputs with the variable names is a good idea. I assume in case of hierarchical models you will substitute the "." with a "_".

I don't exactly know what you mean, right now i don't know of such a substitution. I am currently using the component reference name string (the one which is dumped in all of the dumping functions and is unique)

However, I'm not sure what you mean by "recovery variables". Are we computing all the internal variables besides the outputs? For models of real-life size this may actually blow up the linearized model. Our average power plant models has 4-8 inputs, 6-12 outputs, 100-300 states, and maybe over 20,000 variables. Matrices Cz and Dz will be huge.

At least I would make the output of such variables optional and off by default.

I also don't know what you mean by that, i only implemented state, input and output. As far as i can see a duplicate model is generated which is called something like data_recovery_frame. Unfortunately i did not implement this module and don't really know what this is about. I guess the first one might be corrupted due to wrong outputs, it is written right on time during c-runtime from simulation results so the output is not really known to be safe a priori. So another safer method which can have worse outputs creates the second file. The recovery variables could be for this.

But that is just a guess from what i have seen. If you have any issues with these recovery variables i can have a deeper look!

comment:8 Changed 4 years ago by Karim.Abdelhak

Ah sorry i just realized what you mean, yes for alias names the "." get replaced by "_" due to c syntax

comment:9 in reply to: ↑ 7 Changed 4 years ago by casella

Replying to Karim.Abdelhak:

I also don't know what you mean by that, i only implemented state, input and output. As far as i can see a duplicate model is generated which is called something like data_recovery_frame. Unfortunately i did not implement this module and don't really know what this is about. I guess the first one might be corrupted due to wrong outputs, it is written right on time during c-runtime from simulation results so the output is not really known to be safe a priori. So another safer method which can have worse outputs creates the second file. The recovery variables could be for this.

I am completely lost, sorry. The linearized model is a mathematical model, why should it be corrupted and why should some variables need to be recovered?

But that is just a guess from what i have seen. If you have any issues with these recovery variables i can have a deeper look!

Please do. We shouldn't have outputs of our tools that we don't understand :)

comment:10 Changed 4 years ago by casella

In any case, I tried to linearize the following model

model foo
  input Real input1,input2;
  output Real output1;
  Real state1, state2;
  Real z;
  Real w;
equation
  z = input1 - state1;
  w = 2*state2 + input2;
  der(state1) = z;
  der(state2) = state1;
  output1 = w;
end foo;

By default, I get this linear dump in Modelica

model linear_foo "
  states[2] = [state1, state2]
  inputs[2] = [input1, input2]
  outputs[1] = [output1]"
  parameter Integer n = 2 "number of states";
  parameter Integer p = 2 "number of inputs";
  parameter Integer q = 1 "number of outputs";

  parameter Real x0[n] = {0, 0};
  parameter Real u0[p] = {0, 0};

  parameter Real A[n, n] = [-1, 0; 1, 0];
  parameter Real B[n, p] = [1, 0; 0, 0];
  parameter Real C[q, n] = [0, 2];
  parameter Real D[q, p] = [0, 1];

  Real x[n](start=x0);
  input Real u[p](start=u0);
  output Real y[q];

  Real 'x_state1' = x[1];
  Real 'x_state2' = x[2];
  Real 'u_input1' = u[1];
  Real 'u_input2' = u[2];
  Real 'y_output1' = y[1];
equation
  der(x) = A * x + B * u;
  y = C * x + D * u;
end linear_foo;

the comment actually repeats the same information that was already provided in the alias variable names that are defined at the end of the model. I'm not sure we really need that, I would stick to the old format for the Modelica linearized model, since it provides information that is accessible from the browser and not only in a comment.

The goal of this ticket, which was not formulated clearly (sorry about that) was to have the same information also with the other languages. For example, the matlab output is

function [n, p, q, x0, u0, A, B, C, D] = foo_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u 
  n = 2; % number of states 
  p = 2; % number of inputs 
  q = 1; % number of outputs 

  x0 = {0, 0};
  u0 = {0, 0};

  A = [-1, 0; 1, 0];
  B = [1, 0; 0, 0];
  C = [0, 2];
  D = [0, 1];
end}}}
which still doesn't give me any information about who the inputs, states, and outputs are.

Can we add some comments to the Matlab, Julia, and Python output?

comment:11 Changed 4 years ago by casella

  • Milestone changed from Future to 1.16.0
  • Priority changed from low to high

comment:12 follow-up: Changed 4 years ago by Karim.Abdelhak

I removed the old comments and added the following (including artificial hierarchy to showcase the replacement):

function [n, p, q, x0, u0, A, B, C, D] = simple_test_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u
  n = 2; % number of states
  p = 1; % number of inputs
  q = 1; % number of outputs

  x0 = {1, 2};
  u0 = {0};

  A = [-2.665600749850023e-07, -1.000000001354202; -2.000000002708404, -1.00000010131423];
  B = [0; 0];
  C = [0, 0];
  D = [3];

  x_x1_x = x[1];
  x_x2 = x[2];
  u_u = u[1];
  y_y = y[1];
  z_u = z[1];
  z_y = z[2];

It is tested right now, but i guess there will be no problems with it. I also added a testcase for the three dumping options.

Did you mean the z variables with recovery variables? I don't see a reason for them to be there so i could remove them.

comment:13 follow-up: Changed 4 years ago by Karim.Abdelhak

Everything runs fine according to Jenkins, but i think the z variables really do not belong there. They are some kind of inner variables i suppose. I will remove them, everything should be fine then!

I will also try to implement the fix for #5938 in the same PR since it is strongly related.

EDIT: you were right about the recovery variables, they are inner variables to recover the simulation from any point i guess. This however does only matter for the Modelica output and not for all other targets. I will just remove them.

Last edited 4 years ago by Karim.Abdelhak (previous) (diff)

comment:14 in reply to: ↑ 13 Changed 4 years ago by casella

Replying to Karim.Abdelhak:

I will also try to implement the fix for #5938 in the same PR since it is strongly related.

That would be nice :)

EDIT: you were right about the recovery variables, they are inner variables to recover the simulation from any point i guess. This however does only matter for the Modelica output and not for all other targets.

Absolutely!

I will just remove them.

OK.

comment:15 in reply to: ↑ 12 ; follow-up: Changed 4 years ago by casella

Replying to Karim.Abdelhak:

I removed the old comments and added the following (including artificial hierarchy to showcase the replacement):

function [n, p, q, x0, u0, A, B, C, D] = simple_test_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u
  n = 2; % number of states
  p = 1; % number of inputs
  q = 1; % number of outputs

  x0 = {1, 2};
  u0 = {0};

  A = [-2.665600749850023e-07, -1.000000001354202; -2.000000002708404, -1.00000010131423];
  B = [0; 0];
  C = [0, 0];
  D = [3];

  x_x1_x = x[1];
  x_x2 = x[2];
  u_u = u[1];
  y_y = y[1];
  z_u = z[1];
  z_y = z[2];

This looks fine. I'm only a bit worried about the translation rules from Modelica variable names to Matlab variable names. Are we sure we always get legal Matlab variable names, particularly when hierarchically structured models are involved?

comment:16 in reply to: ↑ 15 Changed 4 years ago by Karim.Abdelhak

Replying to casella:

Are we sure we always get legal Matlab variable names, particularly when hierarchically structured models are involved?

I implemented a replacement of "." with "_", as you can see with x_x1_x. That originally was an attribute x of a class object x1, the first prepended x just shows that it is a state. Prepended y is an output and u an input. (zis inner but i removed that).

I need to test some stuff for indexed variables, bit i think x_x(1) = without allocating beforehand and just expanding on the fly afterwards is allowed in matlab. Even for any amount of dimensions as long as the first call of that array has the same amount of dimensions. Inefficient but allowed. I need to replace square brackets with round brackets though.

I don't know if that works for python or julia... any comments on that?

comment:17 Changed 4 years ago by Karim.Abdelhak

I also just realized a fundamental flaw in this ... the generated code will never be actually executable since there is nothing inbetween. The variable x(1) that gets aliased by something else just does not exist without the state space system. I guess the solver could be added by hand but there is something matlab internal that can be used. I fiddled a little bit with it and this is what i came up with:

Original model:

model simple_test
 number num(x(start={1, 2}));
 parameter Real a=6,b=2,c=4;
 input Real u = sin(0);
 output Real y;
 class number
   Real x[2];
 end number;
equation
 der(num.x[1]) = num.x[1]*(a-b*num.x[1]-num.x[2]);
 der(num.x[2]) = num.x[2]*(c-num.x[1]-num.x[2]);
 y = num.x[1] * u + num.x[2] * u;
end simple_test;

generated linear matlab model:

function [sys] = simple_test_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u 
  Ts = 1.0; % stop time
  n = 2; % number of states 
  p = 1; % number of inputs 
  q = 1; % number of outputs 

  x0 = {1, 2};
  u0 = {0};

  A = [-2.665600749850023e-07, -1.000000001354202; -2.000000002708404, -1.00000010131423];
  B = [0; 0];
  C = [0, 0];
  D = [3];

  sys = ss(A,B,C,D,Ts,'StateName',{'x_num_x1','x_num_x2'}, 'InputName',{'u_u'}, 'OutputName', {'y_y'});
end

This creates a state space model, that can be run as is. It returns the following in matlab:

>> sys = linear_simple_test()

sys =
 
  A = 
               x_num_x1    x_num_x2
   x_num_x1  -2.666e-07          -1
   x_num_x2          -2          -1
 
  B = 
             u_u
   x_num_x1    0
   x_num_x2    0
 
  C = 
        x_num_x1  x_num_x2
   y_y         0         0
 
  D = 
        u_u
   y_y    3
 
Sample time: 1 seconds
Discrete-time state-space model.

This can just be used. I don't know how many variables this thing can handle, but i guess this is a good solution. What do you think? (solutions for python and julia are still to be found i guess)

Note that i also turned all spaces between dimension into underscores such that it still is unique. variable num.x[1,3] would be x_num_x1_3.

Last edited 4 years ago by Karim.Abdelhak (previous) (diff)

comment:18 Changed 4 years ago by casella

I'm in favour of turning everything which is not a number or letter into an underscore, and so be it. It just needs to be recognizable, I don't think there's any need to waste too much time here

comment:19 Changed 4 years ago by casella

But in fact, I would go further than that.

In the Modelica model, we defined those aliases, because you may want to simulate the linearized model and look at their transients in the variables browser when the simulation is finished.

In Matlab, there is no such thing, the goal of the function is only to return the ABCD matrices. The system will be simulated starting from those matrices. So, there is no point defining alias variables in that function at all

What is needed is just to know who the inputs, states and outputs are, so that you can interpret the ABCD matrices correctly. This is achieved easily by putting the Modelica names of inputs, outputs, and states in a comment, using dot notation and without any change, since a comment can contain any character.

For example:

function [n, p, q, x0, u0, A, B, C, D] = simple_test_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u
  n = 2; % number of states
  p = 1; % number of inputs
  q = 1; % number of outputs

  x0 = {1, 2};
  u0 = {0};

  A = [-2.665600749850023e-07, -1.000000001354202; -2.000000002708404, -1.00000010131423];
  B = [0; 0];
  C = [0, 0];
  D = [3];

% u(1) = myInput
% y(1) = myOutput
% x(1) = mass1.v
% x(2) = mass2.v
Last edited 4 years ago by casella (previous) (diff)

comment:20 Changed 4 years ago by casella

Also, the standard notation in system theory is n = number of states, m = number of inputs, p = number of outputs. Can you please update that as well?

Thanks!

Last edited 4 years ago by casella (previous) (diff)

comment:21 Changed 4 years ago by casella

If that is not too difficult, I would also recommend to break the lines of matrix A to avoid getting too long lines with systems having hundreds of states, e.g.

  A = [-2.665600749850023e-07, -1.000000001354202; 
       -2.000000002708404, -1.00000010131423];

no need to do that for B, C, D, since normally you don't have more than a handful of inputs and outputs.

comment:22 Changed 4 years ago by casella

I guess the last comments hold also for Julia and Python

comment:23 Changed 4 years ago by Karim.Abdelhak

Sure, i can do all those rather small updates, but what about the system returning?

  sys = ss(A,B,C,D,Ts,'StateName',{'x_num_x1','x_num_x2'}, 'InputName',{'u_u'}, 'OutputName', {'y_y'});

I want to create this thing because it has all the wanted information dense and structured. It looks like this:

sys =
 
  A = 
               x_num_x1    x_num_x2
   x_num_x1  -2.666e-07          -1
   x_num_x2          -2          -1
 
  B = 
             u_u
   x_num_x1    0
   x_num_x2    0
 
  C = 
        x_num_x1  x_num_x2
   y_y         0         0
 
  D = 
        u_u
   y_y    3
 
Sample time: 1 seconds
Discrete-time state-space model.

No mapping needed for this, and you can access stuff inside with e.g. sys.A and also get the state names in the correct order with sys.StateNames which is a cell array.

I would like to use that structure!

comment:24 follow-up: Changed 4 years ago by casella

Then you need the function to return those cell arrays of strings, and define them in the function by putting the full, unabridged Modelica path in there.

function [n, p, q, x0, u0, A, B, C, D, input_names, output_names, state_names] = simple_test_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u
  n = 2; % number of states
  p = 1; % number of inputs
  q = 1; % number of outputs

  x0 = {1, 2};
  u0 = {0};

  A = [-2.665600749850023e-07, -1.000000001354202; -2.000000002708404, -1.00000010131423];
  B = [0; 0];
  C = [0, 0];
  D = [3];

  inputNames = {'myInput'};
  outputNames = {'myOutput'};
  stateNames = {
    'mass1.v',
    'mass2.v'};
[n, p, q, x0, u0, A, B, C, D, input_names, output_names, state_names] = 
  simple_test_GetLinearModel()
sys = ss(A,B,C,D,Ts,'StateName',state_names, 'InputName', input_names,
         'OutputName', output_names);

comment:25 in reply to: ↑ 24 Changed 4 years ago by Karim.Abdelhak

Replying to casella:

Then you need the function to return those cell arrays of strings, and define them in the function by putting the full, unabridged Modelica path in there.

i don't really get it, why can't my generated function just return the sys object with everything else inside? I guess i have to also return the start values and if you insist the number of states, inputs and outputs, but i guess they are implicitely given by the sizes of the structures inside.

If you take a look again at comment:17 this just works as is and everything is already inside the sys object.

comment:26 Changed 4 years ago by casella

Aha, I missed that, sorry. This looks good to me, assuming everybody who uses Matlab has the Control System Toolbox installed. I guess this is probably going to be the case.

Let's go for it!

comment:27 Changed 4 years ago by Karim.Abdelhak

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

This should be fixed. If anything is not as expected, please reopen the ticket.

Changed 4 years ago by casella

comment:28 follow-up: Changed 4 years ago by casella

  • Resolution fixed deleted
  • Status changed from closed to reopened

Unfortunately it doesn't still work well. If I set the language to Matlab, I get the following error

messages = "Succeeding building the linearized executable, but failed to run the 
linearize command: \"D:/Temp/OMEdit/simple_test.exe\" -l=1.0

Without setting the language it works fine.

I'd suggest to add a test case to the testsuite where you generate a linearized model with each of the options, just in case.

comment:29 in reply to: ↑ 28 Changed 4 years ago by Karim.Abdelhak

Replying to casella:

Unfortunately it doesn't still work well. If I set the language to Matlab, I get the following error

messages = "Succeeding building the linearized executable, but failed to run the 
linearize command: \"D:/Temp/OMEdit/simple_test.exe\" -l=1.0

Without setting the language it works fine.

I'd suggest to add a test case to the testsuite where you generate a linearized model with each of the options, just in case.

I did and it works. Could you provide a failing example?

Edit: Ah i just saw the attachment, thanks!

Last edited 4 years ago by Karim.Abdelhak (previous) (diff)

comment:30 Changed 4 years ago by casella

Any idea what is going wrong?

comment:31 Changed 4 years ago by Karim.Abdelhak

Running the file you provided with version 1.16.0~dev-420-gc007a39 works fine for me. This is the output:

true
""
true
""
record SimulationResult
    resultFile = "/home/kab/Repo/FH/promotion/modelicaLibs/linearization/simple_test_res.mat",
    simulationOptions = "startTime = 0.0, stopTime = 1.0, numberOfIntervals = 500, tolerance = 1e-06, method = 'dassl', fileNamePrefix = 'simple_test', options = '', outputFormat = 'mat', variableFilter = '.*', cflags = '', simflags = ''",
    messages = "stdout            | info    | Linearization will be performed at point of time: 1.000000
LOG_SUCCESS       | info    | The initialization finished successfully without homotopy method.
LOG_SUCCESS       | info    | The simulation finished successfully.
stdout            | info    | Linear model is created at /home/kab/Repo/FH/promotion/modelicaLibs/linearization/linear_simple_test.m
stdout            | info    | The output format can be changed with the command line option --linearizationDumpLanguage.
stdout            | info    | The options are: --linearizationDumpLanguage=modelica, matlab, julia, python.
",
    timeFrontend = 0.008347984000000001,
    timeBackend = 0.016981782,
    timeSimCode = 0.002271676,
    timeTemplates = 0.006069824000000001,
    timeCompile = 1.047986665,
    timeSimulation = 0.024659261,
    timeTotal = 1.10679228
end SimulationResult;
"Warning: The initial conditions are not fully specified. For more information set -d=initialization. In OMEdit Tools->Options->Simulation->OMCFlags, in OMNotebook call setCommandLineOptions("-d=initialization").
"

And this is the matlab file linear_simple_test.m:

function [sys, x0, u0, n, m, p] = simple_test_GetLinearModel()
% der(x) = A * x + B * u
% y = C * x + D * u 
  n = 2; % number of states 
  m = 1; % number of inputs 
  p = 1; % number of outputs 

  x0 = {1.797168953922605, 2.295669049250046};
  u0 = {0};
  Ts = 1; % stop time

  A =	[-3.484345076096534, -1.797168955525673;
	-2.295669056449934, -2.388507150275357];

  B =	[0;
	0];

  C =	[0, 0];

  D =	[4.092838003172653];

  % The Control System Toolbox is required for this. Alternatively just return the matrices A,B,C,D instead.
  sys = ss(A,B,C,D,Ts,'StateName',{'x_num_x(1)','x_num_x(2)'}, 'InputName',{'u_u'}, 'OutputName', {'y_y'});
end

Which i verified with Matlab. It works fine, is it maybe a windows problem? Can you verify your fail on a linux machine? Unfortunately i do not have a windows system available right now, i could check this weekend.

Last edited 4 years ago by Karim.Abdelhak (previous) (diff)

comment:32 Changed 4 years ago by casella

  • Component changed from Backend to Run-time
  • Owner changed from Karim.Abdelhak to adrpo
  • Priority changed from high to blocker
  • Status changed from reopened to assigned
  • Type changed from enhancement to defect

I confirm that it works on Linux. I tired your same version of OMC under windows and it failed due to some issues with cvode linking by mingw. On the next day @adrpo committed a fix for that, and then is started to fail as reported in comment:28.

I guess we need to fix this for 1.16.0. I'd ask @adrpo to check as soon as he's done with replaceable.

comment:33 Changed 4 years ago by adrpo

This doesn't work because the C++ code that populates the string does not agree with the variables that are sent to it and then it segfaults. Is more than weird that this works in Linux!

The string in fprintf is from CodegenC.tpl:

    const char *<%symbolName(modelNamePrefix,"linear_model_frame")%>()
    {
      return "function [sys, x0, u0, n, m, p] = <%symbolName(modelNamePrefix,"GetLinearModel")%>()\n"
      "%% der(x) = A * x + B * u\n%% y = C * x + D * u \n"
      "  n = <%varInfo.numStateVars%>; %% number of states \n  m = <%varInfo.numInVars%>; %% number of inputs \n  p = <%varInfo.numOutVars%>; %% number of outputs \n"
      "\n"
      "  x0 = %s;\n"
      "  u0 = %s;\n"
      "  Ts = %g; %% stop time\n"
      "\n"
      <%matrixA%>
      <%matrixB%>
      <%matrixC%>
      <%matrixD%>
      "  %% The Control System Toolbox is required for this. Alternatively just return the matrices A,B,C,D instead.\n"
      "  sys = ss(A,B,C,D,Ts,'StateName',{<%getVarNameMatlab(vars.stateVars, "x")%>}, 'InputName',{<%getVarNameMatlab(vars.inputVars, "u")%>}, 'OutputName', {<%getVarNameMatlab(vars.outputVars, "y")%>});\n"
      "end";
    }

This is what happens in linearize.cpp:

fprintf(
  fout, 
  data->callback->linear_model_frame(),  // this is the string above, see where %s and %g are in there!
  strX.c_str(), // x0 = %s; <-- OK
  strU.c_str(), // u0 = %s; <-- OK
  strA.c_str(), // Ts = %g; %% stop time <-------- Oh, NOOO, are we sending a string to %g?!
  strB.c_str(), 
  strC.c_str(), 
  strD.c_str(), 
  (double) data->simulationInfo->stopTime); // <--- this should be for %g.

One needs to be very careful to make the strings for the linear models in the CodegenC.tpl agree to the values one sends to it from linearize.cpp.

comment:34 Changed 4 years ago by adrpo

The CodegenC.tpl template generates this string (seen in gdb). Is easier to see the %s and %g in here:

(gdb) set print elements 0
(gdb) p $4
$5 = 0x1f269c0 <equationIndexes.56564+232> "function [sys, x0, u0, n, m, p] = simple_test_GetLinearModel()\n
%% der(x) = A * x + B * u\n
%% y = C * x + D * u \n
  n = 2; %% number of states \n
  m = 1; %% number of inputs \n
  p = 1; %% number of outputs \n\n
  x0 = %s;\n
  u0 = %s;\n
  Ts = %g; %% stop time\n\n
  A =\t[%s];\n\n
  B =\t[%s];\n\n
  C =\t[%s];\n\n
  D =\t[%s];\n\n
  %% The Control System Toolbox is required for this. Alternatively just return the matrices A,B,C,D instead.\n
  sys = ss(A,B,C,D,Ts,'StateName',{'x_num_x(1)','x_num_x(2)'}, 'InputName',{'u_u'}, 'OutputName', {'y_y'});\n
end"

comment:35 Changed 4 years ago by adrpo

I guess the easy fix is to move: Ts = %g; %% stop time\n\n after the D matrix.

comment:36 Changed 4 years ago by adrpo

Indeed, that fixed it: https://github.com/OpenModelica/OpenModelica/pull/939

adrpo33@ida-0030 MINGW64 /c/home/adrpo33/dev/OMTesting/linearize
$ ~/dev/OpenModelica/build/bin/omc linearize_test.mos
true
""
true
""
record SimulationResult
    resultFile = "C:/home/adrpo33/dev/OMTesting/linearize/simple_test_res.mat",
    simulationOptions = "startTime = 0.0, stopTime = 1.0, numberOfIntervals = 500, tolerance = 1e-006, method = 'dassl', fileNamePrefix = 'simple_test', options = '', outputFormat = 'mat', variableFilter = '.*', cflags = '', simflags = ''",
    messages = "stdout            | info    | Linearization will be performed at point of time: 1.000000
LOG_SUCCESS       | info    | The initialization finished successfully without homotopy method.
LOG_SUCCESS       | info    | The simulation finished successfully.
stdout            | info    | Linear model is created at C:\\home\\adrpo33\\dev\\OMTesting\\linearize/linear_simple_test.m
stdout            | info    | The output format can be changed with the command line option --linearizationDumpLanguage.
stdout            | info    | The options are: --linearizationDumpLanguage=modelica, matlab, julia, python.
",
    timeFrontend = 0.0076419,
    timeBackend = 0.0114338,
    timeSimCode = 0.0014409,
    timeTemplates = 0.06275409999999999,
    timeCompile = 9.3955176,
    timeSimulation = 0.3017144,
    timeTotal = 9.7812082
end SimulationResult;
"Warning: The initial conditions are not fully specified. For more information set -d=initialization. In OMEdit Tools->Options->Simulation->OMCFlags, in OMNotebook call setCommandLineOptions("-d=initialization").
"

comment:37 Changed 4 years ago by casella

Thank @adrpo! I'll check that out on the nightly tomorrow

comment:38 Changed 4 years ago by casella

OK, now the .m file is successfully generated, thanks @adrpo for the fix!

I would kindly ask @Karim to fix some last details.

If I linearize model P.simple_test, I get a file linear_P.simple_test.m. Besides this filename potentially being very, very long if the model is found within several nested packages with long names, MATLAB file names can only contain letters, numbers and underscores, so this file name is illegal for all models that are not defined in the top global scope.

Furthermore, the function name does not coincide with the file name, which is also not recommended.

I would suggest to simply name the function get_linearized_model() and the corresponding file get_linearized_model.m; this would solve all these problems at once.

A similar argument applies when generating the linearized model in Modelica. If I linearize model P.simple, I get a file named linear_P.simple_test.mo which contains a model named linear_simple__test, where the dot becomes an underscore and the underscore becomes a double underscore. The Modelica Specification Sect. 13.2.2.2 says:

When mapping a package or class-hierarchy to a file [e.g. the file A.mo], that file shall only define a single class [A] with a name matching the name of the nonstructured entity. In a file hierarchy the files shall have the extension ".mo"

Hence, this is also technically illegal, even though OMC loads it anyway.

I would suggest to call the model Linearized and the file Linearized.mo instead.

I'm not sure about the Python and Julia model, because I do not have experience using those languages, but I'd stick to the same guidelines.

In general, one can generate multiple different linearized models from the same Modelica model, by changing the point of linearization, so trying to make them unique by using the original model name in the linearized model doesn't work anyway. I think it is then much simpler and easier to understand if apply the KISS principle and avoid using that name at all.

Last edited 4 years ago by casella (previous) (diff)

comment:39 Changed 4 years ago by casella

  • Owner changed from adrpo to Karim.Abdelhak

@Karim, there is also a small but absolutey critical bug to fix, besides the naming. In the Matlab output you wrote:

sys = ss(A,B,C,D,Ts,...)

If you add Ts, Matlab interprets this a the sampling time of a discrete-time system, and thus returns a discrete-time object. Which is totally wrong, because of course the meaning of the A,B,C,D matrices is completely different in that case. I was trying to figure out why I was getting unstable responses when computing step(sys), since the eigenvalues were all in the left-hand plane, but then I plotted the Bode plots, and they were limited to Nyquist's frequency, so I understood: the eigenvalues were outside the unit disk...

This really needs to be fixed. I would suggest you also fix the file and function naming, see comment:38, otherwise this won't work out of the box, but would require some file editing, which is not nice.

Thanks!

Last edited 4 years ago by casella (previous) (diff)

comment:40 Changed 4 years ago by Karim.Abdelhak

I removed the stopTime and just return it individually now. Probably not really relevant for the linearized model, but why not. It also states that it is a continuou-time state-space model now.

Furthermore, i changed the stuff you recommended about the file name, in general all files and all functions are now called linearized_model and have a comment which is the original model name.

I have to change some test files because of this, but should be done today. If you have some last requests or objections, please tell me! :)

comment:41 Changed 4 years ago by casella

Sounds good :)

comment:42 Changed 4 years ago by Karim.Abdelhak

The changes are in PR984, one model somehow fails... it has different results from mine locally. As far as i can see my results are correct and the ones jenkins computes seem to be wrong, maybe you can give some intel here.

Everything else seems to work fine.

comment:43 Changed 4 years ago by Karim.Abdelhak

We figured out the problem, hopefully it is easy to fix.

comment:44 Changed 4 years ago by Karim.Abdelhak

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

For later reference: It was a problem regarding parallelization of test files. Since we now generate multiple files with the same name, the linking in parallel interfered with one another.

We modified an old ugly hack to be the new ugly hack but we made nothing worse. There should be a better solution to this problem but that is another problem in general.

This problem is fixed! (for good hopefully)

Last edited 4 years ago by Karim.Abdelhak (previous) (diff)
Note: See TracTickets for help on using tickets.