Reference documentation for deal.II version Git 0277a0d32d 2021-02-27 15:20:51 +0100
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
arkode.cc
Go to the documentation of this file.
1 //-----------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 //-----------------------------------------------------------
15 
16 
17 #include <deal.II/base/config.h>
18 
20 
21 #ifdef DEAL_II_WITH_SUNDIALS
22 
24 # include <deal.II/base/utilities.h>
25 
29 # ifdef DEAL_II_WITH_TRILINOS
32 # endif
33 # ifdef DEAL_II_WITH_PETSC
36 # endif
37 
40 
41 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
42 # include <arkode/arkode_impl.h>
43 # include <sundials/sundials_config.h>
44 # else
45 # include <arkode/arkode_arkstep.h>
46 # include <sunlinsol/sunlinsol_spgmr.h>
47 # include <sunnonlinsol/sunnonlinsol_fixedpoint.h>
48 # if DEAL_II_SUNDIALS_VERSION_LT(5, 0, 0)
50 # endif
51 # endif
52 
53 # include <iostream>
54 
55 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
56 // Make sure we know how to call sundials own ARKode() function
57 const auto &SundialsARKode = ARKode;
58 # endif
59 
61 
62 namespace SUNDIALS
63 {
64  using namespace internal;
65 
66  namespace
67  {
68  template <typename VectorType>
69  int
70  t_arkode_explicit_function(realtype tt,
71  N_Vector yy,
72  N_Vector yp,
73  void * user_data)
74  {
75  Assert(user_data != nullptr, ExcInternalError());
76  ARKode<VectorType> &solver =
77  *static_cast<ARKode<VectorType> *>(user_data);
78 
79  auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
80  auto *dst_yp = internal::unwrap_nvector<VectorType>(yp);
81 
82  return solver.explicit_function(tt, *src_yy, *dst_yp);
83  }
84 
85 
86 
87  template <typename VectorType>
88  int
89  t_arkode_implicit_function(realtype tt,
90  N_Vector yy,
91  N_Vector yp,
92  void * user_data)
93  {
94  Assert(user_data != nullptr, ExcInternalError());
95  ARKode<VectorType> &solver =
96  *static_cast<ARKode<VectorType> *>(user_data);
97 
98  auto *src_yy = internal::unwrap_nvector_const<VectorType>(yy);
99  auto *dst_yp = internal::unwrap_nvector<VectorType>(yp);
100 
101  return solver.implicit_function(tt, *src_yy, *dst_yp);
102  }
103 
104 
105 
106 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
107  template <typename VectorType>
108  int
109  t_arkode_setup_jacobian(ARKodeMem arkode_mem,
110  int convfail,
111  N_Vector ypred,
112  N_Vector fpred,
113  booleantype *jcurPtr,
114  N_Vector,
115  N_Vector,
116  N_Vector)
117  {
118  Assert(arkode_mem->ark_user_data != nullptr, ExcInternalError());
119  ARKode<VectorType> &solver =
120  *static_cast<ARKode<VectorType> *>(arkode_mem->ark_user_data);
121 
122  auto *src_ypred = internal::unwrap_nvector_const<VectorType>(ypred);
123  auto *src_fpred = internal::unwrap_nvector_const<VectorType>(fpred);
124  // avoid reinterpret_cast
125  bool jcurPtr_tmp = false;
126  int err = solver.setup_jacobian(convfail,
127  arkode_mem->ark_tn,
128  arkode_mem->ark_gamma,
129  *src_ypred,
130  *src_fpred,
131  jcurPtr_tmp);
132 # if DEAL_II_SUNDIALS_VERSION_GTE(2, 0, 0)
133  *jcurPtr = jcurPtr_tmp ? SUNTRUE : SUNFALSE;
134 # else
135  *jcurPtr = jcurPtr_tmp ? TRUE : FALSE;
136 # endif
137 
138  return err;
139  }
140 
141 
142 
143  template <typename VectorType>
144  int
145  t_arkode_solve_jacobian(ARKodeMem arkode_mem,
146  N_Vector b,
147 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
148  N_Vector,
149 # endif
150  N_Vector ycur,
151  N_Vector fcur)
152  {
153  Assert(arkode_mem->ark_user_data != nullptr, ExcInternalError());
154  ARKode<VectorType> &solver =
155  *static_cast<ARKode<VectorType> *>(arkode_mem->ark_user_data);
156 
157  auto *dst = internal::unwrap_nvector<VectorType>(b);
158  auto *src_ycur = internal::unwrap_nvector_const<VectorType>(ycur);
159  auto *src_fcur = internal::unwrap_nvector_const<VectorType>(fcur);
160 
161  // make a temporary copy to work on in the user call
162  VectorType src = *dst;
163 
164  int err = solver.solve_jacobian_system(arkode_mem->ark_tn,
165  arkode_mem->ark_gamma,
166  *src_ycur,
167  *src_fcur,
168  src,
169  *dst);
170 
171 
172  return err;
173  }
174 
175 
176 
177  template <typename VectorType>
178  int
179  t_arkode_setup_mass(ARKodeMem arkode_mem, N_Vector, N_Vector, N_Vector)
180  {
181  Assert(arkode_mem->ark_user_data != nullptr, ExcInternalError());
182  ARKode<VectorType> &solver =
183  *static_cast<ARKode<VectorType> *>(arkode_mem->ark_user_data);
184  int err = solver.setup_mass(arkode_mem->ark_tn);
185  return err;
186  }
187 
188 
189 
190  template <typename VectorType>
191  int
192  t_arkode_solve_mass(ARKodeMem arkode_mem,
193 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
194  N_Vector b,
195  N_Vector
196 # else
197  N_Vector b
198 # endif
199  )
200  {
201  Assert(arkode_mem->ark_user_data != nullptr, ExcInternalError());
202  ARKode<VectorType> &solver =
203  *static_cast<ARKode<VectorType> *>(arkode_mem->ark_user_data);
204 
205  auto *dst = internal::unwrap_nvector<VectorType>(b);
206 
207  // make a temporary copy to work on in the user call
208  VectorType src = *dst;
209 
210  int err = solver.solve_mass_system(src, *dst);
211 
212  return err;
213  }
214 # else
215 
216  template <typename VectorType>
217  int
218  t_arkode_jac_times_vec_function(N_Vector v,
219  N_Vector Jv,
220  realtype t,
221  N_Vector y,
222  N_Vector fy,
223  void * user_data,
224  N_Vector)
225  {
226  Assert(user_data != nullptr, ExcInternalError());
227  ARKode<VectorType> &solver =
228  *static_cast<ARKode<VectorType> *>(user_data);
230 
231  auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
232  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
233  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
234 
235  auto *dst_Jv = internal::unwrap_nvector<VectorType>(Jv);
236 
237  return solver.jacobian_times_vector(*src_v, *dst_Jv, t, *src_y, *src_fy);
238  }
239 
240 
241 
242  template <typename VectorType>
243  int
244  t_arkode_jac_times_setup_function(realtype t,
245  N_Vector y,
246  N_Vector fy,
247  void * user_data)
248  {
249  Assert(user_data != nullptr, ExcInternalError());
250  ARKode<VectorType> &solver =
251  *static_cast<ARKode<VectorType> *>(user_data);
252 
253  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
254  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
255 
256  return solver.jacobian_times_setup(t, *src_y, *src_fy);
257  }
258 
259 
260 
261  template <typename VectorType>
262  int
263  t_arkode_prec_solve_function(realtype t,
264  N_Vector y,
265  N_Vector fy,
266  N_Vector r,
267  N_Vector z,
268  realtype gamma,
269  realtype delta,
270  int lr,
271  void * user_data)
272  {
273  Assert(user_data != nullptr, ExcInternalError());
274  ARKode<VectorType> &solver =
275  *static_cast<ARKode<VectorType> *>(user_data);
276 
277  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
278  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
279  auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
280 
281  auto *dst_z = internal::unwrap_nvector<VectorType>(z);
282 
283  return solver.jacobian_preconditioner_solve(
284  t, *src_y, *src_fy, *src_r, *dst_z, gamma, delta, lr);
285  }
286 
287 
288 
289  template <typename VectorType>
290  int
291  t_arkode_prec_setup_function(realtype t,
292  N_Vector y,
293  N_Vector fy,
294  booleantype jok,
295  booleantype *jcurPtr,
296  realtype gamma,
297  void * user_data)
298  {
299  Assert(user_data != nullptr, ExcInternalError());
300  ARKode<VectorType> &solver =
301  *static_cast<ARKode<VectorType> *>(user_data);
302 
303  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
304  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
305 
306  return solver.jacobian_preconditioner_setup(
307  t, *src_y, *src_fy, jok, *jcurPtr, gamma);
308  }
309 
310 
311 
312  template <typename VectorType>
313  int
314  t_arkode_mass_times_vec_function(N_Vector v,
315  N_Vector Mv,
316  realtype t,
317  void * mtimes_data)
318  {
319  Assert(mtimes_data != nullptr, ExcInternalError());
320  ARKode<VectorType> &solver =
321  *static_cast<ARKode<VectorType> *>(mtimes_data);
322 
323  auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
324  auto *dst_Mv = internal::unwrap_nvector<VectorType>(Mv);
325 
326  return solver.mass_times_vector(t, *src_v, *dst_Mv);
327  }
328 
329 
330 
331  template <typename VectorType>
332  int
333  t_arkode_mass_times_setup_function(realtype t, void *mtimes_data)
334  {
335  Assert(mtimes_data != nullptr, ExcInternalError());
336  ARKode<VectorType> &solver =
337  *static_cast<ARKode<VectorType> *>(mtimes_data);
338 
339  return solver.mass_times_setup(t);
340  }
341 
342 
343 
344  template <typename VectorType>
345  int
346  t_arkode_mass_prec_solve_function(realtype t,
347  N_Vector r,
348  N_Vector z,
349  realtype delta,
350  int lr,
351  void * user_data)
352  {
353  Assert(user_data != nullptr, ExcInternalError());
354  ARKode<VectorType> &solver =
355  *static_cast<ARKode<VectorType> *>(user_data);
356 
357  auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
358  auto *dst_z = internal::unwrap_nvector<VectorType>(z);
359 
360  return solver.mass_preconditioner_solve(t, *src_r, *dst_z, delta, lr);
361  }
362 
363 
364 
365  template <typename VectorType>
366  int
367  t_arkode_mass_prec_setup_function(realtype t, void *user_data)
368  {
369  Assert(user_data != nullptr, ExcInternalError());
370  ARKode<VectorType> &solver =
371  *static_cast<ARKode<VectorType> *>(user_data);
372 
373  return solver.mass_preconditioner_setup(t);
374  }
375 # endif
376 
377  } // namespace
378 
379 
380  template <typename VectorType>
382  const MPI_Comm & mpi_comm)
383  : data(data)
384  , arkode_mem(nullptr)
385  , communicator(is_serial_vector<VectorType>::value ?
386  MPI_COMM_SELF :
387  Utilities::MPI::duplicate_communicator(mpi_comm))
388  {
390  }
391 
392  template <typename VectorType>
394  {
395  if (arkode_mem)
396 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
397  ARKodeFree(&arkode_mem);
398 # else
399  ARKStepFree(&arkode_mem);
400 # endif
401 
402 # ifdef DEAL_II_WITH_MPI
404  {
405  const int ierr = MPI_Comm_free(&communicator);
406  (void)ierr;
407  AssertNothrow(ierr == MPI_SUCCESS, ExcMPI(ierr));
408  }
409 # endif
410  }
411 
412 
413 
414  template <typename VectorType>
415  unsigned int
417  {
421  reset(time.get_current_time(), time.get_next_step_size(), solution);
422 
423  if (output_step)
424  output_step(time.get_current_time(), solution, time.get_step_number());
425 
426  while (time.is_at_end() == false)
427  {
429  double actual_next_time;
430 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
431  const auto status = SundialsARKode(
432  arkode_mem, time.get_next_time(), yy, &actual_next_time, ARK_NORMAL);
433 # else
434  const auto status = ARKStepEvolve(
435  arkode_mem, time.get_next_time(), yy, &actual_next_time, ARK_NORMAL);
436 # endif
437  (void)status;
438  AssertARKode(status);
439 
440  time.set_next_step_size(actual_next_time - time.get_current_time());
441  time.advance_time();
442 
443  while (solver_should_restart(time.get_current_time(), solution))
444  reset(time.get_current_time(),
445  time.get_previous_step_size(),
446  solution);
447 
448  if (output_step)
450  solution,
451  time.get_step_number());
452  }
453 
454  return time.get_step_number();
455  }
456 
457 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
458  template <typename VectorType>
459  void
460  ARKode<VectorType>::reset(const double current_time,
461  const double current_time_step,
462  VectorType & solution)
463  {
464  if (arkode_mem)
465  ARKodeFree(&arkode_mem);
466 
467  arkode_mem = ARKodeCreate();
468 
469  int status;
470  (void)status;
471 
473  ExcFunctionNotProvided("explicit_function || implicit_function"));
474 
475  // just a view on the memory in solution, all write operations on yy by
476  // ARKODE will automatically be mirrored to solution
477  yy = internal::make_nvector_view(solution);
478 
479  status = ARKodeInit(
480  arkode_mem,
481  explicit_function ? &t_arkode_explicit_function<VectorType> : nullptr,
482  implicit_function ? &t_arkode_implicit_function<VectorType> : nullptr,
483  current_time,
484  yy);
485  AssertARKode(status);
486 
488  {
489  abs_tolls = make_nvector_view(get_local_tolerances());
490  status =
491  ARKodeSVtolerances(arkode_mem, data.relative_tolerance, abs_tolls);
492  AssertARKode(status);
493  }
494  else
495  {
496  status = ARKodeSStolerances(arkode_mem,
499  AssertARKode(status);
500  }
501 
502  status = ARKodeSetInitStep(arkode_mem, current_time_step);
503  AssertARKode(status);
504 
505  status = ARKodeSetUserData(arkode_mem, this);
506  AssertARKode(status);
507 
508  status = ARKodeSetStopTime(arkode_mem, data.final_time);
509  AssertARKode(status);
510 
511  status =
512  ARKodeSetMaxNonlinIters(arkode_mem, data.maximum_non_linear_iterations);
513  AssertARKode(status);
514 
515  // Initialize solver
516  auto ARKode_mem = static_cast<ARKodeMem>(arkode_mem);
517 
519  {
520  status = ARKodeSetNewton(arkode_mem);
521  AssertARKode(status);
523  {
524  status = ARKodeSetLinear(
526  AssertARKode(status);
527  }
528 
529 
530  ARKode_mem->ark_lsolve = t_arkode_solve_jacobian<VectorType>;
531  if (setup_jacobian)
532  {
533  ARKode_mem->ark_lsetup = t_arkode_setup_jacobian<VectorType>;
534 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
535  ARKode_mem->ark_setupNonNull = true;
536 # endif
537  }
538  }
539  else
540  {
541  status =
542  ARKodeSetFixedPoint(arkode_mem, data.maximum_non_linear_iterations);
543  AssertARKode(status);
544  }
545 
546 
547  if (solve_mass_system)
548  {
549  ARKode_mem->ark_msolve = t_arkode_solve_mass<VectorType>;
550 
551  if (setup_mass)
552  {
553  ARKode_mem->ark_msetup = t_arkode_setup_mass<VectorType>;
554 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
555  ARKode_mem->ark_MassSetupNonNull = true;
556 # endif
557  }
558  }
559 
560  status = ARKodeSetOrder(arkode_mem, data.maximum_order);
561  AssertARKode(status);
562 
563  if (custom_setup)
565  }
566 
567 # else
568 
569  template <typename VectorType>
570  void
571  ARKode<VectorType>::reset(const double current_time,
572  const double current_time_step,
573  VectorType &solution)
574  {
575  if (arkode_mem)
576  ARKStepFree(&arkode_mem);
577 
578  int status;
579  (void)status;
580 
581  // just a view on the memory in solution, all write operations on yy by
582  // ARKODE will automatically be mirrored to solution
583  yy = internal::make_nvector_view(solution);
584 
586  ExcFunctionNotProvided("explicit_function || implicit_function"));
587 
588  arkode_mem = ARKStepCreate(
589  explicit_function ? &t_arkode_explicit_function<VectorType> : nullptr,
590  implicit_function ? &t_arkode_implicit_function<VectorType> : nullptr,
591  current_time,
592  yy);
593 
594  Assert(arkode_mem != nullptr, ExcInternalError());
595 
597  {
598  abs_tolls = internal::make_nvector_view(get_local_tolerances());
599  status =
600  ARKStepSVtolerances(arkode_mem, data.relative_tolerance, abs_tolls);
601  AssertARKode(status);
602  }
603  else
604  {
605  status = ARKStepSStolerances(arkode_mem,
608  AssertARKode(status);
609  }
610 
611  // for version 4.0.0 this call must be made before solver settings
612  status = ARKStepSetUserData(arkode_mem, this);
613  AssertARKode(status);
614 
615  setup_system_solver(solution);
616 
618 
619  status = ARKStepSetInitStep(arkode_mem, current_time_step);
620  AssertARKode(status);
621 
622  status = ARKStepSetStopTime(arkode_mem, data.final_time);
623  AssertARKode(status);
624 
625  status = ARKStepSetOrder(arkode_mem, data.maximum_order);
626  AssertARKode(status);
627 
628  if (custom_setup)
630  }
631 
632 
633 
634  template <typename VectorType>
635  void
637  {
638  // no point in setting up a solver when there is no linear system
639  if (!implicit_function)
640  return;
641 
642  int status;
643  (void)status;
644  // Initialize solver
645  // currently only iterative linear solvers are supported
647  {
648  SUNLinearSolver sun_linear_solver;
650  {
651  linear_solver =
652  std::make_unique<internal::LinearSolverWrapper<VectorType>>(
654  sun_linear_solver = *linear_solver;
655  }
656  else
657  {
658  // use default solver from SUNDIALS
659  // TODO give user options
660  sun_linear_solver =
661  SUNLinSol_SPGMR(yy,
662  PREC_NONE,
663  0 /*krylov subvectors, 0 uses default*/);
664  }
665  status = ARKStepSetLinearSolver(arkode_mem, sun_linear_solver, nullptr);
666  AssertARKode(status);
667  status =
668  ARKStepSetJacTimes(arkode_mem,
670  t_arkode_jac_times_setup_function<VectorType> :
671  nullptr,
672  t_arkode_jac_times_vec_function<VectorType>);
673  AssertARKode(status);
675  {
676  status = ARKStepSetPreconditioner(
677  arkode_mem,
679  t_arkode_prec_setup_function<VectorType> :
680  nullptr,
681  t_arkode_prec_solve_function<VectorType>);
682  AssertARKode(status);
683  }
685  {
686  status = ARKStepSetLinear(
688  AssertARKode(status);
689  }
690  }
691  else
692  {
693  auto y_template = internal::make_nvector_view(solution);
694 
695  SUNNonlinearSolver fixed_point_solver =
696  SUNNonlinSol_FixedPoint(y_template,
698 
699  status = ARKStepSetNonlinearSolver(arkode_mem, fixed_point_solver);
700  AssertARKode(status);
701  }
702 
703  status =
704  ARKStepSetMaxNonlinIters(arkode_mem, data.maximum_non_linear_iterations);
705  AssertARKode(status);
706  }
707 
708 
709 
710  template <typename VectorType>
711  void
713  {
714  int status;
715  (void)status;
716 
717  if (mass_times_vector)
718  {
719  SUNLinearSolver sun_mass_linear_solver;
720  if (solve_mass)
721  {
722  mass_solver =
723  std::make_unique<internal::LinearSolverWrapper<VectorType>>(
724  solve_mass);
725  sun_mass_linear_solver = *mass_solver;
726  }
727  else
728  {
729  sun_mass_linear_solver =
730  SUNLinSol_SPGMR(yy,
731  PREC_NONE,
732  0 /*krylov subvectors, 0 uses default*/);
733  }
734  booleantype mass_time_dependent =
735  data.mass_is_time_independent ? SUNFALSE : SUNTRUE;
736  status = ARKStepSetMassLinearSolver(arkode_mem,
737  sun_mass_linear_solver,
738  nullptr,
739  mass_time_dependent);
740  AssertARKode(status);
741  status =
742  ARKStepSetMassTimes(arkode_mem,
744  t_arkode_mass_times_setup_function<VectorType> :
745  nullptr,
746  t_arkode_mass_times_vec_function<VectorType>,
747  this);
748  AssertARKode(status);
749 
751  {
752  status = ARKStepSetMassPreconditioner(
753  arkode_mem,
755  t_arkode_mass_prec_setup_function<VectorType> :
756  nullptr,
757  t_arkode_mass_prec_solve_function<VectorType>);
758  AssertARKode(status);
759  }
760  }
761  }
762 # endif
763 
764 
765 
766  template <typename VectorType>
767  void
769  {
770  reinit_vector = [](VectorType &) {
771  AssertThrow(false, ExcFunctionNotProvided("reinit_vector"));
772  };
773 
774  solver_should_restart = [](const double, VectorType &) -> bool {
775  return false;
776  };
777  }
778 
779 
780 
781  template <typename VectorType>
782  void *
784  {
785  return arkode_mem;
786  }
787 
788  template class ARKode<Vector<double>>;
789  template class ARKode<BlockVector<double>>;
790 
793 
794 # ifdef DEAL_II_WITH_MPI
795 
796 # ifdef DEAL_II_WITH_TRILINOS
799 
800 # endif // DEAL_II_WITH_TRILINOS
801 
802 # ifdef DEAL_II_WITH_PETSC
803 # ifndef PETSC_USE_COMPLEX
804  template class ARKode<PETScWrappers::MPI::Vector>;
806 # endif // PETSC_USE_COMPLEX
807 # endif // DEAL_II_WITH_PETSC
808 
809 # endif // DEAL_II_WITH_MPI
810 
811 } // namespace SUNDIALS
812 
814 
815 #endif
#define AssertNothrow(cond, exc)
Definition: exceptions.h:1529
const auto & SundialsARKode
Definition: arkode.cc:57
void * arkode_mem
Definition: arkode.h:1263
#define DEAL_II_SUNDIALS_VERSION_LT(major, minor, patch)
Definition: config.h:305
void set_desired_next_step_size(const double time_step_size)
double get_next_step_size() const
void advance_time()
LinearSolveFunction< VectorType > solve_linearized_system
Definition: arkode.h:976
void * get_arkode_memory() const
Definition: arkode.cc:783
double get_current_time() const
internal::NVectorView< VectorType > abs_tolls
Definition: arkode.h:1273
std::function< int(const VectorType &rhs, VectorType &dst)> solve_mass_system
Definition: arkode.h:828
DEAL_II_DEPRECATED std::function< void(VectorType &)> reinit_vector
Definition: arkode.h:587
#define AssertThrow(cond, exc)
Definition: exceptions.h:1576
double get_previous_step_size() const
MPI_Comm communicator
Definition: arkode.h:1280
std::function< void(const double t, const VectorType &sol, const unsigned int step_number)> output_step
Definition: arkode.h:1162
void set_next_step_size(const double time_step_size)
std::function< VectorType &()> get_local_tolerances
Definition: arkode.h:1189
unsigned int maximum_non_linear_iterations
Definition: arkode.h:490
std::function< int(const double t)> setup_mass
Definition: arkode.h:807
std::function< bool(const double t, VectorType &sol)> solver_should_restart
Definition: arkode.h:1181
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > mass_solver
Definition: arkode.h:1289
ARKode(const AdditionalData &data=AdditionalData(), const MPI_Comm &mpi_comm=MPI_COMM_WORLD)
Definition: arkode.cc:381
std::function< int(const double t, const double gamma, const VectorType &ycur, const VectorType &fcur, const VectorType &rhs, VectorType &dst)> solve_jacobian_system
Definition: arkode.h:771
std::function< int(const double t)> mass_times_setup
Definition: arkode.h:885
AdditionalData data
Definition: arkode.h:1258
bool is_at_end() const
#define Assert(cond, exc)
Definition: exceptions.h:1466
std::function< int(double t, const VectorType &y, const VectorType &fy, const VectorType &r, VectorType &z, double gamma, double tol, int lr)> jacobian_preconditioner_solve
Definition: arkode.h:1035
void setup_system_solver(const VectorType &solution)
Definition: arkode.cc:636
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:394
std::function< int(const double t, const VectorType &y, VectorType &res)> implicit_function
Definition: arkode.h:626
std::function< int(realtype t, const VectorType &y, const VectorType &fy)> jacobian_times_setup
Definition: arkode.h:955
void setup_mass_solver()
Definition: arkode.cc:712
LinearSolveFunction< VectorType > solve_mass
Definition: arkode.h:993
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
std::function< int(double t, const VectorType &v, VectorType &Mv)> mass_times_vector
Definition: arkode.h:849
std::function< int(double t, const VectorType &y, const VectorType &fy, int jok, int &jcur, double gamma)> jacobian_preconditioner_setup
Definition: arkode.h:1085
MPI_Comm duplicate_communicator(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:159
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
std::function< int(const double t, const VectorType &y, VectorType &explicit_f)> explicit_function
Definition: arkode.h:607
unsigned int get_step_number() const
std::function< int(const VectorType &v, VectorType &Jv, double t, const VectorType &y, const VectorType &fy)> jacobian_times_vector
Definition: arkode.h:918
std::function< int(double t)> mass_preconditioner_setup
Definition: arkode.h:1142
std::function< void(void *arkode_mem)> custom_setup
Definition: arkode.h:1215
unsigned int solve_ode(VectorType &solution)
Definition: arkode.cc:416
std::function< int(const int convfail, const double t, const double gamma, const VectorType &ypred, const VectorType &fpred, bool &j_is_current)> setup_jacobian
Definition: arkode.h:719
#define AssertARKode(code)
Definition: arkode.h:60
internal::NVectorView< VectorType > yy
Definition: arkode.h:1268
std::function< int(double t, const VectorType &r, VectorType &z, double tol, int lr)> mass_preconditioner_solve
Definition: arkode.h:1116
long double gamma(const unsigned int n)
void set_functions_to_trigger_an_assert()
Definition: arkode.cc:768
double get_next_time() const
void reset(const double t, const double h, VectorType &y)
Definition: arkode.cc:460
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > linear_solver
Definition: arkode.h:1288
static ::ExceptionBase & ExcInternalError()