Alexandria 2.31.0
SDC-CH common library for the Euclid project
Loading...
Searching...
No Matches
SecantMethod.cpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2022 Euclid Science Ground Segment
3 *
4 * This library is free software; you can redistribute it and/or modify it under
5 * the terms of the GNU Lesser General Public License as published by the Free
6 * Software Foundation; either version 3.0 of the License, or (at your option)
7 * any later version.
8 *
9 * This library is distributed in the hope that it will be useful, but WITHOUT
10 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11 * FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
12 * details.
13 *
14 * You should have received a copy of the GNU Lesser General Public License
15 * along with this library; if not, write to the Free Software Foundation, Inc.,
16 * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
17 */
18
20#include <cmath>
21
22namespace Euclid {
23namespace MathUtils {
24
25SecantReturn secantMethod(const Function& func, double x0, double x1, const SecantParams& params) {
27
28 double y0 = func(x0);
29
30 if (!std::isfinite(y0)) {
32 }
33
34 for (ret.iterations = 0; ret.iterations < params.max_iter; ++ret.iterations) {
35 if (std::abs(y0) <= params.atol) {
36 break;
37 }
38
39 double y1 = func(x1);
40 if (!std::isfinite(y1)) {
42 break;
43 }
44 if (y1 == y0) {
46 break;
47 }
48
49 double x2 = x1 - y1 * ((x1 - x0) / (y1 - y0));
50
51 if (x2 <= params.min) {
52 ret.root = params.min;
54 break;
55 } else if (x2 >= params.max) {
56 ret.root = params.max;
58 break;
59 }
60 x0 = x1;
61 x1 = x2;
62 y0 = y1;
63 ret.root = x0;
64 }
65 return ret;
66}
67
68} // namespace MathUtils
69} // namespace Euclid
Interface class representing a function with an arbitrary number of parameters.
Definition Function.h:104
T isfinite(T... args)
SecantReturn secantMethod(const Function &func, double x0, double x1, const SecantParams &params=SecantParams{})
std::array< std::vector< double >, N > Coordinates
Used to pass the grid coordinates to interpn. Internally will make a copy of the required values.