test/root_finder.cc

Code
Comments
Other
Rev Date Author Line
4237 17 Sep 22 peter 1 // $Id$
4237 17 Sep 22 peter 2
4237 17 Sep 22 peter 3 /*
4237 17 Sep 22 peter 4   Copyright (C) 2022 Peter Johansson
4237 17 Sep 22 peter 5
4237 17 Sep 22 peter 6   This file is part of the yat library, https://dev.thep.lu.se/yat
4237 17 Sep 22 peter 7
4237 17 Sep 22 peter 8   The yat library is free software; you can redistribute it and/or
4237 17 Sep 22 peter 9   modify it under the terms of the GNU General Public License as
4237 17 Sep 22 peter 10   published by the Free Software Foundation; either version 3 of the
4237 17 Sep 22 peter 11   License, or (at your option) any later version.
4237 17 Sep 22 peter 12
4237 17 Sep 22 peter 13   The yat library is distributed in the hope that it will be useful,
4237 17 Sep 22 peter 14   but WITHOUT ANY WARRANTY; without even the implied warranty of
4237 17 Sep 22 peter 15   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
4237 17 Sep 22 peter 16   General Public License for more details.
4237 17 Sep 22 peter 17
4237 17 Sep 22 peter 18   You should have received a copy of the GNU General Public License
4237 17 Sep 22 peter 19   along with yat. If not, see <https://www.gnu.org/licenses/>.
4237 17 Sep 22 peter 20 */
4237 17 Sep 22 peter 21
4237 17 Sep 22 peter 22 #include <config.h>
4237 17 Sep 22 peter 23
4237 17 Sep 22 peter 24 #include "Suite.h"
4237 17 Sep 22 peter 25
4237 17 Sep 22 peter 26 #include "yat/utility/Bisection.h"
4237 17 Sep 22 peter 27 #include "yat/utility/Brent.h"
4237 17 Sep 22 peter 28 #include "yat/utility/FalsePosition.h"
4237 17 Sep 22 peter 29
4237 17 Sep 22 peter 30 #include <gsl/gsl_math.h>
4237 17 Sep 22 peter 31
4237 17 Sep 22 peter 32 #include <cmath>
4237 17 Sep 22 peter 33
4237 17 Sep 22 peter 34 using namespace theplu::yat;
4237 17 Sep 22 peter 35 using namespace theplu::yat::utility;
4237 17 Sep 22 peter 36
4237 17 Sep 22 peter 37 struct MyFunction
4237 17 Sep 22 peter 38 {
4237 17 Sep 22 peter 39   double operator()(double x) const
4237 17 Sep 22 peter 40   {
4237 17 Sep 22 peter 41     return std::exp(x) - 2;
4237 17 Sep 22 peter 42   }
4237 17 Sep 22 peter 43 };
4237 17 Sep 22 peter 44
4237 17 Sep 22 peter 45
4237 17 Sep 22 peter 46 void run_test(RootFinder& finder, test::Suite& suite)
4237 17 Sep 22 peter 47 {
4237 17 Sep 22 peter 48   MyFunction func;
4237 17 Sep 22 peter 49   suite.out() << "name: " << finder.name() << "\n";
4237 17 Sep 22 peter 50
4237 17 Sep 22 peter 51   double margin = 1e-3;
4237 17 Sep 22 peter 52   double x = finder(func, -10, 10, RootFinder::Interval(margin, 0));
4237 17 Sep 22 peter 53   suite.out() << "x: " << x << "\n";
4237 17 Sep 22 peter 54   if (!suite.equal_fix(x, M_LN2, margin)) {
4237 17 Sep 22 peter 55     suite.err() << "error using Interval as stopper\n";
4237 17 Sep 22 peter 56     suite.add(false);
4237 17 Sep 22 peter 57   }
4237 17 Sep 22 peter 58
4237 17 Sep 22 peter 59   x = finder(func, -10, 10, RootFinder::Delta(0, 1e-5));
4237 17 Sep 22 peter 60   suite.out() << "x: " << x << "\n";
4237 17 Sep 22 peter 61
4237 17 Sep 22 peter 62   margin = 1e-3;
4237 17 Sep 22 peter 63   x = finder(func, -10, 10,
4237 17 Sep 22 peter 64              RootFinder::Residual<MyFunction>(func, margin));
4237 17 Sep 22 peter 65   suite.out() << "x: " << x << "\n";
4237 17 Sep 22 peter 66   if (!suite.equal_fix(func(x), 0, margin)) {
4237 17 Sep 22 peter 67     suite.err() << "error using Residual as stopper\n";
4237 17 Sep 22 peter 68     suite.add(false);
4237 17 Sep 22 peter 69   }
4237 17 Sep 22 peter 70 }
4237 17 Sep 22 peter 71
4237 17 Sep 22 peter 72
4237 17 Sep 22 peter 73 template<typename T>
4237 17 Sep 22 peter 74 void run_test(test::Suite& suite)
4237 17 Sep 22 peter 75 {
4237 17 Sep 22 peter 76   T finder;
4237 17 Sep 22 peter 77   run_test(finder, suite);
4237 17 Sep 22 peter 78 }
4237 17 Sep 22 peter 79
4237 17 Sep 22 peter 80
4237 17 Sep 22 peter 81 int main(int argc, char* argv[])
4237 17 Sep 22 peter 82 {
4237 17 Sep 22 peter 83   test::Suite suite(argc, argv);
4237 17 Sep 22 peter 84
4237 17 Sep 22 peter 85   suite.out() << "test Bisection\n";
4237 17 Sep 22 peter 86   run_test<Bisection>(suite);
4237 17 Sep 22 peter 87   suite.out() << "test Brent\n";
4237 17 Sep 22 peter 88   run_test<Brent>(suite);
4237 17 Sep 22 peter 89   suite.out() << "test FalsePosition\n";
4237 17 Sep 22 peter 90   run_test<FalsePosition>(suite);
4237 17 Sep 22 peter 91   return suite.return_value();
4237 17 Sep 22 peter 92 }