001 /* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017 018 package org.apache.commons.math.ode; 019 020 import junit.framework.*; 021 import java.util.Random; 022 023 import org.apache.commons.math.ode.ContinuousOutputModel; 024 import org.apache.commons.math.ode.DerivativeException; 025 import org.apache.commons.math.ode.FirstOrderIntegrator; 026 import org.apache.commons.math.ode.IntegratorException; 027 import org.apache.commons.math.ode.nonstiff.DormandPrince54Integrator; 028 import org.apache.commons.math.ode.nonstiff.DormandPrince853Integrator; 029 import org.apache.commons.math.ode.sampling.DummyStepInterpolator; 030 import org.apache.commons.math.ode.sampling.StepInterpolator; 031 032 public class ContinuousOutputModelTest 033 extends TestCase { 034 035 public ContinuousOutputModelTest(String name) { 036 super(name); 037 pb = null; 038 integ = null; 039 } 040 041 public void testBoundaries() 042 throws DerivativeException, IntegratorException { 043 integ.addStepHandler(new ContinuousOutputModel()); 044 integ.integrate(pb, 045 pb.getInitialTime(), pb.getInitialState(), 046 pb.getFinalTime(), new double[pb.getDimension()]); 047 ContinuousOutputModel cm = (ContinuousOutputModel) integ.getStepHandlers().iterator().next(); 048 cm.setInterpolatedTime(2.0 * pb.getInitialTime() - pb.getFinalTime()); 049 cm.setInterpolatedTime(2.0 * pb.getFinalTime() - pb.getInitialTime()); 050 cm.setInterpolatedTime(0.5 * (pb.getFinalTime() + pb.getInitialTime())); 051 } 052 053 public void testRandomAccess() 054 throws DerivativeException, IntegratorException { 055 056 ContinuousOutputModel cm = new ContinuousOutputModel(); 057 integ.addStepHandler(cm); 058 integ.integrate(pb, 059 pb.getInitialTime(), pb.getInitialState(), 060 pb.getFinalTime(), new double[pb.getDimension()]); 061 062 Random random = new Random(347588535632l); 063 double maxError = 0.0; 064 for (int i = 0; i < 1000; ++i) { 065 double r = random.nextDouble(); 066 double time = r * pb.getInitialTime() + (1.0 - r) * pb.getFinalTime(); 067 cm.setInterpolatedTime(time); 068 double[] interpolatedY = cm.getInterpolatedState (); 069 double[] theoreticalY = pb.computeTheoreticalState(time); 070 double dx = interpolatedY[0] - theoreticalY[0]; 071 double dy = interpolatedY[1] - theoreticalY[1]; 072 double error = dx * dx + dy * dy; 073 if (error > maxError) { 074 maxError = error; 075 } 076 } 077 078 assertTrue(maxError < 1.0e-9); 079 080 } 081 082 public void testModelsMerging() 083 throws DerivativeException, IntegratorException { 084 085 // theoretical solution: y[0] = cos(t), y[1] = sin(t) 086 FirstOrderDifferentialEquations problem = 087 new FirstOrderDifferentialEquations() { 088 private static final long serialVersionUID = 2472449657345878299L; 089 public void computeDerivatives(double t, double[] y, double[] dot) 090 throws DerivativeException { 091 dot[0] = -y[1]; 092 dot[1] = y[0]; 093 } 094 public int getDimension() { 095 return 2; 096 } 097 }; 098 099 // integrate backward from π to 0; 100 ContinuousOutputModel cm1 = new ContinuousOutputModel(); 101 FirstOrderIntegrator integ1 = 102 new DormandPrince853Integrator(0, 1.0, 1.0e-8, 1.0e-8); 103 integ1.addStepHandler(cm1); 104 integ1.integrate(problem, Math.PI, new double[] { -1.0, 0.0 }, 105 0, new double[2]); 106 107 // integrate backward from 2π to π 108 ContinuousOutputModel cm2 = new ContinuousOutputModel(); 109 FirstOrderIntegrator integ2 = 110 new DormandPrince853Integrator(0, 0.1, 1.0e-12, 1.0e-12); 111 integ2.addStepHandler(cm2); 112 integ2.integrate(problem, 2.0 * Math.PI, new double[] { 1.0, 0.0 }, 113 Math.PI, new double[2]); 114 115 // merge the two half circles 116 ContinuousOutputModel cm = new ContinuousOutputModel(); 117 cm.append(cm2); 118 cm.append(new ContinuousOutputModel()); 119 cm.append(cm1); 120 121 // check circle 122 assertEquals(2.0 * Math.PI, cm.getInitialTime(), 1.0e-12); 123 assertEquals(0, cm.getFinalTime(), 1.0e-12); 124 assertEquals(cm.getFinalTime(), cm.getInterpolatedTime(), 1.0e-12); 125 for (double t = 0; t < 2.0 * Math.PI; t += 0.1) { 126 cm.setInterpolatedTime(t); 127 double[] y = cm.getInterpolatedState(); 128 assertEquals(Math.cos(t), y[0], 1.0e-7); 129 assertEquals(Math.sin(t), y[1], 1.0e-7); 130 } 131 132 } 133 134 public void testErrorConditions() 135 throws DerivativeException { 136 137 ContinuousOutputModel cm = new ContinuousOutputModel(); 138 cm.handleStep(buildInterpolator(0, new double[] { 0.0, 1.0, -2.0 }, 1), true); 139 140 // dimension mismatch 141 assertTrue(checkAppendError(cm, 1.0, new double[] { 0.0, 1.0 }, 2.0)); 142 143 // hole between time ranges 144 assertTrue(checkAppendError(cm, 10.0, new double[] { 0.0, 1.0, -2.0 }, 20.0)); 145 146 // propagation direction mismatch 147 assertTrue(checkAppendError(cm, 1.0, new double[] { 0.0, 1.0, -2.0 }, 0.0)); 148 149 // no errors 150 assertFalse(checkAppendError(cm, 1.0, new double[] { 0.0, 1.0, -2.0 }, 2.0)); 151 152 } 153 154 private boolean checkAppendError(ContinuousOutputModel cm, 155 double t0, double[] y0, double t1) 156 throws DerivativeException { 157 try { 158 ContinuousOutputModel otherCm = new ContinuousOutputModel(); 159 otherCm.handleStep(buildInterpolator(t0, y0, t1), true); 160 cm.append(otherCm); 161 } catch(IllegalArgumentException iae) { 162 //expected behavior 163 return true; 164 } 165 return false; 166 } 167 168 private StepInterpolator buildInterpolator(double t0, double[] y0, double t1) { 169 DummyStepInterpolator interpolator = new DummyStepInterpolator(y0, t1 >= t0); 170 interpolator.storeTime(t0); 171 interpolator.shift(); 172 interpolator.storeTime(t1); 173 return interpolator; 174 } 175 176 public void checkValue(double value, double reference) { 177 assertTrue(Math.abs(value - reference) < 1.0e-10); 178 } 179 180 public static Test suite() { 181 return new TestSuite(ContinuousOutputModelTest.class); 182 } 183 184 @Override 185 public void setUp() { 186 pb = new TestProblem3(0.9); 187 double minStep = 0; 188 double maxStep = pb.getFinalTime() - pb.getInitialTime(); 189 integ = new DormandPrince54Integrator(minStep, maxStep, 1.0e-8, 1.0e-8); 190 } 191 192 @Override 193 public void tearDown() { 194 pb = null; 195 integ = null; 196 } 197 198 TestProblem3 pb; 199 FirstOrderIntegrator integ; 200 201 }