From 1032025f3a076f97fc42e43fbd927278b6aec327 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Mon, 17 May 2021 17:40:55 +0200 Subject: [PATCH] Fix Kinsol interface. --- include/deal.II/sundials/kinsol.h | 43 ++------ source/sundials/kinsol.cc | 103 ++++++++++-------- tests/sundials/kinsol_01.cc | 47 +++++--- ...ol_01.with_sundials=true.expect=run.output | 3 - .../kinsol_01.with_sundials=true.output | 3 + tests/sundials/kinsol_02.cc | 62 ++++++----- ...ol_02.with_sundials=true.expect=run.output | 3 - .../kinsol_02.with_sundials=true.output | 3 + tests/sundials/kinsol_03.cc | 2 +- tests/sundials/kinsol_03_new_interface.cc | 2 +- tests/sundials/kinsol_04.cc | 2 +- tests/sundials/kinsol_04_new_interface.cc | 2 +- tests/sundials/kinsol_05.cc | 2 +- tests/sundials/kinsol_05_new_interface.cc | 2 +- tests/sundials/kinsol_fixed_point.prm | 16 +++ .../{kinsol_01.prm => kinsol_linesearch.prm} | 0 tests/sundials/kinsol_newton.prm | 16 +++ tests/sundials/kinsol_picard.prm | 16 +++ 18 files changed, 190 insertions(+), 137 deletions(-) delete mode 100644 tests/sundials/kinsol_01.with_sundials=true.expect=run.output create mode 100644 tests/sundials/kinsol_01.with_sundials=true.output delete mode 100644 tests/sundials/kinsol_02.with_sundials=true.expect=run.output create mode 100644 tests/sundials/kinsol_02.with_sundials=true.output create mode 100644 tests/sundials/kinsol_fixed_point.prm rename tests/sundials/{kinsol_01.prm => kinsol_linesearch.prm} (100%) create mode 100644 tests/sundials/kinsol_newton.prm create mode 100644 tests/sundials/kinsol_picard.prm diff --git a/include/deal.II/sundials/kinsol.h b/include/deal.II/sundials/kinsol.h index 939215b05f..9bf6c55180 100644 --- a/include/deal.II/sundials/kinsol.h +++ b/include/deal.II/sundials/kinsol.h @@ -66,10 +66,7 @@ namespace SUNDIALS * * KINSOL's Newton solver employs the inexact Newton method. As this solver * is intended mainly for large systems, the user is required to provide - * their own solver function. If a solver function is not provided, the - * internal dense solver of KINSOL is used. Be warned that this solver - * computes the Jacobian approximately, and may be efficient only for small - * systems. + * their own solver function. * * At the highest level, KINSOL implements the following iteration * scheme: @@ -136,22 +133,6 @@ namespace SUNDIALS * algorithm is that the full Newton step tends to be taken close to the * solution. * - * As a user option, KINSOL permits the application of inequality - * constraints, $u_i > 0$ and $u_i < 0$, as well as $u_i \geq 0$ and $u_i - * \leq 0$, where $u_i$ is the $i$-th component of $u$. Any such constraint, - * or no constraint, may be imposed on each component by providing the - * optional functions - * - get_lower_than_zero_constrained_entries() - * - get_greater_than_zero_constrained_entries() - * - get_lower_equal_than_zero_constrained_entries() - * - get_greater_or_equal_than_zero_constrained_entries() - * - * KINSOL will reduce step lengths in order to ensure that no constraint is - * violated. Specifically, if a new Newton iterate will violate a constraint, - * the maximum step length along the Newton direction that will satisfy all - * constraints is found, and $\delta_n$ is scaled to take a step of that - * length. - * * The basic fixed-point iteration scheme implemented in KINSOL is given by: * - Set $u_0 =$ an initial guess * - For $n = 0, 1, 2, \dots$ until convergence do: @@ -180,20 +161,17 @@ namespace SUNDIALS * or * - iteration_function; * - * Specifying residual() allows the user to use Newton strategies (i.e., - * $F(u)=0$ will be solved), while specifying iteration_function(), fixed - * point iteration or Picard iteration will be used (i.e., $G(u)=u$ will be - * solved). + * Specifying residual() allows the user to use Newton and Picard strategies + * (i.e., $F(u)=0$ will be solved), while specifying iteration_function(), a + * fixed point iteration will be used (i.e., $G(u)=u$ will be solved). * - * If the use of a Newton method is desired, then the user should also supply - * - solve_jacobian_system; + * If the use of a Newton or Picard method is desired, then the user should + * also supply + * - solve_jacobian_system or solve_with_jacobian; * and optionally * - setup_jacobian; * - * If the solve_jacobian_system() function is not supplied, then KINSOL will - * use its internal dense solver for Newton methods, with approximate - * Jacobian. This may be very expensive for large systems. Fixed point - * iteration does not require the solution of any linear system. + * Fixed point iteration does not require the solution of any linear system. * * Also the following functions could be rewritten, to provide additional * scaling factors for both the solution and the residual evaluation during @@ -711,11 +689,6 @@ namespace SUNDIALS */ void *kinsol_mem; - /** - * MPI communicator. SUNDIALS solver runs happily in parallel. - */ - MPI_Comm communicator; - /** * Memory pool of vectors. */ diff --git a/source/sundials/kinsol.cc b/source/sundials/kinsol.cc index cb12781b0e..bdbac29035 100644 --- a/source/sundials/kinsol.cc +++ b/source/sundials/kinsol.cc @@ -133,7 +133,7 @@ namespace SUNDIALS { template int - residual_or_iteration_callback(N_Vector yy, N_Vector FF, void *user_data) + residual_callback(N_Vector yy, N_Vector FF, void *user_data) { KINSOL &solver = *static_cast *>(user_data); @@ -144,7 +144,26 @@ namespace SUNDIALS int err = 0; if (solver.residual) err = solver.residual(*src_yy, *dst_FF); - else if (solver.iteration_function) + else + Assert(false, ExcInternalError()); + + return err; + } + + + + template + int + iteration_callback(N_Vector yy, N_Vector FF, void *user_data) + { + KINSOL &solver = + *static_cast *>(user_data); + + auto *src_yy = internal::unwrap_nvector_const(yy); + auto *dst_FF = internal::unwrap_nvector(FF); + + int err = 0; + if (solver.iteration_function) err = solver.iteration_function(*src_yy, *dst_FF); else Assert(false, ExcInternalError()); @@ -281,20 +300,15 @@ namespace SUNDIALS return err; } } - # endif } // namespace template - KINSOL::KINSOL(const AdditionalData &data, - const MPI_Comm & mpi_comm) + KINSOL::KINSOL(const AdditionalData &data, const MPI_Comm &) : data(data) , kinsol_mem(nullptr) - , communicator(is_serial_vector::value ? - MPI_COMM_SELF : - Utilities::MPI::duplicate_communicator(mpi_comm)) { set_functions_to_trigger_an_assert(); } @@ -306,15 +320,6 @@ namespace SUNDIALS { if (kinsol_mem) KINFree(&kinsol_mem); - -# ifdef DEAL_II_WITH_MPI - if (is_serial_vector::value == false) - { - const int ierr = MPI_Comm_free(&communicator); - (void)ierr; - AssertNothrow(ierr == MPI_SUCCESS, ExcMPI(ierr)); - } -# endif } @@ -323,8 +328,6 @@ namespace SUNDIALS unsigned int KINSOL::solve(VectorType &initial_guess_and_solution) { - unsigned int system_size = initial_guess_and_solution.size(); - NVectorView u_scale, f_scale; VectorType u_scale_temp, f_scale_temp; @@ -347,6 +350,20 @@ namespace SUNDIALS f_scale = internal::make_nvector_view(f_scale_temp); } + // Make sure we have what we need + if (data.strategy == AdditionalData::fixed_point) + { + Assert(iteration_function, + ExcFunctionNotProvided("iteration_function")); + } + else + { + Assert(residual, ExcFunctionNotProvided("residual")); + Assert(solve_jacobian_system || solve_with_jacobian, + ExcFunctionNotProvided( + "solve_jacobian_system || solve_with_jacobian")); + } + auto solution = internal::make_nvector_view(initial_guess_and_solution); if (kinsol_mem) @@ -354,17 +371,26 @@ namespace SUNDIALS kinsol_mem = KINCreate(); - int status = - KINInit(kinsol_mem, residual_or_iteration_callback, solution); + int status = 0; (void)status; - AssertKINSOL(status); status = KINSetUserData(kinsol_mem, static_cast(this)); AssertKINSOL(status); + // This must be called before KINSetMAA status = KINSetNumMaxIters(kinsol_mem, data.maximum_non_linear_iterations); AssertKINSOL(status); + // From the manual: this must be called BEFORE KINInit + status = KINSetMAA(kinsol_mem, data.anderson_subspace_size); + AssertKINSOL(status); + + if (data.strategy == AdditionalData::fixed_point) + status = KINInit(kinsol_mem, iteration_callback, solution); + else + status = KINInit(kinsol_mem, residual_callback, solution); + AssertKINSOL(status); + status = KINSetFuncNormTol(kinsol_mem, data.function_tolerance); AssertKINSOL(status); @@ -383,9 +409,6 @@ namespace SUNDIALS status = KINSetMaxBetaFails(kinsol_mem, data.maximum_beta_failures); AssertKINSOL(status); - status = KINSetMAA(kinsol_mem, data.anderson_subspace_size); - AssertKINSOL(status); - status = KINSetRelErrFunc(kinsol_mem, data.dq_relative_error); AssertKINSOL(status); @@ -398,8 +421,11 @@ namespace SUNDIALS { /* interface up to and including 4.0 */ # if DEAL_II_SUNDIALS_VERSION_LT(4, 1, 0) - auto KIN_mem = static_cast(kinsol_mem); - KIN_mem->kin_lsolve = solve_with_jacobian_callback; + auto KIN_mem = static_cast(kinsol_mem); + // Old version only works with solve_jacobian_system + Assert(solve_jacobian_system, + ExcFunctionNotProvided("solve_jacobian_system")) + KIN_mem->kin_lsolve = solve_with_jacobian_callback; if (setup_jacobian) // user assigned a function object to the Jacobian // set-up slot KIN_mem->kin_lsetup = setup_jacobian_callback; @@ -486,21 +512,6 @@ namespace SUNDIALS } # endif } - else - { - J = SUNDenseMatrix(system_size, system_size); - LS = SUNDenseLinearSolver(u_scale, J); - status = KINDlsSetLinearSolver(kinsol_mem, LS, J); - AssertKINSOL(status); - } - - if (data.strategy == AdditionalData::newton || - data.strategy == AdditionalData::linesearch) - Assert(residual, ExcFunctionNotProvided("residual")); - - if (data.strategy == AdditionalData::fixed_point || - data.strategy == AdditionalData::picard) - Assert(iteration_function, ExcFunctionNotProvided("iteration_function")); // call to KINSol status = KINSol(kinsol_mem, solution, data.strategy, u_scale, f_scale); @@ -510,8 +521,10 @@ namespace SUNDIALS status = KINGetNumNonlinSolvIters(kinsol_mem, &nniters); AssertKINSOL(status); - SUNMatDestroy(J); - SUNLinSolFree(LS); + if (J != nullptr) + SUNMatDestroy(J); + if (LS != nullptr) + SUNLinSolFree(LS); KINFree(&kinsol_mem); return static_cast(nniters); @@ -524,6 +537,8 @@ namespace SUNDIALS reinit_vector = [](VectorType &) { AssertThrow(false, ExcFunctionNotProvided("reinit_vector")); }; + + setup_jacobian = [](const VectorType &, const VectorType &) { return 0; }; } template class KINSOL>; diff --git a/tests/sundials/kinsol_01.cc b/tests/sundials/kinsol_01.cc index 15e9005bd7..028491fcea 100644 --- a/tests/sundials/kinsol_01.cc +++ b/tests/sundials/kinsol_01.cc @@ -22,15 +22,22 @@ #include "../tests.h" -// Solve a nonlinear system but provide only residual function. KINSOL -// then uses its internal solvers which are based on a -// finite-difference approximation to the Jacobian and a direct -// solver. +// Solve a nonlinear system using fixed point iteration, and Anderson +// acceleration /** - * Solve the non linear problem - * - * F(u) = 0 , where f_i(u) = u_i^2 - i^2, 0 <= i < N + * The following is a simple example problem, with the coding + * needed for its solution by the accelerated fixed point solver in + * KINSOL. + * The problem is from chemical kinetics, and consists of solving + * the first time step in a Backward Euler solution for the + * following three rate equations: + * dy1/dt = -.04*y1 + 1.e4*y2*y3 + * dy2/dt = .04*y1 - 1.e4*y2*y3 - 3.e2*(y2)^2 + * dy3/dt = 3.e2*(y2)^2 + * on the interval from t = 0.0 to t = 0.1, with initial + * conditions: y1 = 1.0, y2 = y3 = 0. The problem is stiff. + * Run statistics (optional outputs) are printed at the end. */ int main(int argc, char **argv) @@ -46,31 +53,35 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_fixed_point.prm"); prm.parse_input(ifile); // Size of the problem - unsigned int N = 10; + unsigned int N = 3; SUNDIALS::KINSOL kinsol(data); kinsol.reinit_vector = [N](VectorType &v) { v.reinit(N); }; - kinsol.residual = [](const VectorType &u, VectorType &F) -> int { - for (unsigned int i = 0; i < u.size(); ++i) - F[i] = u[i] * u[i] - (i + 1) * (i + 1); - return 0; - }; + // Robert example + kinsol.iteration_function = [](const VectorType &u, VectorType &F) -> int { + const double dstep = 0.1; + const double y10 = 1.0; + const double y20 = 0.0; + const double y30 = 0.0; + const double yd1 = dstep * (-0.04 * u[0] + 1.0e4 * u[1] * u[2]); + const double yd3 = dstep * 3.0e2 * u[1] * u[1]; + + F[0] = yd1 + y10; + F[1] = -yd1 - yd3 + y20; + F[2] = yd3 + y30; - kinsol.iteration_function = [](const VectorType &u, VectorType &F) -> int { - for (unsigned int i = 0; i < u.size(); ++i) - F[i] = u[i] * u[i] - i * i - u[i]; return 0; }; VectorType v(N); - v = 1.0; + v[0] = 1; auto niter = kinsol.solve(v); v.print(deallog.get_file_stream()); deallog << "Converged in " << niter << " iterations." << std::endl; diff --git a/tests/sundials/kinsol_01.with_sundials=true.expect=run.output b/tests/sundials/kinsol_01.with_sundials=true.expect=run.output deleted file mode 100644 index 2a13a0ef6e..0000000000 --- a/tests/sundials/kinsol_01.with_sundials=true.expect=run.output +++ /dev/null @@ -1,3 +0,0 @@ - -1.00000 2.00000 3.00000 4.00000 5.00000 6.00000 7.00000 8.00000 9.00000 10.0000 -DEAL::Converged in 23 iterations. diff --git a/tests/sundials/kinsol_01.with_sundials=true.output b/tests/sundials/kinsol_01.with_sundials=true.output new file mode 100644 index 0000000000..61ccde6ace --- /dev/null +++ b/tests/sundials/kinsol_01.with_sundials=true.output @@ -0,0 +1,3 @@ + +9.968e-01 2.953e-03 2.616e-04 +DEAL::Converged in 8 iterations. diff --git a/tests/sundials/kinsol_02.cc b/tests/sundials/kinsol_02.cc index 5a01ab35c6..dc56917859 100644 --- a/tests/sundials/kinsol_02.cc +++ b/tests/sundials/kinsol_02.cc @@ -22,23 +22,18 @@ #include "../tests.h" -// Solve a nonlinear system but provide only residual function. KINSOL -// then uses its internal solvers which are based on a -// finite-difference approximation to the Jacobian and a direct -// solver. +// Solve a nonlinear system in the form accepted by Picard iteration. // -// Compared to the _01 test, this is simply a more complicated function: // We solve the nonlinear problem // -// F(u) = 0 +// F(u) = L u + N(u) = 0 // -// with a 2-dimensional vector u and where +// where L is a constant matrix, and N(u) is non linear. // -// F(u) = [ cos(u1 + u2) - 1 ] -> u1=-u2 -// [ sin(u1 - u2) ] -> u1=u2 +// We set L = id and +// +// N_i(u) = .1*u_i^2 - i - 1 // -// In other words, we need to find the solution u1=u2=0. - int main(int argc, char **argv) { @@ -49,43 +44,54 @@ main(int argc, char **argv) using VectorType = Vector; + // Size of the problem + unsigned int N = 2; + + FullMatrix L(N, N); + L(0, 0) = 1; + L(1, 1) = 1; + L(0, 1) = 1; + + FullMatrix Linv(N, N); + Linv.invert(L); + SUNDIALS::KINSOL::AdditionalData data; ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_picard.prm"); prm.parse_input(ifile); - // Size of the problem - unsigned int N = 2; - SUNDIALS::KINSOL kinsol(data); kinsol.reinit_vector = [N](VectorType &v) { v.reinit(N); }; - kinsol.residual = [](const VectorType &u, VectorType &F) -> int { - F(0) = std::cos(u[0] + u[1]) - 1; - F(1) = std::sin(u[0] - u[1]); + kinsol.residual = [&](const VectorType &u, VectorType &F) -> int { + F = u; + + F[0] += .1 * u[0] * u[0] - 1; + F[1] += .1 * u[1] * u[1] - 2; return 0; }; + kinsol.solve_with_jacobian = + [&](const VectorType &rhs, VectorType &dst, double) -> int { + dst = rhs; + return 0; + }; - kinsol.iteration_function = [](const VectorType &u, VectorType &F) -> int { - // We want a Newton-type scheme, not a fixed point iteration. So we - // shouldn't get into this function. - std::abort(); - - // But if anyone wanted to see how it would look like: - F(0) = std::cos(u[0] + u[1]) - 1 - u[0]; - F(1) = std::sin(u[0] - u[1]) - u[1]; + kinsol.solve_jacobian_system = [&](const VectorType &, + const VectorType &, + const VectorType &rhs, + VectorType & dst) -> int { + dst = rhs; return 0; }; VectorType v(N); - v(0) = 0.5; - v(1) = 1.234; auto niter = kinsol.solve(v); + v.print(deallog.get_file_stream()); deallog << "Converged in " << niter << " iterations." << std::endl; } diff --git a/tests/sundials/kinsol_02.with_sundials=true.expect=run.output b/tests/sundials/kinsol_02.with_sundials=true.expect=run.output deleted file mode 100644 index ee3d8fab69..0000000000 --- a/tests/sundials/kinsol_02.with_sundials=true.expect=run.output +++ /dev/null @@ -1,3 +0,0 @@ - -9.761e-04 9.761e-04 -DEAL::Converged in 27 iterations. diff --git a/tests/sundials/kinsol_02.with_sundials=true.output b/tests/sundials/kinsol_02.with_sundials=true.output new file mode 100644 index 0000000000..38ba409016 --- /dev/null +++ b/tests/sundials/kinsol_02.with_sundials=true.output @@ -0,0 +1,3 @@ + +9.161e-01 1.708e+00 +DEAL::Converged in 7 iterations. diff --git a/tests/sundials/kinsol_03.cc b/tests/sundials/kinsol_03.cc index d5368565cd..1c2c5ebee6 100644 --- a/tests/sundials/kinsol_03.cc +++ b/tests/sundials/kinsol_03.cc @@ -58,7 +58,7 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_newton.prm"); prm.parse_input(ifile); // Size of the problem diff --git a/tests/sundials/kinsol_03_new_interface.cc b/tests/sundials/kinsol_03_new_interface.cc index f473627aa8..a0d79e3993 100644 --- a/tests/sundials/kinsol_03_new_interface.cc +++ b/tests/sundials/kinsol_03_new_interface.cc @@ -58,7 +58,7 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_linesearch.prm"); prm.parse_input(ifile); // Size of the problem diff --git a/tests/sundials/kinsol_04.cc b/tests/sundials/kinsol_04.cc index f7c8baaf41..b0f079685d 100644 --- a/tests/sundials/kinsol_04.cc +++ b/tests/sundials/kinsol_04.cc @@ -70,7 +70,7 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_linesearch.prm"); prm.parse_input(ifile); // Size of the problem diff --git a/tests/sundials/kinsol_04_new_interface.cc b/tests/sundials/kinsol_04_new_interface.cc index 73ab99866a..da932551ed 100644 --- a/tests/sundials/kinsol_04_new_interface.cc +++ b/tests/sundials/kinsol_04_new_interface.cc @@ -70,7 +70,7 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_linesearch.prm"); prm.parse_input(ifile); // Size of the problem diff --git a/tests/sundials/kinsol_05.cc b/tests/sundials/kinsol_05.cc index 5c5c3464f0..a5dfbed054 100644 --- a/tests/sundials/kinsol_05.cc +++ b/tests/sundials/kinsol_05.cc @@ -64,7 +64,7 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_linesearch.prm"); prm.parse_input(ifile); // Update the Jacobian in each iteration: diff --git a/tests/sundials/kinsol_05_new_interface.cc b/tests/sundials/kinsol_05_new_interface.cc index 83648603e1..4d2e05ecf5 100644 --- a/tests/sundials/kinsol_05_new_interface.cc +++ b/tests/sundials/kinsol_05_new_interface.cc @@ -64,7 +64,7 @@ main(int argc, char **argv) ParameterHandler prm; data.add_parameters(prm); - std::ifstream ifile(SOURCE_DIR "/kinsol_01.prm"); + std::ifstream ifile(SOURCE_DIR "/kinsol_linesearch.prm"); prm.parse_input(ifile); // Update the Jacobian in each iteration: diff --git a/tests/sundials/kinsol_fixed_point.prm b/tests/sundials/kinsol_fixed_point.prm new file mode 100644 index 0000000000..a64b3cac28 --- /dev/null +++ b/tests/sundials/kinsol_fixed_point.prm @@ -0,0 +1,16 @@ +set Function norm stopping tolerance = 1e-10 +set Maximum number of nonlinear iterations = 200 +set Scaled step stopping tolerance = 1e-10 +set Solution strategy = fixed_point +subsection Fixed point and Picard parameters + set Anderson acceleration subspace size = 2 +end +subsection Linesearch parameters + set Maximum number of beta-condition failures = 0 +end +subsection Newton parameters + set Maximum allowable scaled length of the Newton step = 0.000000 + set Maximum iterations without matrix setup = 0 + set No initial matrix setup = false + set Relative error for different quotient computation = 0.000000 +end diff --git a/tests/sundials/kinsol_01.prm b/tests/sundials/kinsol_linesearch.prm similarity index 100% rename from tests/sundials/kinsol_01.prm rename to tests/sundials/kinsol_linesearch.prm diff --git a/tests/sundials/kinsol_newton.prm b/tests/sundials/kinsol_newton.prm new file mode 100644 index 0000000000..3e5ce272a8 --- /dev/null +++ b/tests/sundials/kinsol_newton.prm @@ -0,0 +1,16 @@ +set Function norm stopping tolerance = 0.000000 +set Maximum number of nonlinear iterations = 200 +set Scaled step stopping tolerance = 0.000000 +set Solution strategy = newton +subsection Fixed point and Picard parameters + set Anderson acceleration subspace size = 5 +end +subsection Linesearch parameters + set Maximum number of beta-condition failures = 0 +end +subsection Newton parameters + set Maximum allowable scaled length of the Newton step = 0.000000 + set Maximum iterations without matrix setup = 0 + set No initial matrix setup = false + set Relative error for different quotient computation = 0.000000 +end diff --git a/tests/sundials/kinsol_picard.prm b/tests/sundials/kinsol_picard.prm new file mode 100644 index 0000000000..f684cc3e2b --- /dev/null +++ b/tests/sundials/kinsol_picard.prm @@ -0,0 +1,16 @@ +set Function norm stopping tolerance = 1e-10 +set Maximum number of nonlinear iterations = 200 +set Scaled step stopping tolerance = 1e-10 +set Solution strategy = picard +subsection Fixed point and Picard parameters + set Anderson acceleration subspace size = 2 +end +subsection Linesearch parameters + set Maximum number of beta-condition failures = 0 +end +subsection Newton parameters + set Maximum allowable scaled length of the Newton step = 0.000000 + set Maximum iterations without matrix setup = 0 + set No initial matrix setup = false + set Relative error for different quotient computation = 0.000000 +end -- 2.39.5