In the secant function, a 1-D Newton-Raphson solver is implemented. An initial guess for the solution is provided. Note that this is different than the Secant function because if something goes out of bounds, it will just make its best guess. @param f A pointer to an instance of the FuncWrapper1D class that implements the call() function @param x0 The initial guess for the solutionh @param dx Th
| 427 | @returns If no errors are found, the solution, otherwise the value _HUGE, the value for infinity |
| 428 | */ |
| 429 | double ExtrapolatingSecant(FuncWrapper1D* f, double x0, double dx, double tol, int maxiter) { |
| 430 | #if defined(COOLPROP_DEEP_DEBUG) |
| 431 | static std::vector<double> xlog, flog; |
| 432 | xlog.clear(); |
| 433 | flog.clear(); |
| 434 | #endif |
| 435 | |
| 436 | // Initialization |
| 437 | double x1 = 0, x2 = 0, x3 = 0, y0 = 0, y1 = 0, y2 = 0, x = x0, fval = 999; |
| 438 | f->iter = 1; |
| 439 | f->errstring.clear(); |
| 440 | |
| 441 | // The relaxation factor (less than 1 for smaller steps) |
| 442 | double omega = f->options.get_double("omega", 1.0); |
| 443 | |
| 444 | if (std::abs(dx) == 0) { |
| 445 | f->errstring = "dx cannot be zero"; |
| 446 | return _HUGE; |
| 447 | } |
| 448 | while (f->iter <= 2 || std::abs(fval) > tol) { |
| 449 | if (f->iter == 1) { |
| 450 | x1 = x0; |
| 451 | x = x1; |
| 452 | } |
| 453 | if (f->iter == 2) { |
| 454 | x2 = x0 + dx; |
| 455 | x = x2; |
| 456 | } |
| 457 | if (f->iter > 2) { |
| 458 | x = x2; |
| 459 | } |
| 460 | |
| 461 | if (f->input_not_in_range(x)) { |
| 462 | throw ValueError(format("Input [%g] is out of range", x)); |
| 463 | } |
| 464 | |
| 465 | fval = f->call(x); |
| 466 | |
| 467 | #if defined(COOLPROP_DEEP_DEBUG) |
| 468 | xlog.push_back(x); |
| 469 | flog.push_back(fval); |
| 470 | #endif |
| 471 | |
| 472 | if (!ValidNumber(fval)) { |
| 473 | if (f->iter == 1) { |
| 474 | return x; |
| 475 | } else { |
| 476 | return x2 - omega * y1 / (y1 - y0) * (x2 - x1); |
| 477 | } |
| 478 | }; |
| 479 | if (f->iter == 1) { |
| 480 | y1 = fval; |
| 481 | } |
| 482 | if (f->iter > 1) { |
| 483 | double deltax = x2 - x1; |
| 484 | if (std::abs(deltax) < 1e-14) { |
| 485 | return x; |
| 486 | } |
no test coverage detected