use of org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction in project Orekit by CS-SI.
the class PoissonSeriesParserTest method testDerivativesFromDoubleAPI.
@Test
public void testDerivativesFromDoubleAPI() throws OrekitException {
Utils.setDataRoot("regular-data");
String directory = "/assets/org/orekit/IERS-conventions/";
PoissonSeriesParser parser = new PoissonSeriesParser(17).withPolynomialPart('t', PolynomialParser.Unit.NO_UNITS).withFirstDelaunay(4).withFirstPlanetary(9).withSinCos(0, 2, 1.0, 3, 1.0);
InputStream xStream = getClass().getResourceAsStream(directory + "2010/tab5.2a.txt");
PoissonSeries xSeries = parser.parse(xStream, "2010/tab5.2a.txt");
InputStream yStream = getClass().getResourceAsStream(directory + "2010/tab5.2b.txt");
PoissonSeries ySeries = parser.parse(yStream, "2010/tab5.2b.txt");
InputStream zStream = getClass().getResourceAsStream(directory + "2010/tab5.2d.txt");
PoissonSeries zSeries = parser.parse(zStream, "2010/tab5.2d.txt");
final PoissonSeries.CompiledSeries compiled = PoissonSeries.compile(xSeries, ySeries, zSeries);
TimeScale ut1 = TimeScalesFactory.getUT1(FramesFactory.getEOPHistory(IERSConventions.IERS_2010, true));
final FundamentalNutationArguments arguments = IERSConventions.IERS_2010.getNutationArguments(ut1);
UnivariateDifferentiableVectorFunction finite = new FiniteDifferencesDifferentiator(4, 0.4).differentiate((double t) -> compiled.value(arguments.evaluateAll(AbsoluteDate.J2000_EPOCH.shiftedBy(t))));
DSFactory factory = new DSFactory(1, 1);
for (double t = 0; t < Constants.JULIAN_DAY; t += 120) {
// computation of derivatives from API
double[] dAPI = compiled.derivative(arguments.evaluateAll(AbsoluteDate.J2000_EPOCH.shiftedBy(t)));
// finite differences computation of derivatives
DerivativeStructure[] d = finite.value(factory.variable(0, t));
Assert.assertEquals(d.length, dAPI.length);
for (int i = 0; i < d.length; ++i) {
Assert.assertEquals(d[i].getPartialDerivative(1), dAPI[i], FastMath.abs(2.0e-7 * d[i].getPartialDerivative(1)));
}
}
}
use of org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction in project Orekit by CS-SI.
the class IERSConventionsTest method checkVectorFunctionConsistency.
private void checkVectorFunctionConsistency(final TimeVectorFunction function, final int dimension, final AbsoluteDate date, final double span, final double sampleStep, final double h, final double valueTolerance, final double derivativeTolerance) {
UnivariateDifferentiableVectorFunction differentiated = new FiniteDifferencesDifferentiator(4, h).differentiate(new UnivariateVectorFunction() {
@Override
public double[] value(final double dt) {
return function.value(date.shiftedBy(dt));
}
});
DSFactory factory = new DSFactory(1, 1);
FieldAbsoluteDate<DerivativeStructure> dsDate = new FieldAbsoluteDate<>(date, factory.constant(0.0));
double maxValueError = 0;
double maxDerivativeError = 0;
for (double dt = 0; dt < span; dt += sampleStep) {
DerivativeStructure dsdt = factory.variable(0, dt);
DerivativeStructure[] yRef = differentiated.value(dsdt);
DerivativeStructure[] y = function.value(dsDate.shiftedBy(dsdt));
Assert.assertEquals(dimension, yRef.length);
Assert.assertEquals(dimension, y.length);
for (int i = 0; i < dimension; ++i) {
maxValueError = FastMath.max(maxValueError, FastMath.abs(yRef[i].getValue() - y[i].getValue()));
maxDerivativeError = FastMath.max(maxDerivativeError, FastMath.abs(yRef[i].getPartialDerivative(1) - y[i].getPartialDerivative(1)));
}
}
Assert.assertEquals(0, maxValueError, valueTolerance);
Assert.assertEquals(0, maxDerivativeError, derivativeTolerance);
}
use of org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction in project Orekit by CS-SI.
the class GroundStationTest method doTestCartesianDerivatives.
private void doTestCartesianDerivatives(double latitude, double longitude, double altitude, double stepFactor, double relativeTolerancePositionValue, double relativeTolerancePositionDerivative, double relativeToleranceVelocityValue, double relativeToleranceVelocityDerivative, String... parameterPattern) throws OrekitException {
Utils.setDataRoot("regular-data");
final Frame eme2000 = FramesFactory.getEME2000();
final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
final AbsoluteDate date0 = date.shiftedBy(50000);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final GroundStation station = new GroundStation(new TopocentricFrame(earth, new GeodeticPoint(latitude, longitude, altitude), "dummy"));
final DSFactory factory = new DSFactory(parameterPattern.length, 1);
final FieldAbsoluteDate<DerivativeStructure> dateDS = new FieldAbsoluteDate<>(factory.getDerivativeField(), date);
ParameterDriver[] selectedDrivers = new ParameterDriver[parameterPattern.length];
UnivariateDifferentiableVectorFunction[] dFCartesian = new UnivariateDifferentiableVectorFunction[parameterPattern.length];
final ParameterDriver[] allDrivers = selectAllDrivers(station);
for (ParameterDriver driver : allDrivers) {
driver.setReferenceDate(date0);
}
Map<String, Integer> indices = new HashMap<>();
for (int k = 0; k < dFCartesian.length; ++k) {
for (int i = 0; i < allDrivers.length; ++i) {
if (allDrivers[i].getName().matches(parameterPattern[k])) {
selectedDrivers[k] = allDrivers[i];
dFCartesian[k] = differentiatedStationPV(station, eme2000, date, selectedDrivers[k], stepFactor);
indices.put(selectedDrivers[k].getName(), k);
}
}
}
;
DSFactory factory11 = new DSFactory(1, 1);
RandomGenerator generator = new Well19937a(0x084d58a19c498a54l);
double maxPositionValueRelativeError = 0;
double maxPositionDerivativeRelativeError = 0;
double maxVelocityValueRelativeError = 0;
double maxVelocityDerivativeRelativeError = 0;
for (int i = 0; i < 1000; ++i) {
// randomly change one parameter
ParameterDriver changed = allDrivers[generator.nextInt(allDrivers.length)];
changed.setNormalizedValue(2 * generator.nextDouble() - 1);
// transform to check
FieldTransform<DerivativeStructure> t = station.getOffsetToInertial(eme2000, dateDS, factory, indices);
FieldPVCoordinates<DerivativeStructure> pv = t.transformPVCoordinates(FieldPVCoordinates.getZero(factory.getDerivativeField()));
for (int k = 0; k < dFCartesian.length; ++k) {
// reference values and derivatives computed using finite differences
DerivativeStructure[] refCartesian = dFCartesian[k].value(factory11.variable(0, selectedDrivers[k].getValue()));
// position
final Vector3D refP = new Vector3D(refCartesian[0].getValue(), refCartesian[1].getValue(), refCartesian[2].getValue());
final Vector3D resP = new Vector3D(pv.getPosition().getX().getValue(), pv.getPosition().getY().getValue(), pv.getPosition().getZ().getValue());
maxPositionValueRelativeError = FastMath.max(maxPositionValueRelativeError, Vector3D.distance(refP, resP) / refP.getNorm());
final Vector3D refPD = new Vector3D(refCartesian[0].getPartialDerivative(1), refCartesian[1].getPartialDerivative(1), refCartesian[2].getPartialDerivative(1));
final Vector3D resPD = new Vector3D(pv.getPosition().getX().getAllDerivatives()[k + 1], pv.getPosition().getY().getAllDerivatives()[k + 1], pv.getPosition().getZ().getAllDerivatives()[k + 1]);
maxPositionDerivativeRelativeError = FastMath.max(maxPositionDerivativeRelativeError, Vector3D.distance(refPD, resPD) / refPD.getNorm());
// velocity
final Vector3D refV = new Vector3D(refCartesian[3].getValue(), refCartesian[4].getValue(), refCartesian[5].getValue());
final Vector3D resV = new Vector3D(pv.getVelocity().getX().getValue(), pv.getVelocity().getY().getValue(), pv.getVelocity().getZ().getValue());
maxVelocityValueRelativeError = FastMath.max(maxVelocityValueRelativeError, Vector3D.distance(refV, resV) / refV.getNorm());
final Vector3D refVD = new Vector3D(refCartesian[3].getPartialDerivative(1), refCartesian[4].getPartialDerivative(1), refCartesian[5].getPartialDerivative(1));
final Vector3D resVD = new Vector3D(pv.getVelocity().getX().getAllDerivatives()[k + 1], pv.getVelocity().getY().getAllDerivatives()[k + 1], pv.getVelocity().getZ().getAllDerivatives()[k + 1]);
maxVelocityDerivativeRelativeError = FastMath.max(maxVelocityDerivativeRelativeError, Vector3D.distance(refVD, resVD) / refVD.getNorm());
}
}
Assert.assertEquals(0.0, maxPositionValueRelativeError, relativeTolerancePositionValue);
Assert.assertEquals(0.0, maxPositionDerivativeRelativeError, relativeTolerancePositionDerivative);
Assert.assertEquals(0.0, maxVelocityValueRelativeError, relativeToleranceVelocityValue);
Assert.assertEquals(0.0, maxVelocityDerivativeRelativeError, relativeToleranceVelocityDerivative);
}
use of org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction in project Orekit by CS-SI.
the class GroundStationTest method doTestAngularDerivatives.
private void doTestAngularDerivatives(double latitude, double longitude, double altitude, double stepFactor, double toleranceRotationValue, double toleranceRotationDerivative, double toleranceRotationRateValue, double toleranceRotationRateDerivative, String... parameterPattern) throws OrekitException {
Utils.setDataRoot("regular-data");
final Frame eme2000 = FramesFactory.getEME2000();
final AbsoluteDate date = AbsoluteDate.J2000_EPOCH;
final AbsoluteDate date0 = date.shiftedBy(50000);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, FramesFactory.getITRF(IERSConventions.IERS_2010, true));
final GroundStation station = new GroundStation(new TopocentricFrame(earth, new GeodeticPoint(latitude, longitude, altitude), "dummy"));
final DSFactory factory = new DSFactory(parameterPattern.length, 1);
final FieldAbsoluteDate<DerivativeStructure> dateDS = new FieldAbsoluteDate<>(factory.getDerivativeField(), date);
ParameterDriver[] selectedDrivers = new ParameterDriver[parameterPattern.length];
UnivariateDifferentiableVectorFunction[] dFAngular = new UnivariateDifferentiableVectorFunction[parameterPattern.length];
final ParameterDriver[] allDrivers = selectAllDrivers(station);
for (ParameterDriver driver : allDrivers) {
driver.setReferenceDate(date0);
}
Map<String, Integer> indices = new HashMap<>();
for (int k = 0; k < dFAngular.length; ++k) {
for (int i = 0; i < allDrivers.length; ++i) {
if (allDrivers[i].getName().matches(parameterPattern[k])) {
selectedDrivers[k] = allDrivers[i];
dFAngular[k] = differentiatedTransformAngular(station, eme2000, date, selectedDrivers[k], stepFactor);
indices.put(selectedDrivers[k].getName(), k);
}
}
}
;
DSFactory factory11 = new DSFactory(1, 1);
RandomGenerator generator = new Well19937a(0xa01a1d8fe5d80af7l);
double maxRotationValueError = 0;
double maxRotationDerivativeError = 0;
double maxRotationRateValueError = 0;
double maxRotationRateDerivativeError = 0;
for (int i = 0; i < 1000; ++i) {
// randomly change one parameter
ParameterDriver changed = allDrivers[generator.nextInt(allDrivers.length)];
changed.setNormalizedValue(2 * generator.nextDouble() - 1);
// transform to check
FieldTransform<DerivativeStructure> t = station.getOffsetToInertial(eme2000, dateDS, factory, indices);
for (int k = 0; k < dFAngular.length; ++k) {
// reference values and derivatives computed using finite differences
DerivativeStructure[] refAngular = dFAngular[k].value(factory11.variable(0, selectedDrivers[k].getValue()));
// rotation
final Rotation refQ = new Rotation(refAngular[0].getValue(), refAngular[1].getValue(), refAngular[2].getValue(), refAngular[3].getValue(), true);
final Rotation resQ = t.getRotation().toRotation();
maxRotationValueError = FastMath.max(maxRotationValueError, Rotation.distance(refQ, resQ));
double sign = FastMath.copySign(1.0, refAngular[0].getValue() * t.getRotation().getQ0().getValue() + refAngular[1].getValue() * t.getRotation().getQ1().getValue() + refAngular[2].getValue() * t.getRotation().getQ2().getValue() + refAngular[3].getValue() * t.getRotation().getQ3().getValue());
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[0].getPartialDerivative(1) - t.getRotation().getQ0().getAllDerivatives()[k + 1]));
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[1].getPartialDerivative(1) - t.getRotation().getQ1().getAllDerivatives()[k + 1]));
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[2].getPartialDerivative(1) - t.getRotation().getQ2().getAllDerivatives()[k + 1]));
maxRotationDerivativeError = FastMath.max(maxRotationDerivativeError, FastMath.abs(sign * refAngular[3].getPartialDerivative(1) - t.getRotation().getQ3().getAllDerivatives()[k + 1]));
// rotation rate
final Vector3D refRate = new Vector3D(refAngular[4].getValue(), refAngular[5].getValue(), refAngular[6].getValue());
final Vector3D resRate = t.getRotationRate().toVector3D();
final Vector3D refRateD = new Vector3D(refAngular[4].getPartialDerivative(1), refAngular[5].getPartialDerivative(1), refAngular[6].getPartialDerivative(1));
final Vector3D resRateD = new Vector3D(t.getRotationRate().getX().getAllDerivatives()[k + 1], t.getRotationRate().getY().getAllDerivatives()[k + 1], t.getRotationRate().getZ().getAllDerivatives()[k + 1]);
maxRotationRateValueError = FastMath.max(maxRotationRateValueError, Vector3D.distance(refRate, resRate));
maxRotationRateDerivativeError = FastMath.max(maxRotationRateDerivativeError, Vector3D.distance(refRateD, resRateD));
}
}
Assert.assertEquals(0.0, maxRotationValueError, toleranceRotationValue);
Assert.assertEquals(0.0, maxRotationDerivativeError, toleranceRotationDerivative);
Assert.assertEquals(0.0, maxRotationRateValueError, toleranceRotationRateValue);
Assert.assertEquals(0.0, maxRotationRateDerivativeError, toleranceRotationRateDerivative);
}
use of org.hipparchus.analysis.differentiation.UnivariateDifferentiableVectorFunction in project Orekit by CS-SI.
the class PredefinedIAUPolesTest method testDerivatives.
@Test
public void testDerivatives() {
final DSFactory factory = new DSFactory(1, 1);
final AbsoluteDate ref = AbsoluteDate.J2000_EPOCH;
final FieldAbsoluteDate<DerivativeStructure> refDS = new FieldAbsoluteDate<>(factory.getDerivativeField(), ref);
FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(8, 60.0);
for (final IAUPole iaupole : PredefinedIAUPoles.values()) {
UnivariateDifferentiableVectorFunction dPole = differentiator.differentiate(new UnivariateVectorFunction() {
@Override
public double[] value(double t) {
return iaupole.getPole(ref.shiftedBy(t)).toArray();
}
});
UnivariateDifferentiableFunction dMeridian = differentiator.differentiate(new UnivariateFunction() {
@Override
public double value(double t) {
return iaupole.getPrimeMeridianAngle(ref.shiftedBy(t));
}
});
for (double dt = 0; dt < Constants.JULIAN_YEAR; dt += 3600) {
final DerivativeStructure dtDS = factory.variable(0, dt);
final DerivativeStructure[] refPole = dPole.value(dtDS);
final DerivativeStructure[] fieldPole = iaupole.getPole(refDS.shiftedBy(dtDS)).toArray();
for (int i = 0; i < 3; ++i) {
Assert.assertEquals(refPole[i].getValue(), fieldPole[i].getValue(), 2.0e-15);
Assert.assertEquals(refPole[i].getPartialDerivative(1), fieldPole[i].getPartialDerivative(1), 4.0e-17);
}
final DerivativeStructure refMeridian = dMeridian.value(dtDS);
final DerivativeStructure fieldMeridian = iaupole.getPrimeMeridianAngle(refDS.shiftedBy(dtDS));
Assert.assertEquals(refMeridian.getValue(), fieldMeridian.getValue(), 4.0e-12);
Assert.assertEquals(refMeridian.getPartialDerivative(1), fieldMeridian.getPartialDerivative(1), 9.0e-14);
}
}
}
Aggregations