-1

Following this question, I'm trying to generate two time-dependent random functions omega1 and tau using this example. The difference is that I need to have two different sample periods of 0.05 and 0.17 for omega1 and tau respectively. I just duplicated the parts I thought would do the job:

model testData

  extends Modelica.Icons.Example;
  import Modelica.Math.Random.Generators;
  import Modelica.Math.Random.Utilities;

  parameter Real k = 50.0;
  parameter Real J = 0.001;
  Real theta1;
  Real theta2;
  Real omega2;
  
  parameter Modelica.SIunits.Period samplePeriod1 = 0.05;
  parameter Integer globalSeed1 = 30020;
  parameter Integer localSeed1 = 614657;
  output Real omega1;
  
  parameter Modelica.SIunits.Period samplePeriod2 = 0.17;
  parameter Integer globalSeed2 = 30020;
  parameter Integer localSeed2 = 614657;
  output Real tau;

protected
  discrete Integer state1024[33](each start=0, each fixed = true);
    
algorithm
  when initial() then
    state1024 := Generators.Xorshift1024star.initialState(localSeed1, globalSeed1);
    omega1 := 0;
  elsewhen sample(0, samplePeriod1) then
    (omega1, state1024) := Generators.Xorshift1024star.random(pre(state1024));
    omega1 := (omega1 - 0.5) * 13;
  end when;
  
  when initial() then
    state1024 := Generators.Xorshift1024star.initialState(localSeed2, globalSeed2);
    omega1 := 0;
  elsewhen sample(0, samplePeriod2) then
    (tau, state1024) := Generators.Xorshift1024star.random(pre(state1024));
    tau := (tau - 0.5) * 3;
  end when;
  
public
  parameter Integer id1 = Utilities.initializeImpureRandom(globalSeed1);
  discrete Real rImpure1;
  Integer iImpure1;
  
  parameter Integer id2 = Utilities.initializeImpureRandom(globalSeed2);
  discrete Real rImpure2;
  Integer iImpure2;
  
algorithm
  when initial() then
    rImpure1 := 0;
    iImpure1 := 0;
  elsewhen sample(0, samplePeriod1) then
    rImpure1 := Utilities.impureRandom(id=id1);
    iImpure1 := Utilities.impureRandomInteger(
          id=id1,
          imin=-1234,
          imax=2345);
  end when;
  
  when initial() then
    rImpure2 := 0;
    iImpure2 := 0;
  elsewhen sample(0, samplePeriod2) then
    rImpure2 := Utilities.impureRandom(id=id2);
    iImpure2 := Utilities.impureRandomInteger(
          id=id2,
          imin=-1234,
          imax=2345);
  end when;

initial equation
  theta1 = 0;
  theta2 = 0;
  der(theta2) = 0;

equation
  der(theta1) = omega1;
  der(theta2) = omega2;
  J * der(omega2) = tau + k * (theta1 - theta2);

annotation(experiment(StartTime = 0, StopTime = 10, Tolerance = 1e-6, Interval = 0.02));

end testData;

however I get the error messages:

Symbolic Error

The given system is mixed-determined. [index > 3]

Please checkout the option "--maxMixedDeterminedIndex".

Translation Error

No system for the symbolic initialization was generated

I would appreciate if you could help me know what is the problem and how I can solve it.

P.S. considering that this code is apparantly compiling fine on Dymola, this could be a problem with OpenModelica. So I'm adding th JModelica tag in the case those guys can help me know if this compiles over there or not.

Community
  • 1
  • 1
Foad S. Farimani
  • 12,396
  • 15
  • 78
  • 193
  • 1
    The code works fine in Dymola. So it seems to be something tool-related. Looking at the code a bit closer the nearly identical `when(initial)` sections seem strange. – Markus A. Sep 11 '19 at 09:58
  • I don't have a license for Dymola. any idea how I can solve it in either OpenModelica or JMoledlica? – Foad S. Farimani Sep 11 '19 at 10:00

1 Answers1

3

You have omega1 := 0; in two when initial()statements. Replace it by tau := 0; in the second one and the example will work.

