--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#ifndef dealii_sundials_utilities_h
+#define dealii_sundials_utilities_h
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_SUNDIALS
+# include <exception>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace SUNDIALS
+{
+ namespace Utilities
+ {
+ /**
+ * A function that calls the function object given by its first argument
+ * with the set of arguments following at the end. If the call returns
+ * regularly, the current function returns zero to indicate success. If the
+ * call fails with an exception of type RecoverableUserCallbackError, then
+ * the current function returns 1 to indicate that the called function
+ * object thought the error it encountered is recoverable. If the call fails
+ * with any other exception, then the current function returns with an error
+ * code of -1. In each of the last two cases, the exception thrown by `f`
+ * is captured and `eptr` is set to the exception. In case of success,
+ * `eptr` is set to `nullptr`.
+ */
+ template <typename F, typename... Args>
+ int
+ call_and_possibly_capture_exception(const F & f,
+ std::exception_ptr &eptr,
+ Args &&...args)
+ {
+ // See whether there is already something in the exception pointer
+ // variable. This can only happen if we had previously had
+ // a recoverable exception, and the underlying library actually
+ // did recover successfully. In that case, we can abandon the
+ // exception previously thrown. If eptr contains anything other,
+ // then we really don't know how that could have happened, and
+ // should probably bail out:
+ if (eptr)
+ {
+ try
+ {
+ std::rethrow_exception(eptr);
+ }
+ catch (const RecoverableUserCallbackError &)
+ {
+ // ok, ignore, but reset the pointer
+ eptr = nullptr;
+ }
+ catch (...)
+ {
+ // uh oh:
+ AssertThrow(false, ExcInternalError());
+ }
+ }
+
+ // Call the function and if that succeeds, return zero:
+ try
+ {
+ f(std::forward<Args>(args)...);
+ eptr = nullptr;
+ return 0;
+ }
+ // If the call failed with a recoverable error, then
+ // ignore the exception for now (but store a pointer to it)
+ // and return a positive return value (+1). If the underlying
+ // implementation manages to recover
+ catch (const RecoverableUserCallbackError &)
+ {
+ eptr = std::current_exception();
+ return 1;
+ }
+ // For any other exception, capture the exception and
+ // return -1:
+ catch (const std::exception &)
+ {
+ eptr = std::current_exception();
+ return -1;
+ }
+ }
+ } // namespace Utilities
+} // namespace SUNDIALS
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_WITH_SUNDIALS
+
+#endif
# include <deal.II/sundials/n_vector.h>
# include <deal.II/sundials/sunlinsol_wrapper.h>
+# include <deal.II/sundials/utilities.h>
# include <arkode/arkode_arkstep.h>
# include <sunlinsol/sunlinsol_spgmr.h>
namespace SUNDIALS
{
- namespace
- {
- /**
- * A function that calls the function object given by its first argument
- * with the set of arguments following at the end. If the call returns
- * regularly, the current function returns zero to indicate success. If the
- * call fails with an exception of type RecoverableUserCallbackError, then
- * the current function returns 1 to indicate that the called function
- * object thought the error it encountered is recoverable. If the call fails
- * with any other exception, then the current function returns with an error
- * code of -1. In each of the last two cases, the exception thrown by `f`
- * is captured and `eptr` is set to the exception. In case of success,
- * `eptr` is set to `nullptr`.
- */
- template <typename F, typename... Args>
- int
- call_and_possibly_capture_exception(const F & f,
- std::exception_ptr &eptr,
- Args &&...args)
- {
- // See whether there is already something in the exception pointer
- // variable. This can only happen if we had previously had
- // a recoverable exception, and the underlying library actually
- // did recover successfully. In that case, we can abandon the
- // exception previously thrown. If eptr contains anything other,
- // then we really don't know how that could have happened, and
- // should probably bail out:
- if (eptr)
- {
- try
- {
- std::rethrow_exception(eptr);
- }
- catch (const RecoverableUserCallbackError &)
- {
- // ok, ignore, but reset the pointer
- eptr = nullptr;
- }
- catch (...)
- {
- // uh oh:
- AssertThrow(false, ExcInternalError());
- }
- }
-
- // Call the function and if that succeeds, return zero:
- try
- {
- f(std::forward<Args>(args)...);
- eptr = nullptr;
- return 0;
- }
- // If the call failed with a recoverable error, then
- // ignore the exception for now (but store a pointer to it)
- // and return a positive return value (+1). If the underlying
- // implementation manages to recover
- catch (const RecoverableUserCallbackError &)
- {
- eptr = std::current_exception();
- return 1;
- }
- // For any other exception, capture the exception and
- // return -1:
- catch (const std::exception &)
- {
- eptr = std::current_exception();
- return -1;
- }
- }
- } // namespace
-
-
-
template <typename VectorType>
ARKode<VectorType>::ARKode(const AdditionalData &data)
: ARKode(data, MPI_COMM_SELF)
auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
auto *dst_yp = internal::unwrap_nvector<VectorType>(yp);
- return call_and_possibly_capture_exception(solver.explicit_function,
- solver.pending_exception,
- tt,
- *src_yy,
- *dst_yp);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.explicit_function,
+ solver.pending_exception,
+ tt,
+ *src_yy,
+ *dst_yp);
};
auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
auto *dst_yp = internal::unwrap_nvector<VectorType>(yp);
- return call_and_possibly_capture_exception(solver.implicit_function,
- solver.pending_exception,
- tt,
- *src_yy,
- *dst_yp);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.implicit_function,
+ solver.pending_exception,
+ tt,
+ *src_yy,
+ *dst_yp);
};
arkode_mem = ARKStepCreate(explicit_function ? explicit_function_callback :
auto *dst_Jv = internal::unwrap_nvector<VectorType>(Jv);
- return call_and_possibly_capture_exception(
+ return Utilities::call_and_possibly_capture_exception(
solver.jacobian_times_vector,
solver.pending_exception,
*src_v,
auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
- return call_and_possibly_capture_exception(
+ return Utilities::call_and_possibly_capture_exception(
solver.jacobian_times_setup,
solver.pending_exception,
t,
auto *dst_z = internal::unwrap_nvector<VectorType>(z);
- return call_and_possibly_capture_exception(
+ return Utilities::call_and_possibly_capture_exception(
solver.jacobian_preconditioner_solve,
solver.pending_exception,
t,
auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
- return call_and_possibly_capture_exception(
+ return Utilities::call_and_possibly_capture_exception(
solver.jacobian_preconditioner_setup,
solver.pending_exception,
t,
ARKode<VectorType> &solver =
*static_cast<ARKode<VectorType> *>(mtimes_data);
- return call_and_possibly_capture_exception(solver.mass_times_setup,
- solver.pending_exception,
- t);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.mass_times_setup, solver.pending_exception, t);
};
auto mass_matrix_times_vector_callback =
auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
auto *dst_Mv = internal::unwrap_nvector<VectorType>(Mv);
- return call_and_possibly_capture_exception(solver.mass_times_vector,
- solver.pending_exception,
- t,
- *src_v,
- *dst_Mv);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.mass_times_vector,
+ solver.pending_exception,
+ t,
+ *src_v,
+ *dst_Mv);
};
status = ARKStepSetMassTimes(arkode_mem,
ARKode<VectorType> &solver =
*static_cast<ARKode<VectorType> *>(user_data);
- return call_and_possibly_capture_exception(
+ return Utilities::call_and_possibly_capture_exception(
solver.mass_preconditioner_setup, solver.pending_exception, t);
};
auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
auto *dst_z = internal::unwrap_nvector<VectorType>(z);
- return call_and_possibly_capture_exception(
+ return Utilities::call_and_possibly_capture_exception(
solver.mass_preconditioner_solve,
solver.pending_exception,
t,
# endif
# include <deal.II/sundials/n_vector.h>
+# include <deal.II/sundials/utilities.h>
# include <iomanip>
# include <iostream>
namespace SUNDIALS
{
- namespace
- {
- /**
- * A function that calls the function object given by its first argument
- * with the set of arguments following at the end. If the call returns
- * regularly, the current function returns zero to indicate success. If the
- * call fails with an exception of type RecoverableUserCallbackError, then
- * the current function returns 1 to indicate that the called function
- * object thought the error it encountered is recoverable. If the call fails
- * with any other exception, then the current function returns with an error
- * code of -1. In each of the last two cases, the exception thrown by `f`
- * is captured and `eptr` is set to the exception. In case of success,
- * `eptr` is set to `nullptr`.
- */
- template <typename F, typename... Args>
- int
- call_and_possibly_capture_exception(const F & f,
- std::exception_ptr &eptr,
- Args &&...args)
- {
- // See whether there is already something in the exception pointer
- // variable. This can only happen if we had previously had
- // a recoverable exception, and the underlying library actually
- // did recover successfully. In that case, we can abandon the
- // exception previously thrown. If eptr contains anything other,
- // then we really don't know how that could have happened, and
- // should probably bail out:
- if (eptr)
- {
- try
- {
- std::rethrow_exception(eptr);
- }
- catch (const RecoverableUserCallbackError &)
- {
- // ok, ignore, but reset the pointer
- eptr = nullptr;
- }
- catch (...)
- {
- // uh oh:
- AssertThrow(false, ExcInternalError());
- }
- }
-
- // Call the function and if that succeeds, return zero:
- try
- {
- f(std::forward<Args>(args)...);
- eptr = nullptr;
- return 0;
- }
- // If the call failed with a recoverable error, then
- // ignore the exception for now (but store a pointer to it)
- // and return a positive return value (+1). If the underlying
- // implementation manages to recover
- catch (const RecoverableUserCallbackError &)
- {
- eptr = std::current_exception();
- return 1;
- }
- // For any other exception, capture the exception and
- // return -1:
- catch (const std::exception &)
- {
- eptr = std::current_exception();
- return -1;
- }
- }
- } // namespace
-
-
template <typename VectorType>
IDA<VectorType>::IDA(const AdditionalData &data)
: IDA(data, MPI_COMM_SELF)
auto *src_yp = internal::unwrap_nvector_const<VectorType>(yp);
auto *residual = internal::unwrap_nvector<VectorType>(rr);
- return call_and_possibly_capture_exception(solver.residual,
- solver.pending_exception,
- tt,
- *src_yy,
- *src_yp,
- *residual);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.residual,
+ solver.pending_exception,
+ tt,
+ *src_yy,
+ *src_yp,
+ *residual);
},
current_time,
yy,
auto *src_b = internal::unwrap_nvector_const<VectorType>(b);
auto *dst_x = internal::unwrap_nvector<VectorType>(x);
if (solver.solve_with_jacobian)
- return call_and_possibly_capture_exception(solver.solve_with_jacobian,
- solver.pending_exception,
- *src_b,
- *dst_x,
- tol);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.solve_with_jacobian,
+ solver.pending_exception,
+ *src_b,
+ *dst_x,
+ tol);
else if (solver.solve_jacobian_system)
- return call_and_possibly_capture_exception(solver.solve_jacobian_system,
- solver.pending_exception,
- *src_b,
- *dst_x);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.solve_jacobian_system,
+ solver.pending_exception,
+ *src_b,
+ *dst_x);
else
{
// We have already checked this outside, so we should never get here.
auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
auto *src_yp = internal::unwrap_nvector_const<VectorType>(yp);
- return call_and_possibly_capture_exception(solver.setup_jacobian,
- solver.pending_exception,
- tt,
- *src_yy,
- *src_yp,
- cj);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.setup_jacobian,
+ solver.pending_exception,
+ tt,
+ *src_yy,
+ *src_yp,
+ cj);
});
AssertIDA(status);
status = IDASetMaxOrd(ida_mem, data.maximum_order);
# endif
# include <deal.II/sundials/n_vector.h>
+# include <deal.II/sundials/utilities.h>
// Make sure we #include the SUNDIALS config file...
# include <sundials/sundials_config.h>
namespace SUNDIALS
{
- namespace
- {
- /**
- * A function that calls the function object given by its first argument
- * with the set of arguments following at the end. If the call returns
- * regularly, the current function returns zero to indicate success. If the
- * call fails with an exception of type RecoverableUserCallbackError, then
- * the current function returns 1 to indicate that the called function
- * object thought the error it encountered is recoverable. If the call fails
- * with any other exception, then the current function returns with an error
- * code of -1. In each of the last two cases, the exception thrown by `f`
- * is captured and `eptr` is set to the exception. In case of success,
- * `eptr` is set to `nullptr`.
- */
- template <typename F, typename... Args>
- int
- call_and_possibly_capture_exception(const F & f,
- std::exception_ptr &eptr,
- Args &&...args)
- {
- // See whether there is already something in the exception pointer
- // variable. This can only happen if we had previously had
- // a recoverable exception, and the underlying library actually
- // did recover successfully. In that case, we can abandon the
- // exception previously thrown. If eptr contains anything other,
- // then we really don't know how that could have happened, and
- // should probably bail out:
- if (eptr)
- {
- try
- {
- std::rethrow_exception(eptr);
- }
- catch (const RecoverableUserCallbackError &)
- {
- // ok, ignore, but reset the pointer
- eptr = nullptr;
- }
- catch (...)
- {
- // uh oh:
- AssertThrow(false, ExcInternalError());
- }
- }
-
- // Call the function and if that succeeds, return zero:
- try
- {
- f(std::forward<Args>(args)...);
- eptr = nullptr;
- return 0;
- }
- // If the call failed with a recoverable error, then
- // ignore the exception for now (but store a pointer to it)
- // and return a positive return value (+1). If the underlying
- // implementation manages to recover
- catch (const RecoverableUserCallbackError &)
- {
- eptr = std::current_exception();
- return 1;
- }
- // For any other exception, capture the exception and
- // return -1:
- catch (...)
- {
- eptr = std::current_exception();
- return -1;
- }
- }
- } // namespace
-
template <typename VectorType>
KINSOL<VectorType>::AdditionalData::AdditionalData(
const SolutionStrategy &strategy,
Assert(solver.iteration_function, ExcInternalError());
- const int err =
- call_and_possibly_capture_exception(solver.iteration_function,
- solver.pending_exception,
- *src_yy,
- *dst_FF);
+ const int err = Utilities::call_and_possibly_capture_exception(
+ solver.iteration_function,
+ solver.pending_exception,
+ *src_yy,
+ *dst_FF);
return err;
},
Assert(solver.residual, ExcInternalError());
- const int err = call_and_possibly_capture_exception(
+ const int err = Utilities::call_and_possibly_capture_exception(
solver.residual, solver.pending_exception, *src_yy, *dst_FF);
return err;
auto src_b = internal::unwrap_nvector_const<VectorType>(b);
auto dst_x = internal::unwrap_nvector<VectorType>(x);
- const int err =
- call_and_possibly_capture_exception(solver.solve_with_jacobian,
- solver.pending_exception,
- *src_b,
- *dst_x,
- tol);
+ const int err = Utilities::call_and_possibly_capture_exception(
+ solver.solve_with_jacobian,
+ solver.pending_exception,
+ *src_b,
+ *dst_x,
+ tol);
return err;
}
// src_ycur and src_fcur, and so we simply pass dummy vector in.
// These vectors will have zero lengths because we don't reinit
// them above.
- const int err = call_and_possibly_capture_exception(
+ const int err = Utilities::call_and_possibly_capture_exception(
solver.solve_jacobian_system,
solver.pending_exception,
*src_ycur,
auto fcur = internal::unwrap_nvector<VectorType>(f);
// Call the user-provided setup function with these arguments:
- return call_and_possibly_capture_exception(solver.setup_jacobian,
- solver.pending_exception,
- *ycur,
- *fcur);
+ return Utilities::call_and_possibly_capture_exception(
+ solver.setup_jacobian, solver.pending_exception, *ycur, *fcur);
});
AssertKINSOL(status);
}
projects / dealii.git / commitdiff