Apache
/
commons-math
mirror of https://github.com/apache/commons-math.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Fixed Butcher arrays constants for Luther integrator.

This commit is contained in:
Luc Maisonobe 2016-01-06 12:41:02 +01:00
parent 85b941ea1f
commit ba1ce9b9aa
2 changed files with 19 additions and 21 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
src/main/java/org/apache/commons/math4/ode/nonstiff/LutherFieldIntegrator.java
View File

@ -60,9 +60,6 @@ import org.apache.commons.math4.util.MathArrays;
public class LutherFieldIntegrator<T extends RealFieldElement<T>> public class LutherFieldIntegrator<T extends RealFieldElement<T>>
extends RungeKuttaFieldIntegrator<T> { extends RungeKuttaFieldIntegrator<T> {
/** Square root. */
private final T q;
/** Simple constructor. /** Simple constructor.
* Build a fourth-order Luther integrator with the given step. * Build a fourth-order Luther integrator with the given step.
* @param field field to which the time and state vector elements belong * @param field field to which the time and state vector elements belong
@ -70,12 +67,12 @@ public class LutherFieldIntegrator<T extends RealFieldElement<T>>
*/ */
public LutherFieldIntegrator(final Field<T> field, final T step) { public LutherFieldIntegrator(final Field<T> field, final T step) {
super(field, "Luther", step); super(field, "Luther", step);
q = getField().getOne().multiply(21).sqrt();
} }
/** {@inheritDoc} */ /** {@inheritDoc} */
@Override @Override
protected T[] getC() { protected T[] getC() {
final T q = getField().getZero().add(21).sqrt();
final T[] c = MathArrays.buildArray(getField(), 6); final T[] c = MathArrays.buildArray(getField(), 6);
c[0] = getField().getOne(); c[0] = getField().getOne();
c[1] = fraction(1, 2); c[1] = fraction(1, 2);
@ -89,6 +86,7 @@ public class LutherFieldIntegrator<T extends RealFieldElement<T>>
/** {@inheritDoc} */ /** {@inheritDoc} */
@Override @Override
protected T[][] getA() { protected T[][] getA() {
final T q = getField().getZero().add(21).sqrt();
final T[][] a = MathArrays.buildArray(getField(), 6, -1); final T[][] a = MathArrays.buildArray(getField(), 6, -1);
for (int i = 0; i < a.length; ++i) { for (int i = 0; i < a.length; ++i) {
a[i] = MathArrays.buildArray(getField(), i + 1); a[i] = MathArrays.buildArray(getField(), i + 1);
@ -98,22 +96,22 @@ public class LutherFieldIntegrator<T extends RealFieldElement<T>>
a[1][1] = fraction(1, 8); a[1][1] = fraction(1, 8);
a[2][0] = fraction(8, 27); a[2][0] = fraction(8, 27);
a[2][1] = fraction(2, 27); a[2][1] = fraction(2, 27);
a[2][2] = fraction(8, 27); a[2][2] = a[2][0];
a[3][0] = q.multiply( 9).subtract( 21).divide( 392); a[3][0] = q.multiply( 9).add( -21).divide( 392);
a[3][1] = q.multiply( 8).subtract( 56).divide( 392); a[3][1] = q.multiply( 8).add( -56).divide( 392);
a[3][2] = q.multiply( 48).subtract(336).divide(-392); a[3][2] = q.multiply( -48).add( 336).divide( 392);
a[3][3] = q.multiply( 3).subtract( 63).divide( 392); a[3][3] = q.multiply( 3).add( -63).divide( 392);
a[4][0] = q.multiply(255).add(1155).divide(-1960); a[4][0] = q.multiply(-255).add(-1155).divide(1960);
a[4][1] = q.multiply( 40).add( 280).divide(-1960); a[4][1] = q.multiply( -40).add( -280).divide(1960);
a[4][2] = q.multiply(20) .divide(-1960); a[4][2] = q.multiply(-320) .divide(1960);
a[4][3] = q.multiply(363).add( 63).divide( 1960); a[4][3] = q.multiply( 363).add( 63).divide(1960);
a[4][4] = q.multiply(392).add(2352).divide( 1960); a[4][4] = q.multiply( 392).add( 2352).divide(1960);
a[5][0] = q.multiply(105).add( 330).divide( 180); a[5][0] = q.multiply( 105).add( 330).divide( 180);
a[5][1] = fraction(120, 180); a[5][1] = fraction(2, 3);
a[5][2] = q.multiply(280).add(-200).divide( 180); a[5][2] = q.multiply( 280).add( -200).divide( 180);
a[5][3] = q.multiply(189).add(-126).divide(-180); a[5][3] = q.multiply(-189).add( 126).divide( 180);
a[5][4] = q.multiply(126).add( 686).divide(-180); a[5][4] = q.multiply(-126).add( -686).divide( 180);
a[5][5] = q.multiply( 70).add(-490).divide(-180); a[5][5] = q.multiply( -70).add( 490).divide( 180);
return a; return a;
} }

2
src/main/java/org/apache/commons/math4/ode/nonstiff/RungeKuttaFieldIntegrator.java
View File

@ -93,7 +93,7 @@ public abstract class RungeKuttaFieldIntegrator<T extends RealFieldElement<T>>
* @return p/q computed in the instance field * @return p/q computed in the instance field
*/ */
protected T fraction(final int p, final int q) { protected T fraction(final int p, final int q) {
return getField().getOne().multiply(p).divide(q); return getField().getZero().add(p).divide(q);
} }
/** Get the time steps from Butcher array (without the first zero). /** Get the time steps from Butcher array (without the first zero).