I recommend to cleanup your code a bit. I found various smaller issues and needless code lines.

  • everything related to the impure random numbers can be removed
  • localSeed2 and globalSeed2 are useless when they are initialized like the other seed variables
  • state1024 is initialized at 3 different places (even though it works with OpenModelica): with start values and fixed=true and in two different when initial() statements
  • omega2 and tau2 don't need to be outputs. The Tool determines by itself what it has to compute.
  • And finally: Modelica models are a lot easier to debug and understand if existing blocks and physical components are used instead of writing lengthy code in a single class. Your model can also be built graphically with blocks from Modelica.Blocks.Noise and components from Modelica.Mechanics.Rotational.

Below is an updated version of your code with units, only one section for initialization and removed algorithm section (not necessary anymore due to the additional variables rand_omega and rand_tau).

model testData2

  extends Modelica.Icons.Example;
  import Modelica.Math.Random.Generators;
  import Modelica.Math.Random.Utilities;
  import SI = Modelica.SIunits;

  parameter SI.RotationalSpringConstant k = 50.0;
  parameter SI.Inertia J = 0.001;

  parameter SI.Period samplePeriod_tau = 0.17;
  parameter SI.Period samplePeriod_omega = 0.05;

  parameter Integer globalSeed = 30020;
  parameter Integer localSeed_tau = 614657;
  parameter Integer localSeed_omega = 45613;

  SI.Angle theta1, theta2;
  SI.AngularVelocity omega1, omega2, rand_omega;
  SI.Torque tau, rand_tau;

protected 
  discrete Integer state1024_tau[33];
  discrete Integer state1024_omega[33];

initial equation 

  state1024_omega = Generators.Xorshift1024star.initialState(localSeed_omega, globalSeed);
  state1024_tau = Generators.Xorshift1024star.initialState(localSeed_tau, globalSeed);

  theta1 = 0;
  theta2 = 0;
  der(theta2) = 0;

equation 

  when sample(0, samplePeriod_omega) then
    (rand_omega, state1024_omega) = Generators.Xorshift1024star.random(pre(state1024_omega));
  end when;

  when sample(0, samplePeriod_tau) then
    (rand_tau, state1024_tau) = Generators.Xorshift1024star.random(pre(state1024_tau));
  end when;

  der(theta1) = omega1;
  der(theta2) = omega2;

  omega1 = (rand_omega - 0.5) * 13;
  tau = (rand_tau - 0.5) * 3;

  J * der(omega2) = 0 + k * (theta1 - theta2);

annotation(experiment(StartTime = 0, StopTime = 10, Tolerance = 1e-6, Interval = 0.02));
end testData2;
marco
  • 5,944
  • 10
  • 22
  • Thanks a lot for your kind support. some points: 1. will you be kind to include the correct code? I'm not sure how the `state1024` should be implemented. 2. The reason I prefer code is that It is more convenient for me to use, plus later when I want to use JModelica I can bind it to Python... 3. shouldn't I have two different seeds fro the two random functions? they will not be the same with similar seeds? – Foad S. Farimani Sep 11 '19 at 10:48
  • 1
    1. done. 2. Never tried it, but I am pretty sure that JModelica runs also Modelica models which have been built graphically. The graphical information is contained in annotations, which are removed during translation. 3. Two different seeds yes, but not by using the same value twice. – marco Sep 11 '19 at 12:22
  • There is one small issue. The results I get from your code is significantly different from one dited verion I made base on your comments. see my code [here](https://gist.github.com/Foadsf/b0447e8d772b16fe8c78ece99f21e7f5). Do you also observe such a big difference? any idea why is that? – Foad S. Farimani Sep 11 '19 at 12:35
  • 1
    I would say because of different random numbers. Your code uses a single random generator to compute `omega1` and `tau`, whereas mine has two independent random generators which can be initialized with their respective seeds. – marco Sep 11 '19 at 12:47
  • I think the main reason is that the `displayUnit` of `Modelica.SIunits.Angle` is in degrees. but in my code it was radians. – Foad S. Farimani Sep 11 '19 at 13:15
  • 1
    Sure, that makes visually a large difference :) – marco Sep 12 '19 at 05:50
  • I would like also to invite you to the [Modelica Discord](https://discordapp.com/invite/bp2yeYU). Modelica community could benefit if we had better communication methods. – Foad S. Farimani Sep 12 '19 at 06:38