From dbd081ad1808d4e8550dd23971b60b862c7904e0 Mon Sep 17 00:00:00 2001 From: Tavian Barnes Date: Mon, 11 Oct 2010 19:01:41 -0400 Subject: Add a test for a system of vector equations. --- tests/EquationSystem-Vector.cpp | 70 +++++++++++++++++++++++++++++++++++++++++ tests/Makefile.am | 26 ++++++++------- 2 files changed, 84 insertions(+), 12 deletions(-) create mode 100644 tests/EquationSystem-Vector.cpp (limited to 'tests') diff --git a/tests/EquationSystem-Vector.cpp b/tests/EquationSystem-Vector.cpp new file mode 100644 index 0000000..b121026 --- /dev/null +++ b/tests/EquationSystem-Vector.cpp @@ -0,0 +1,70 @@ +/************************************************************************* + * Copyright (C) 2010 Tavian Barnes * + * * + * This file is part of The vZ Test Suite. * + * * + * The vZ Test Suite is free software; you can redistribute it and/or * + * modify it under the terms of the GNU General Public License as * + * published by the Free Software Foundation; either version 3 of the * + * License, or (at your option) any later version. * + * * + * The vZ Test Suite is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * + * General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see . * + *************************************************************************/ + +#include "vZ.hpp" +#include +#include +#include +#include + +typedef vZ::EquationSystem<2, vZ::Vector<2> > Y; + +// y'' = -(|y'|^2/|y|)*(y/|y|) (circular motion) +Y +f(double t, Y y) +{ + Y r; + r[0] = y[1]; + r[1] = -(dot(y[1], y[1])*y[0])/(dot(y[0], y[0])); + return r; +} + +int +main() +{ + Y y; + y[0] = vZ::Vector<2>(0.0, 1.0); + y[1] = vZ::Vector<2>(8*std::atan(1.0), 0.0); + vZ::GenericDP45Integrator integrator(f); + integrator.tol(1e-6) + .y(y) + .x(0.0) + .h(0.06); + + integrator.integrate(2.0); + + vZ::Vector<2> actual = integrator.y()[0]; + vZ::Vector<2> expected(0.0, 1.0); + + std::cout << std::setprecision(10) + << "Numerical: " << actual << std::endl + << "Expected: " << expected << std::endl + << "h: " << integrator.h() << std::endl + << "Iterations: " << integrator.iterations() << std::endl + << "Rejections: " << integrator.rejections() << std::endl; + + double error = norm(expected - actual)/norm(expected); + if (error > 4.0e-6) { + std::cerr << "Error: " << 100.0*error << "%" << std::endl; + return EXIT_FAILURE; + } else { + std::cout << "Error: " << 100.0*error << "%" << std::endl; + return EXIT_SUCCESS; + } +} diff --git a/tests/Makefile.am b/tests/Makefile.am index d072e35..2091b15 100644 --- a/tests/Makefile.am +++ b/tests/Makefile.am @@ -29,17 +29,19 @@ check_PROGRAMS = Euler-test \ CK45-test \ DP45-test \ Vector-test \ - EquationSystem-test + EquationSystem-test \ + EquationSystem-Vector-test TESTS = $(check_PROGRAMS) -Euler_test_SOURCES = Euler.cpp -Midpoint_test_SOURCES = Midpoint.cpp -Heun_test_SOURCES = Heun.cpp -RK4_test_SOURCES = RK4.cpp -HE12_test_SOURCES = HE12.cpp -BS23_test_SOURCES = BS23.cpp -RKF45_test_SOURCES = RKF45.cpp -CK45_test_SOURCES = CK45.cpp -DP45_test_SOURCES = DP45.cpp -Vector_test_SOURCES = Vector.cpp -EquationSystem_test_SOURCES = EquationSystem.cpp +Euler_test_SOURCES = Euler.cpp +Midpoint_test_SOURCES = Midpoint.cpp +Heun_test_SOURCES = Heun.cpp +RK4_test_SOURCES = RK4.cpp +HE12_test_SOURCES = HE12.cpp +BS23_test_SOURCES = BS23.cpp +RKF45_test_SOURCES = RKF45.cpp +CK45_test_SOURCES = CK45.cpp +DP45_test_SOURCES = DP45.cpp +Vector_test_SOURCES = Vector.cpp +EquationSystem_test_SOURCES = EquationSystem.cpp +EquationSystem_Vector_test_SOURCES = EquationSystem-Vector.cpp -- cgit v1.2.3