use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.
the class BoxAndSolarArraySpacecraftTest method testNormalSunAlignedDS.
@Test
public void testNormalSunAlignedDS() throws OrekitException {
BoxAndSolarArraySpacecraft s = new BoxAndSolarArraySpacecraft(0, 0, 0, (date, frame) -> new TimeStampedPVCoordinates(date, new Vector3D(0, 1e6, 0), Vector3D.ZERO), 20.0, Vector3D.PLUS_J, 0.0, 1.0, 0.0);
DSFactory factory = new DSFactory(1, 2);
FieldVector3D<DerivativeStructure> normal = s.getNormal(AbsoluteDate.J2000_EPOCH, FramesFactory.getEME2000(), FieldVector3D.getZero(factory.getDerivativeField()), FieldRotation.getIdentity(factory.getDerivativeField()));
Assert.assertEquals(0, FieldVector3D.dotProduct(normal, Vector3D.PLUS_J).getReal(), 1.0e-16);
}
use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.
the class DragForceTest method accelerationDerivatives.
@Override
protected FieldVector3D<DerivativeStructure> accelerationDerivatives(final ForceModel forceModel, final AbsoluteDate date, final Frame frame, final FieldVector3D<DerivativeStructure> position, final FieldVector3D<DerivativeStructure> velocity, final FieldRotation<DerivativeStructure> rotation, final DerivativeStructure mass) throws OrekitException {
try {
java.lang.reflect.Field atmosphereField = DragForce.class.getDeclaredField("atmosphere");
atmosphereField.setAccessible(true);
Atmosphere atmosphere = (Atmosphere) atmosphereField.get(forceModel);
java.lang.reflect.Field spacecraftField = DragForce.class.getDeclaredField("spacecraft");
spacecraftField.setAccessible(true);
DragSensitive spacecraft = (DragSensitive) spacecraftField.get(forceModel);
// retrieve derivation properties
final DSFactory factory = mass.getFactory();
// get atmosphere properties in atmosphere own frame
final Frame atmFrame = atmosphere.getFrame();
final Transform toBody = frame.getTransformTo(atmFrame, date);
final FieldVector3D<DerivativeStructure> posBodyDS = toBody.transformPosition(position);
final Vector3D posBody = posBodyDS.toVector3D();
final Vector3D vAtmBody = atmosphere.getVelocity(date, posBody, atmFrame);
// to the Atmosphere interface
if (factory.getCompiler().getOrder() > 1) {
throw new OrekitException(OrekitMessages.OUT_OF_RANGE_DERIVATION_ORDER, factory.getCompiler().getOrder());
}
final double delta = 1.0;
final double x = posBody.getX();
final double y = posBody.getY();
final double z = posBody.getZ();
final double rho0 = atmosphere.getDensity(date, posBody, atmFrame);
final double dRhodX = (atmosphere.getDensity(date, new Vector3D(x + delta, y, z), atmFrame) - rho0) / delta;
final double dRhodY = (atmosphere.getDensity(date, new Vector3D(x, y + delta, z), atmFrame) - rho0) / delta;
final double dRhodZ = (atmosphere.getDensity(date, new Vector3D(x, y, z + delta), atmFrame) - rho0) / delta;
final double[] dXdQ = posBodyDS.getX().getAllDerivatives();
final double[] dYdQ = posBodyDS.getY().getAllDerivatives();
final double[] dZdQ = posBodyDS.getZ().getAllDerivatives();
final double[] rhoAll = new double[dXdQ.length];
rhoAll[0] = rho0;
for (int i = 1; i < rhoAll.length; ++i) {
rhoAll[i] = dRhodX * dXdQ[i] + dRhodY * dYdQ[i] + dRhodZ * dZdQ[i];
}
final DerivativeStructure rho = factory.build(rhoAll);
// we consider that at first order the atmosphere velocity in atmosphere frame
// does not depend on local position; however atmosphere velocity in inertial
// frame DOES depend on position since the transform between the frames depends
// on it, due to central body rotation rate and velocity composition.
// So we use the transform to get the correct partial derivatives on vAtm
final FieldVector3D<DerivativeStructure> vAtmBodyDS = new FieldVector3D<>(factory.constant(vAtmBody.getX()), factory.constant(vAtmBody.getY()), factory.constant(vAtmBody.getZ()));
final FieldPVCoordinates<DerivativeStructure> pvAtmBody = new FieldPVCoordinates<>(posBodyDS, vAtmBodyDS);
final FieldPVCoordinates<DerivativeStructure> pvAtm = toBody.getInverse().transformPVCoordinates(pvAtmBody);
// now we can compute relative velocity, it takes into account partial derivatives with respect to position
final FieldVector3D<DerivativeStructure> relativeVelocity = pvAtm.getVelocity().subtract(velocity);
// compute acceleration with all its partial derivatives
return spacecraft.dragAcceleration(new FieldAbsoluteDate<>(factory.getDerivativeField(), date), frame, position, rotation, mass, rho, relativeVelocity, forceModel.getParameters(factory.getDerivativeField()));
} catch (IllegalArgumentException | IllegalAccessException | NoSuchFieldException | SecurityException e) {
return null;
}
}
use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.
the class JB2008Test method testDensityGradient.
@Test
public void testDensityGradient() throws OrekitException {
final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, itrf);
final JB2008 atm = new JB2008(new InputParams(), CelestialBodyFactory.getSun(), earth);
final AbsoluteDate date = InputParams.TC[6];
// Build the position
final double alt = 400.;
final double lat = 60.;
final double lon = -70.;
final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(lat), FastMath.toRadians(lon), alt * 1000.);
final Vector3D pos = earth.transform(point);
// Run
DerivativeStructure zero = new DSFactory(1, 1).variable(0, 0.0);
FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(5, 10.0);
DerivativeStructure rhoX = differentiator.differentiate((double x) -> {
try {
return atm.getDensity(date, new Vector3D(1, pos, x, Vector3D.PLUS_I), itrf);
} catch (OrekitException oe) {
return Double.NaN;
}
}).value(zero);
DerivativeStructure rhoY = differentiator.differentiate((double y) -> {
try {
return atm.getDensity(date, new Vector3D(1, pos, y, Vector3D.PLUS_J), itrf);
} catch (OrekitException oe) {
return Double.NaN;
}
}).value(zero);
DerivativeStructure rhoZ = differentiator.differentiate((double z) -> {
try {
return atm.getDensity(date, new Vector3D(1, pos, z, Vector3D.PLUS_K), itrf);
} catch (OrekitException oe) {
return Double.NaN;
}
}).value(zero);
DSFactory factory3 = new DSFactory(3, 1);
Field<DerivativeStructure> field = factory3.getDerivativeField();
final DerivativeStructure rhoDS = atm.getDensity(new FieldAbsoluteDate<>(field, date), new FieldVector3D<>(factory3.variable(0, pos.getX()), factory3.variable(1, pos.getY()), factory3.variable(2, pos.getZ())), itrf);
Assert.assertEquals(rhoX.getValue(), rhoDS.getReal(), rhoX.getValue() * 2.0e-14);
Assert.assertEquals(rhoY.getValue(), rhoDS.getReal(), rhoY.getValue() * 2.0e-14);
Assert.assertEquals(rhoZ.getValue(), rhoDS.getReal(), rhoZ.getValue() * 2.0e-14);
Assert.assertEquals(rhoX.getPartialDerivative(1), rhoDS.getPartialDerivative(1, 0, 0), FastMath.abs(6.0e-10 * rhoX.getPartialDerivative(1)));
Assert.assertEquals(rhoY.getPartialDerivative(1), rhoDS.getPartialDerivative(0, 1, 0), FastMath.abs(6.0e-10 * rhoY.getPartialDerivative(1)));
Assert.assertEquals(rhoZ.getPartialDerivative(1), rhoDS.getPartialDerivative(0, 0, 1), FastMath.abs(6.0e-10 * rhoY.getPartialDerivative(1)));
}
use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.
the class JB2008Test method testDensityField.
@Test
public void testDensityField() throws OrekitException {
final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING, itrf);
final JB2008 atm = new JB2008(new InputParams(), CelestialBodyFactory.getSun(), earth);
final AbsoluteDate date = InputParams.TC[4];
for (double alt = 100; alt < 1000; alt += 50) {
for (double lat = -1.2; lat < 1.2; lat += 0.4) {
for (double lon = 0; lon < 6.28; lon += 0.8) {
final GeodeticPoint point = new GeodeticPoint(lat, lon, alt * 1000.);
final Vector3D pos = earth.transform(point);
Field<Decimal64> field = Decimal64Field.getInstance();
// Run
final double rho = atm.getDensity(date, pos, itrf);
final Decimal64 rho64 = atm.getDensity(new FieldAbsoluteDate<>(field, date), new FieldVector3D<>(field.getOne(), pos), itrf);
Assert.assertEquals(rho, rho64.getReal(), rho * 4.0e-13);
}
}
}
}
use of org.hipparchus.geometry.euclidean.threed.FieldVector3D in project Orekit by CS-SI.
the class AbstractForceModelTest method checkParameterDerivative.
protected void checkParameterDerivative(SpacecraftState state, ForceModel forceModel, String name, double hFactor, double tol) throws OrekitException {
final DSFactory factory11 = new DSFactory(1, 1);
final Field<DerivativeStructure> field = factory11.getDerivativeField();
final FieldSpacecraftState<DerivativeStructure> stateF = new FieldSpacecraftState<DerivativeStructure>(field, state);
final ParameterDriver[] drivers = forceModel.getParametersDrivers();
final DerivativeStructure[] parametersDS = new DerivativeStructure[drivers.length];
for (int i = 0; i < parametersDS.length; ++i) {
if (drivers[i].getName().equals(name)) {
parametersDS[i] = factory11.variable(0, drivers[i].getValue());
} else {
parametersDS[i] = factory11.constant(drivers[i].getValue());
}
}
FieldVector3D<DerivativeStructure> accDer = forceModel.acceleration(stateF, parametersDS);
Vector3D derivative = new Vector3D(accDer.getX().getPartialDerivative(1), accDer.getY().getPartialDerivative(1), accDer.getZ().getPartialDerivative(1));
int selected = -1;
final double[] parameters = new double[drivers.length];
for (int i = 0; i < drivers.length; ++i) {
parameters[i] = drivers[i].getValue();
if (drivers[i].getName().equals(name)) {
selected = i;
}
}
double p0 = parameters[selected];
double hParam = hFactor * p0;
drivers[selected].setValue(p0 - 1 * hParam);
parameters[selected] = drivers[selected].getValue();
Assert.assertEquals(p0 - 1 * hParam, parameters[selected], 1.0e-10);
final Vector3D gammaM1h = forceModel.acceleration(state, parameters);
drivers[selected].setValue(p0 + 1 * hParam);
parameters[selected] = drivers[selected].getValue();
Assert.assertEquals(p0 + 1 * hParam, parameters[selected], 1.0e-10);
final Vector3D gammaP1h = forceModel.acceleration(state, parameters);
drivers[selected].setValue(p0);
final Vector3D reference = new Vector3D(1 / (2 * hParam), gammaP1h.subtract(gammaM1h));
final Vector3D delta = derivative.subtract(reference);
Assert.assertEquals(0, delta.getNorm(), tol * reference.getNorm());
}
Aggregations