Christophe Le Bris <chris.lebris@gmail.com> a écrit :
Hello,
Hi Christophe,
I tried to reproduce the behavior exposed in my previous email with a
minimal test.
I integrate numerically a trajectory with only central force and an
additional equation that increases the value of an additional state of
1/sec. The "Additional equation" should also increase the SMA of 2m/sec
but, after integration for the duration of one keplerian orbit, my
additional state changed to the value of the initial orbital period but
the
orbit did not changed => the increase of 2m/sec has not been taken into
account.
I confirm the behaviour you get, and also confirm your thorough analysis.
The latest version on Hipparchus indeed does not allow the secondary
equation to modify the main equation.
However, I am on the fence with it. Should we really allow an additional
equation to modify the main equation? If I remember well, the fact that the
derivatives of the main state are provided as an argument is mainly because
they may be needed to compute the secondary derivatives, it was not because
users could modify the primary derivatives.
This was however not obvious in the Apache Commons Math era because the
signature of the method was anyway to provide arrays and have the user
fill them. So the interdace SecondaryEquations did specify the
computeDerivatives method as:
/** Compute the derivatives related to the secondary state parameters.
* @param t current value of the independent <I>time</I> variable
* @param primary array containing the current value of the primary
state vector
* @param primaryDot array containing the derivative of the primary
state vector
* @param secondary array containing the current value of the secondary
state vector
* @param secondaryDot placeholder array where to put the derivative of
the secondary state vector
* @exception MaxCountExceededException if the number of functions
evaluations is exceeded
* @exception DimensionMismatchException if arrays dimensions do not
match equations settings
*/
void computeDerivatives(double t, double[] primary, double[] primaryDot,
double[] secondary, double[] secondaryDot)
throws MaxCountExceededException, DimensionMismatchException;
As you see, primaryDot and secondaryDot are passed the same way, despite
only the
secondaryDot is really intended to be filled by user.
Since Hipparchus, the design is that arrays the user should fill must be
return values (i.e. allocated within the method), so the equivalent
interface,
which is now called SecondaryODE defines:
/** Compute the derivatives related to the secondary state parameters.
* @param t current value of the independent <I>time</I> variable
* @param primary array containing the current value of the primary
state vector
* @param primaryDot array containing the derivative of the primary
state vector
* @param secondary array containing the current value of the
secondary state vector
* @return derivative of the secondary state vector
* @exception MathIllegalStateException if the number of functions
evaluations is exceeded
* @exception MathIllegalArgumentException if arrays dimensions do not
match equations settings
*/
double[] computeDerivatives(double t, double[] primary, double[]
primaryDot, double[] secondary)
throws MathIllegalArgumentException, MathIllegalStateException;
It makes more clear than only the derivatives of the secondary state
vector are expected
from this method. The primaryDot array is still passed because it may be
needed in
to compute the other derivatives, but is was not intended that users
change it.
For the sake of separation of concerns, I would prefer to have realy the
primary
equations completely computed by one part of the code, and secondary
equations
by another part and not messing with the primary. So if you remove the
secondary
equations, you only don't get the secondary state they compute, but you do
not change
at all the premary equations. When the secondary equations do change the
primary
derivatives, they are not secondary anymore: they are needed to compute the
primary state.
I understand that moving the mapper.insertEquationData call corresponding
to the main
state from before the computation of the secondary equations to after the
computation
allows as a side effect changes made by the secondary equations to be
taken into
account, but it is not really clean. It implies that the primary equations
were
indeed incomplete, and it implies the secondary equation returns two
different results,
on in its regular return value, and one by changing a value in an array
passed by the caller,
and expected this array to be used only after the call, not before.
In your use case, would it be possible to split your computation in two
parts, one in
the secondary equation (to compute the return value as in your return new
double[]
{ 2, 0, 0, 0, 0, 0, 0 } example) and another one as a force model (to
compute the
main state derivative as in your pDot[0] = 1 example)?
best regards,
Luc
*Initial orbit = Keplerian parameters: {a: 7800000.0; e: 0.001; i: 28.0;
pa: 0.0; raan: 0.0; v: 0.0;}Initial keplerian period = 6855.717044697212
secFinal value of added state = 6855.717044697215Final orbit = Keplerian
parameters: {a: 7800000.000000003; e: 0.001000000000000334; i:
28.000000000000007; pa: -1.1173044601630712E-25; raan:
2.925205344836198E-29; v: -5.234453515214931E-13;}*
I expect the problem to be in Hipparchus 1.2 (at least) and not Orekit. In
*org.hipparchus.ode.ExpandableODE*, I think the line 137
"mapper.insertEquationData(index, primaryStateDot, yDot);" should be
placed
just before the "return".
If I execute the same program with Orekit 7.0/Commons Math 3.4.1, the
behavior is different and is ok:
*Initial orbit = keplerian parameters: {a: 7800000.0; e: 0.001; i: 28.0;
pa: 0.0; raan: 0.0; v: 0.0;}Initial keplerian period = 6855.717044697212
secFinal value of added state = 6855.717044697215Final orbit = keplerian
parameters: {a: 7813711.434089395; e: 0.0010000000000002728; i:
27.999999999999993; pa: -1.3080031127518713E-13; raan:
3.7263393091775305E-17; v: -0.4748813935419053;}*
Thank you!
Christophe
---------- Forwarded message ----------
From: Christophe Le Bris <chris.lebris@gmail.com>
Date: 2017-11-30 18:21 GMT+01:00
Subject: Problem with additional equations
To: "orekit-users@orekit.org" <orekit-users@orekit.org>
Hello,
I have a problem when integrating "additional equations" that impact the
derivatives of the primary state (ie: value returned by
"computeDerivative"
method is not null).
I have the feeling that these values are not taken into account during the
integration.
I use Orekit 9.1 with Hipparchus 1.2.
It seems that the origin of my problem is located in Hipparchus library
(in
the "computeDerivative" method of the class "ExpandableODE") and not in
Orekit but I would like to know if someone already used this functionality
with the last versions of Orekit/Hipparchus ?
Thank you.
Christophe