~ubuntu-branches/ubuntu/maverick/commons-math/maverick

« back to all changes in this revision

Viewing changes to src/main/java/org/apache/commons/math/ode/nonstiff/GillStepInterpolator.java

  • Committer: Bazaar Package Importer
  • Author(s): Damien Raude-Morvan
  • Date: 2009-08-22 01:13:25 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20090822011325-hi4peq1ua5weguwn
Tags: 2.0-1
* New upstream release.
* Set Maintainer field to Debian Java Team
* Add myself as Uploaders
* Switch to Quilt patch system:
  - Refresh all patchs
  - Remove B-D on dpatch, Add B-D on quilt
  - Include patchsys-quilt.mk in debian/rules
* Bump Standards-Version to 3.8.3:
  - Add a README.source to describe patch system
* Maven POMs:
  - Add a Build-Depends-Indep dependency on maven-repo-helper
  - Use mh_installpom and mh_installjar to install the POM and the jar to the
    Maven repository
* Use default-jdk/jre:
  - Depends on java5-runtime-headless
  - Build-Depends on default-jdk
  - Use /usr/lib/jvm/default-java as JAVA_HOME
* Move api documentation to /usr/share/doc/libcommons-math-java/api
* Build-Depends on junit4 instead of junit

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*
 
2
 * Licensed to the Apache Software Foundation (ASF) under one or more
 
3
 * contributor license agreements.  See the NOTICE file distributed with
 
4
 * this work for additional information regarding copyright ownership.
 
5
 * The ASF licenses this file to You under the Apache License, Version 2.0
 
6
 * (the "License"); you may not use this file except in compliance with
 
7
 * the License.  You may obtain a copy of the License at
 
8
 *
 
9
 *      http://www.apache.org/licenses/LICENSE-2.0
 
10
 *
 
11
 * Unless required by applicable law or agreed to in writing, software
 
12
 * distributed under the License is distributed on an "AS IS" BASIS,
 
13
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 
14
 * See the License for the specific language governing permissions and
 
15
 * limitations under the License.
 
16
 */
 
17
 
 
18
package org.apache.commons.math.ode.nonstiff;
 
19
 
 
20
import org.apache.commons.math.ode.DerivativeException;
 
21
import org.apache.commons.math.ode.sampling.AbstractStepInterpolator;
 
22
import org.apache.commons.math.ode.sampling.StepInterpolator;
 
23
 
 
24
/**
 
25
 * This class implements a step interpolator for the Gill fourth
 
26
 * order Runge-Kutta integrator.
 
27
 *
 
28
 * <p>This interpolator allows to compute dense output inside the last
 
29
 * step computed. The interpolation equation is consistent with the
 
30
 * integration scheme :
 
31
 *
 
32
 * <pre>
 
33
 *   y(t_n + theta h) = y (t_n + h)
 
34
 *                    - (1 - theta) (h/6) [ (1 - theta) (1 - 4 theta) y'_1
 
35
 *                                        + (1 - theta) (1 + 2 theta) ((2-q) y'_2 + (2+q) y'_3)
 
36
 *                                        + (1 + theta + 4 theta^2) y'_4
 
37
 *                                        ]
 
38
 * </pre>
 
39
 * where theta belongs to [0 ; 1], q = sqrt(2) and where y'_1 to y'_4
 
40
 * are the four evaluations of the derivatives already computed during
 
41
 * the step.</p>
 
42
 *
 
43
 * @see GillIntegrator
 
44
 * @version $Revision: 782432 $ $Date: 2009-06-07 15:08:26 -0400 (Sun, 07 Jun 2009) $
 
45
 * @since 1.2
 
46
 */
 
47
 
 
48
class GillStepInterpolator
 
49
  extends RungeKuttaStepInterpolator {
 
50
    
 
51
  /** Simple constructor.
 
52
   * This constructor builds an instance that is not usable yet, the
 
53
   * {@link AbstractStepInterpolator#reinitialize} method should be called
 
54
   * before using the instance in order to initialize the internal arrays. This
 
55
   * constructor is used only in order to delay the initialization in
 
56
   * some cases. The {@link RungeKuttaIntegrator} class uses the
 
57
   * prototyping design pattern to create the step interpolators by
 
58
   * cloning an uninitialized model and latter initializing the copy.
 
59
   */
 
60
  public GillStepInterpolator() {
 
61
  }
 
62
 
 
63
  /** Copy constructor.
 
64
   * @param interpolator interpolator to copy from. The copy is a deep
 
65
   * copy: its arrays are separated from the original arrays of the
 
66
   * instance
 
67
   */
 
68
  public GillStepInterpolator(final GillStepInterpolator interpolator) {
 
69
    super(interpolator);
 
70
  }
 
71
 
 
72
  /** {@inheritDoc} */
 
73
  @Override
 
74
  protected StepInterpolator doCopy() {
 
75
    return new GillStepInterpolator(this);
 
76
  }
 
77
 
 
78
 
 
79
  /** {@inheritDoc} */
 
80
  @Override
 
81
  protected void computeInterpolatedStateAndDerivatives(final double theta,
 
82
                                          final double oneMinusThetaH)
 
83
    throws DerivativeException {
 
84
 
 
85
    final double twoTheta  = 2 * theta;
 
86
    final double fourTheta = 4 * theta;
 
87
    final double s         = oneMinusThetaH / 6.0;
 
88
    final double oMt       = 1 - theta;
 
89
    final double soMt      = s * oMt;
 
90
    final double c23       = soMt * (1 + twoTheta);
 
91
    final double coeff1    = soMt * (1 - fourTheta);
 
92
    final double coeff2    = c23  * tMq;
 
93
    final double coeff3    = c23  * tPq;
 
94
    final double coeff4    = s * (1 + theta * (1 + fourTheta));
 
95
    final double coeffDot1 = theta * (twoTheta - 3) + 1;
 
96
    final double cDot23    = theta * oMt;
 
97
    final double coeffDot2 = cDot23  * tMq;
 
98
    final double coeffDot3 = cDot23  * tPq;
 
99
    final double coeffDot4 = theta * (twoTheta - 1);
 
100
 
 
101
    for (int i = 0; i < interpolatedState.length; ++i) {
 
102
        final double yDot1 = yDotK[0][i];
 
103
        final double yDot2 = yDotK[1][i];
 
104
        final double yDot3 = yDotK[2][i];
 
105
        final double yDot4 = yDotK[3][i];
 
106
        interpolatedState[i] =
 
107
            currentState[i] - coeff1 * yDot1 - coeff2 * yDot2 - coeff3 * yDot3 - coeff4 * yDot4;
 
108
        interpolatedDerivatives[i] =
 
109
            coeffDot1 * yDot1 + coeffDot2 * yDot2 + coeffDot3 * yDot3 + coeffDot4 * yDot4;
 
110
     }
 
111
 
 
112
  }
 
113
 
 
114
  /** First Gill coefficient. */
 
115
  private static final double tMq = 2 - Math.sqrt(2.0);
 
116
 
 
117
  /** Second Gill coefficient. */
 
118
  private static final double tPq = 2 + Math.sqrt(2.0);
 
119
 
 
120
  /** Serializable version identifier */
 
121
  private static final long serialVersionUID = -107804074496313322L;
 
122
 
 
123
}