Hello @all,

I hit the limitation of free evaluation version of dymola with model listed below.
Since I have got a similar model working in 2d, I think it's correct.
As the matter of fact Dymola reports it successfull, but refuses to simulate it,
because "it's too complex for the evaluation version".

I was trying to simulate it with openmodelica 1.4.5.
checkModel(twobody_3d) is successfull => it has 19 equations and 19 variables.
But simulate fails, with no indication of what went wrong.
Can somebody take a look at the model?

My plan was to "print" the variables into a file, which I can open/read from
Octave and show the masses moving around.
I was surprised that the "dymola-way" didn't work.
Doesn't openmodelica have function "modelica.utlities.streams.print"?

I also was surprised to find out that there is just arctan function but no atan2
by the way, arctan is very very unusual to use, it took me also some time
to figure out, that is not called just atan.
But that's not big issue, I wrote the replacement myself.

Thanks in advance, daniel



model twobody_3d

import Modelica.Math.*;
import Modelica.Utilities.*;

parameter Real m1(unit = "kg") = 10;
parameter Real m2(unit = "kg") = 10;

parameter Real p_start[3](unit = "m") = {0,0,0};
parameter Real q_start[3](unit = "m") = {10,0,0};

parameter Real gamma(unit = "N*m^2/kg^2") = 100;

Real dx;
Real dy;
Real dz;
Real phi;
Real theta;
Real d;
Real force;

Real p[3](start = p_start);
Real q[3](start = q_start);
Real vp[3](start = {0,0,0});
Real vq[3](start = {0,0,0});

equation
vp = der(p);
vq = der(q);

dx = q[1] - p[1];
dy = q[2] - p[2];
dz = q[3] - p[3];

d = sqrt(dx^2 + dy^2 + dz^2);

phi = arctan2(dy, dz);
theta = arctan2(dz, sqrt(dx^2 + dy^2));
// phi = atan2(dy, dz);
// theta = atan2(dz, sqrt(dx^2 + dy^2)

force = gamma*m1*m2/d^2;

-m1*der(vp[1]) + force*cos(theta)*cos(phi) = 0;
-m1*der(vp[2]) + force*cos(theta)*sin(phi) = 0;
-m1*der(vp[3]) + force*sin(theta) = 0;

-m2*der(vq[1]) - force*cos(theta)*cos(phi) = 0;
-m2*der(vq[2]) - force*cos(theta)*sin(phi) = 0;
-m2*der(vq[3]) - force*sin(theta) = 0;

//when sample(0,0.01) then
// Streams.print(String(time) + " " +
// String(p[1]) + " " + String(p[2]) + " " + String(p[3]) +
// String(q[1]) + " " + String(q[2]) + " " + String(q[3]),
// "F:/modelica/twobody/data3d");
//end when;

end twobody_3d;

Hi,

You should check all the *.log files when
a simulation fails without any error.
Especially output.log.

We don't support Modelica.Utilities, but is easy to
process the simulation result file which is text for
the same purpose. See twobody_3d_res.plt
For example here is a script that transforms the results in .csv format:
http://article.gmane.org/gmane.comp.misc.openmodelica/517


Here are the files I used:
//----- script twobody_3d.mos
loadModel(Modelica);
getErrorString(); // check if Modelica was loaded OK
loadFile("twobody_3d.mo");
getErrorString(); // check if file was loaded OK
instantiateModel(twobody_3d);
getErrorString(); // are there any error messages?
// tolerance has to be quite large not
// to go into problems with dassl.
simulate(twobody_3d, tolerance=1e-2);
getErrorString(); // are there any error messages?
readFile("output.log"); // read output.log to see if there are any simulation errors
//-----

//---- file twobody_3d.mo
model twobody_3d

import Modelica.Math.*;
import Modelica.Utilities.*;

parameter Real m1(unit = "kg") = 10;
parameter Real m2(unit = "kg") = 10;

parameter Real p_start[3](unit = "m") = {0,0,0};
parameter Real q_start[3](unit = "m") = {10,0,0};

parameter Real gamma(unit = "N*m^2/kg^2") = 100;

Real dx;
Real dy;
Real dz;
Real phi;
Real theta;
Real d;
Real force;

Real p[3](start = p_start);
Real q[3](start = q_start);
Real vp[3](start = {0,0,0});
Real vq[3](start = {0,0,0});

equation
vp = der(p);
vq = der(q);

dx = q[1] - p[1];
dy = q[2] - p[2];
dz = q[3] - p[3];

d = sqrt(dx^2 + dy^2 + dz^2);

// phi = arctan2(dy, dz);
// theta = arctan2(dz, sqrt(dx^2 + dy^2));
phi = atan2(dy, dz);
theta = atan2(dz, sqrt(dx^2 + dy^2));

force = gamma*m1*m2/d^2;

-m1*der(vp[1]) + force*cos(theta)*cos(phi) = 0;
-m1*der(vp[2]) + force*cos(theta)*sin(phi) = 0;
-m1*der(vp[3]) + force*sin(theta) = 0;

-m2*der(vq[1]) - force*cos(theta)*cos(phi) = 0;
-m2*der(vq[2]) - force*cos(theta)*sin(phi) = 0;
-m2*der(vq[3]) - force*sin(theta) = 0;

//when sample(0,0.01) then
// Streams.print(String(time) + " " +
// String(p[1]) + " " + String(p[2]) + " " + String(p[3]) +
// String(q[1]) + " " + String(q[2]) + " " + String(q[3]),
// "F:/modelica/twobody/data3d");
//end when;

end twobody_3d;
//-----

Cheers,
Adrian Pop/

Thank you Adrian for fast helping answer,

everything works fine!