From 5ab0a94abc94d4e2452bbbb3597dd72221107e0e Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 17 May 2023 23:50:16 -0600 Subject: [PATCH] Avoid making a variable 'public'. --- include/deal.II/sundials/arkode.h | 1 - source/sundials/arkode.cc | 502 +++++++++++++----------------- 2 files changed, 225 insertions(+), 278 deletions(-) diff --git a/include/deal.II/sundials/arkode.h b/include/deal.II/sundials/arkode.h index 5ad3d25363..331a8f6bb0 100644 --- a/include/deal.II/sundials/arkode.h +++ b/include/deal.II/sundials/arkode.h @@ -1087,7 +1087,6 @@ namespace SUNDIALS std::unique_ptr> linear_solver; std::unique_ptr> mass_solver; - public: /** * A pointer to any exception that may have been thrown in user-defined * call-backs and that we have to deal after the KINSOL function we call diff --git a/source/sundials/arkode.cc b/source/sundials/arkode.cc index bbd66ae247..dd29cf2188 100644 --- a/source/sundials/arkode.cc +++ b/source/sundials/arkode.cc @@ -121,251 +121,6 @@ namespace SUNDIALS } // namespace - namespace - { - template - int - explicit_function_callback(realtype tt, - N_Vector yy, - N_Vector yp, - void * user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_yy = internal::unwrap_nvector_const(yy); - auto *dst_yp = internal::unwrap_nvector(yp); - - return call_and_possibly_capture_exception(solver.explicit_function, - solver.pending_exception, - tt, - *src_yy, - *dst_yp); - } - - - - template - int - implicit_function_callback(realtype tt, - N_Vector yy, - N_Vector yp, - void * user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_yy = internal::unwrap_nvector_const(yy); - auto *dst_yp = internal::unwrap_nvector(yp); - - return call_and_possibly_capture_exception(solver.implicit_function, - solver.pending_exception, - tt, - *src_yy, - *dst_yp); - } - - - - template - int - jacobian_times_vector_callback(N_Vector v, - N_Vector Jv, - realtype t, - N_Vector y, - N_Vector fy, - void * user_data, - N_Vector) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_v = internal::unwrap_nvector_const(v); - auto *src_y = internal::unwrap_nvector_const(y); - auto *src_fy = internal::unwrap_nvector_const(fy); - - auto *dst_Jv = internal::unwrap_nvector(Jv); - - return call_and_possibly_capture_exception(solver.jacobian_times_vector, - solver.pending_exception, - *src_v, - *dst_Jv, - t, - *src_y, - *src_fy); - } - - - - template - int - jacobian_times_vector_setup_callback(realtype t, - N_Vector y, - N_Vector fy, - void * user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_y = internal::unwrap_nvector_const(y); - auto *src_fy = internal::unwrap_nvector_const(fy); - - return call_and_possibly_capture_exception(solver.jacobian_times_setup, - solver.pending_exception, - t, - *src_y, - *src_fy); - } - - - - template - int - solve_with_jacobian_callback(realtype t, - N_Vector y, - N_Vector fy, - N_Vector r, - N_Vector z, - realtype gamma, - realtype delta, - int lr, - void * user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_y = internal::unwrap_nvector_const(y); - auto *src_fy = internal::unwrap_nvector_const(fy); - auto *src_r = internal::unwrap_nvector_const(r); - - auto *dst_z = internal::unwrap_nvector(z); - - return call_and_possibly_capture_exception( - solver.jacobian_preconditioner_solve, - solver.pending_exception, - t, - *src_y, - *src_fy, - *src_r, - *dst_z, - gamma, - delta, - lr); - } - - - - template - int - jacobian_solver_setup_callback(realtype t, - N_Vector y, - N_Vector fy, - booleantype jok, - booleantype *jcurPtr, - realtype gamma, - void * user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_y = internal::unwrap_nvector_const(y); - auto *src_fy = internal::unwrap_nvector_const(fy); - - return call_and_possibly_capture_exception( - solver.jacobian_preconditioner_setup, - solver.pending_exception, - t, - *src_y, - *src_fy, - jok, - *jcurPtr, - gamma); - } - - - - template - int - mass_matrix_times_vector_callback(N_Vector v, - N_Vector Mv, - realtype t, - void * mtimes_data) - { - Assert(mtimes_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(mtimes_data); - - auto *src_v = internal::unwrap_nvector_const(v); - auto *dst_Mv = internal::unwrap_nvector(Mv); - - return call_and_possibly_capture_exception( - solver.mass_times_vector, solver.pending_exception, t, *src_v, *dst_Mv); - } - - - - template - int - mass_matrix_times_vector_setup_callback(realtype t, void *mtimes_data) - { - Assert(mtimes_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(mtimes_data); - - return call_and_possibly_capture_exception(solver.mass_times_setup, - solver.pending_exception, - t); - } - - - - template - int - solve_with_mass_matrix_callback(realtype t, - N_Vector r, - N_Vector z, - realtype delta, - int lr, - void * user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - auto *src_r = internal::unwrap_nvector_const(r); - auto *dst_z = internal::unwrap_nvector(z); - - return call_and_possibly_capture_exception( - solver.mass_preconditioner_solve, - solver.pending_exception, - t, - *src_r, - *dst_z, - delta, - lr); - } - - - - template - int - mass_matrix_solver_setup_callback(realtype t, void *user_data) - { - Assert(user_data != nullptr, ExcInternalError()); - ARKode &solver = - *static_cast *>(user_data); - - return call_and_possibly_capture_exception( - solver.mass_preconditioner_setup, solver.pending_exception, t); - } - } // namespace - template ARKode::ARKode(const AdditionalData &data) @@ -586,14 +341,48 @@ namespace SUNDIALS Assert(explicit_function || implicit_function, ExcFunctionNotProvided("explicit_function || implicit_function")); - arkode_mem = ARKStepCreate( - explicit_function ? &explicit_function_callback : nullptr, - implicit_function ? &implicit_function_callback : nullptr, - current_time, - initial_condition_nvector + auto explicit_function_callback = + [](realtype tt, N_Vector yy, N_Vector yp, void *user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_yy = internal::unwrap_nvector_const(yy); + auto *dst_yp = internal::unwrap_nvector(yp); + + return call_and_possibly_capture_exception(solver.explicit_function, + solver.pending_exception, + tt, + *src_yy, + *dst_yp); + }; + + + auto implicit_function_callback = + [](realtype tt, N_Vector yy, N_Vector yp, void *user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_yy = internal::unwrap_nvector_const(yy); + auto *dst_yp = internal::unwrap_nvector(yp); + + return call_and_possibly_capture_exception(solver.implicit_function, + solver.pending_exception, + tt, + *src_yy, + *dst_yp); + }; + + arkode_mem = ARKStepCreate(explicit_function ? explicit_function_callback : + ARKRhsFn(nullptr), + implicit_function ? implicit_function_callback : + ARKRhsFn(nullptr), + current_time, + initial_condition_nvector # if DEAL_II_SUNDIALS_VERSION_GTE(6, 0, 0) - , - arkode_ctx + , + arkode_ctx # endif ); Assert(arkode_mem != nullptr, ExcInternalError()); @@ -694,21 +483,120 @@ namespace SUNDIALS } status = ARKStepSetLinearSolver(arkode_mem, sun_linear_solver, nullptr); AssertARKode(status); - status = ARKStepSetJacTimes( - arkode_mem, - jacobian_times_setup ? - jacobian_times_vector_setup_callback : - nullptr, - jacobian_times_vector_callback); + + auto jacobian_times_vector_callback = [](N_Vector v, + N_Vector Jv, + realtype t, + N_Vector y, + N_Vector fy, + void * user_data, + N_Vector) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_v = internal::unwrap_nvector_const(v); + auto *src_y = internal::unwrap_nvector_const(y); + auto *src_fy = internal::unwrap_nvector_const(fy); + + auto *dst_Jv = internal::unwrap_nvector(Jv); + + return call_and_possibly_capture_exception( + solver.jacobian_times_vector, + solver.pending_exception, + *src_v, + *dst_Jv, + t, + *src_y, + *src_fy); + }; + + auto jacobian_times_vector_setup_callback = + [](realtype t, N_Vector y, N_Vector fy, void *user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_y = internal::unwrap_nvector_const(y); + auto *src_fy = internal::unwrap_nvector_const(fy); + + return call_and_possibly_capture_exception( + solver.jacobian_times_setup, + solver.pending_exception, + t, + *src_y, + *src_fy); + }; + status = ARKStepSetJacTimes(arkode_mem, + jacobian_times_setup ? + jacobian_times_vector_setup_callback : + ARKLsJacTimesSetupFn(nullptr), + jacobian_times_vector_callback); AssertARKode(status); if (jacobian_preconditioner_solve) { - status = ARKStepSetPreconditioner( - arkode_mem, - jacobian_preconditioner_setup ? - jacobian_solver_setup_callback : - nullptr, - solve_with_jacobian_callback); + auto solve_with_jacobian_callback = [](realtype t, + N_Vector y, + N_Vector fy, + N_Vector r, + N_Vector z, + realtype gamma, + realtype delta, + int lr, + void * user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_y = internal::unwrap_nvector_const(y); + auto *src_fy = internal::unwrap_nvector_const(fy); + auto *src_r = internal::unwrap_nvector_const(r); + + auto *dst_z = internal::unwrap_nvector(z); + + return call_and_possibly_capture_exception( + solver.jacobian_preconditioner_solve, + solver.pending_exception, + t, + *src_y, + *src_fy, + *src_r, + *dst_z, + gamma, + delta, + lr); + }; + + auto jacobian_solver_setup_callback = [](realtype t, + N_Vector y, + N_Vector fy, + booleantype jok, + booleantype *jcurPtr, + realtype gamma, + void *user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_y = internal::unwrap_nvector_const(y); + auto *src_fy = internal::unwrap_nvector_const(fy); + + return call_and_possibly_capture_exception( + solver.jacobian_preconditioner_setup, + solver.pending_exception, + t, + *src_y, + *src_fy, + jok, + *jcurPtr, + gamma); + }; + + status = ARKStepSetPreconditioner(arkode_mem, + jacobian_preconditioner_setup ? + jacobian_solver_setup_callback : + ARKLsPrecSetupFn(nullptr), + solve_with_jacobian_callback); AssertARKode(status); } if (data.implicit_function_is_linear) @@ -801,23 +689,83 @@ namespace SUNDIALS nullptr, mass_time_dependent); AssertARKode(status); - status = ARKStepSetMassTimes( - arkode_mem, - mass_times_setup ? - mass_matrix_times_vector_setup_callback : - nullptr, - mass_matrix_times_vector_callback, - this); + + auto mass_matrix_times_vector_setup_callback = + [](realtype t, void *mtimes_data) -> int { + Assert(mtimes_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(mtimes_data); + + return call_and_possibly_capture_exception(solver.mass_times_setup, + solver.pending_exception, + t); + }; + + auto mass_matrix_times_vector_callback = + [](N_Vector v, N_Vector Mv, realtype t, void *mtimes_data) -> int { + Assert(mtimes_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(mtimes_data); + + auto *src_v = internal::unwrap_nvector_const(v); + auto *dst_Mv = internal::unwrap_nvector(Mv); + + return call_and_possibly_capture_exception(solver.mass_times_vector, + solver.pending_exception, + t, + *src_v, + *dst_Mv); + }; + + status = ARKStepSetMassTimes(arkode_mem, + mass_times_setup ? + mass_matrix_times_vector_setup_callback : + ARKLsMassTimesSetupFn(nullptr), + mass_matrix_times_vector_callback, + this); AssertARKode(status); if (mass_preconditioner_solve) { - status = ARKStepSetMassPreconditioner( - arkode_mem, - mass_preconditioner_setup ? - mass_matrix_solver_setup_callback : - nullptr, - solve_with_mass_matrix_callback); + auto mass_matrix_solver_setup_callback = + [](realtype t, void *user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + return call_and_possibly_capture_exception( + solver.mass_preconditioner_setup, solver.pending_exception, t); + }; + + auto solve_with_mass_matrix_callback = [](realtype t, + N_Vector r, + N_Vector z, + realtype delta, + int lr, + void *user_data) -> int { + Assert(user_data != nullptr, ExcInternalError()); + ARKode &solver = + *static_cast *>(user_data); + + auto *src_r = internal::unwrap_nvector_const(r); + auto *dst_z = internal::unwrap_nvector(z); + + return call_and_possibly_capture_exception( + solver.mass_preconditioner_solve, + solver.pending_exception, + t, + *src_r, + *dst_z, + delta, + lr); + }; + + status = + ARKStepSetMassPreconditioner(arkode_mem, + mass_preconditioner_setup ? + mass_matrix_solver_setup_callback : + ARKLsMassPrecSetupFn(nullptr), + solve_with_mass_matrix_callback); AssertARKode(status); } } -- 2.39.5