Reference documentation for deal.II version Git 1dc1051882 2021-04-22 23:57:03 +0200
\(\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);
229 
230  auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
231  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
232  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
233 
234  auto *dst_Jv = internal::unwrap_nvector<VectorType>(Jv);
235 
236  return solver.jacobian_times_vector(*src_v, *dst_Jv, t, *src_y, *src_fy);
237  }
238 
239 
240 
241  template <typename VectorType>
242  int
243  t_arkode_jac_times_setup_function(realtype t,
244  N_Vector y,
245  N_Vector fy,
246  void * user_data)
247  {
248  Assert(user_data != nullptr, ExcInternalError());
249  ARKode<VectorType> &solver =
250  *static_cast<ARKode<VectorType> *>(user_data);
251 
252  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
253  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
254 
255  return solver.jacobian_times_setup(t, *src_y, *src_fy);
256  }
257 
258 
259 
260  template <typename VectorType>
261  int
262  t_arkode_prec_solve_function(realtype t,
263  N_Vector y,
264  N_Vector fy,
265  N_Vector r,
266  N_Vector z,
267  realtype gamma,
268  realtype delta,
269  int lr,
270  void * user_data)
271  {
272  Assert(user_data != nullptr, ExcInternalError());
273  ARKode<VectorType> &solver =
274  *static_cast<ARKode<VectorType> *>(user_data);
275 
276  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
277  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
278  auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
279 
280  auto *dst_z = internal::unwrap_nvector<VectorType>(z);
281 
282  return solver.jacobian_preconditioner_solve(
283  t, *src_y, *src_fy, *src_r, *dst_z, gamma, delta, lr);
284  }
285 
286 
287 
288  template <typename VectorType>
289  int
290  t_arkode_prec_setup_function(realtype t,
291  N_Vector y,
292  N_Vector fy,
293  booleantype jok,
294  booleantype *jcurPtr,
295  realtype gamma,
296  void * user_data)
297  {
298  Assert(user_data != nullptr, ExcInternalError());
299  ARKode<VectorType> &solver =
300  *static_cast<ARKode<VectorType> *>(user_data);
301 
302  auto *src_y = internal::unwrap_nvector_const<VectorType>(y);
303  auto *src_fy = internal::unwrap_nvector_const<VectorType>(fy);
304 
305  return solver.jacobian_preconditioner_setup(
306  t, *src_y, *src_fy, jok, *jcurPtr, gamma);
307  }
308 
309 
310 
311  template <typename VectorType>
312  int
313  t_arkode_mass_times_vec_function(N_Vector v,
314  N_Vector Mv,
315  realtype t,
316  void * mtimes_data)
317  {
318  Assert(mtimes_data != nullptr, ExcInternalError());
319  ARKode<VectorType> &solver =
320  *static_cast<ARKode<VectorType> *>(mtimes_data);
321 
322  auto *src_v = internal::unwrap_nvector_const<VectorType>(v);
323  auto *dst_Mv = internal::unwrap_nvector<VectorType>(Mv);
324 
325  return solver.mass_times_vector(t, *src_v, *dst_Mv);
326  }
327 
328 
329 
330  template <typename VectorType>
331  int
332  t_arkode_mass_times_setup_function(realtype t, void *mtimes_data)
333  {
334  Assert(mtimes_data != nullptr, ExcInternalError());
335  ARKode<VectorType> &solver =
336  *static_cast<ARKode<VectorType> *>(mtimes_data);
337 
338  return solver.mass_times_setup(t);
339  }
340 
341 
342 
343  template <typename VectorType>
344  int
345  t_arkode_mass_prec_solve_function(realtype t,
346  N_Vector r,
347  N_Vector z,
348  realtype delta,
349  int lr,
350  void * user_data)
351  {
352  Assert(user_data != nullptr, ExcInternalError());
353  ARKode<VectorType> &solver =
354  *static_cast<ARKode<VectorType> *>(user_data);
355 
356  auto *src_r = internal::unwrap_nvector_const<VectorType>(r);
357  auto *dst_z = internal::unwrap_nvector<VectorType>(z);
358 
359  return solver.mass_preconditioner_solve(t, *src_r, *dst_z, delta, lr);
360  }
361 
362 
363 
364  template <typename VectorType>
365  int
366  t_arkode_mass_prec_setup_function(realtype t, void *user_data)
367  {
368  Assert(user_data != nullptr, ExcInternalError());
369  ARKode<VectorType> &solver =
370  *static_cast<ARKode<VectorType> *>(user_data);
371 
372  return solver.mass_preconditioner_setup(t);
373  }
374 # endif
375 
376  } // namespace
377 
378 
379  template <typename VectorType>
381  const MPI_Comm & mpi_comm)
382  : data(data)
383  , arkode_mem(nullptr)
384  , communicator(is_serial_vector<VectorType>::value ?
385  MPI_COMM_SELF :
386  Utilities::MPI::duplicate_communicator(mpi_comm))
387  , last_end_time(data.initial_time)
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 
422  return do_evolve_time(solution, time, /* force_solver_restart = */ true);
423  }
424 
425 
426 
427  template <typename VectorType>
428  unsigned int
430  const double intermediate_time,
431  const bool reset_solver)
432  {
433  AssertThrow(
434  intermediate_time > last_end_time,
435  ::ExcMessage(
436  "The requested intermediate time is smaller than the last requested "
437  "intermediate time."));
438 
439  const bool do_reset = reset_solver || arkode_mem == nullptr;
440  DiscreteTime time(last_end_time, intermediate_time, data.initial_step_size);
441  return do_evolve_time(solution, time, do_reset);
442  }
443 
444 
445 
446  template <typename VectorType>
447  int
449  DiscreteTime &time,
450  const bool do_reset)
451  {
452  if (do_reset)
453  {
454  reset(time.get_current_time(), time.get_next_step_size(), solution);
455  if (output_step)
457  solution,
458  time.get_step_number());
459  }
460 
461  auto solution_nvector = internal::make_nvector_view(solution);
462 
463  while (!time.is_at_end())
464  {
466  double actual_next_time;
467 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
468  const auto status = SundialsARKode(arkode_mem,
469  time.get_next_time(),
470  solution_nvector,
471  &actual_next_time,
472  ARK_NORMAL);
473 # else
474  const auto status = ARKStepEvolve(arkode_mem,
475  time.get_next_time(),
476  solution_nvector,
477  &actual_next_time,
478  ARK_NORMAL);
479 # endif
480  (void)status;
481  AssertARKode(status);
482 
483  time.set_next_step_size(actual_next_time - time.get_current_time());
484  time.advance_time();
485 
486  while (solver_should_restart(time.get_current_time(), solution))
487  reset(time.get_current_time(),
488  time.get_previous_step_size(),
489  solution);
490 
491  if (output_step)
493  solution,
494  time.get_step_number());
495  }
497  return time.get_step_number();
498  }
499 
500 
501 
502 # if DEAL_II_SUNDIALS_VERSION_LT(4, 0, 0)
503  template <typename VectorType>
504  void
505  ARKode<VectorType>::reset(const double current_time,
506  const double current_time_step,
507  const VectorType &solution)
508  {
509  last_end_time = current_time;
510  if (arkode_mem)
511  ARKodeFree(&arkode_mem);
512 
513  arkode_mem = ARKodeCreate();
514 
515  int status;
516  (void)status;
517 
519  ExcFunctionNotProvided("explicit_function || implicit_function"));
520 
521  // just a view on the memory in solution, all write operations on yy by
522  // ARKODE will automatically be mirrored to solution
523  auto initial_condition_nvector = internal::make_nvector_view(solution);
524 
525  status = ARKodeInit(
526  arkode_mem,
527  explicit_function ? &t_arkode_explicit_function<VectorType> : nullptr,
528  implicit_function ? &t_arkode_implicit_function<VectorType> : nullptr,
529  current_time,
530  initial_condition_nvector);
531  AssertARKode(status);
532 
534  {
535  const auto abs_tols = make_nvector_view(get_local_tolerances());
536  status =
537  ARKodeSVtolerances(arkode_mem, data.relative_tolerance, abs_tols);
538  AssertARKode(status);
539  }
540  else
541  {
542  status = ARKodeSStolerances(arkode_mem,
545  AssertARKode(status);
546  }
547 
548  status = ARKodeSetInitStep(arkode_mem, current_time_step);
549  AssertARKode(status);
550 
551  status = ARKodeSetUserData(arkode_mem, this);
552  AssertARKode(status);
553 
554  status = ARKodeSetStopTime(arkode_mem, data.final_time);
555  AssertARKode(status);
556 
557  status =
558  ARKodeSetMaxNonlinIters(arkode_mem, data.maximum_non_linear_iterations);
559  AssertARKode(status);
560 
561  // Initialize solver
562  auto ARKode_mem = static_cast<ARKodeMem>(arkode_mem);
563 
565  {
566  status = ARKodeSetNewton(arkode_mem);
567  AssertARKode(status);
569  {
570  status = ARKodeSetLinear(
572  AssertARKode(status);
573  }
574 
575 
576  ARKode_mem->ark_lsolve = t_arkode_solve_jacobian<VectorType>;
577  if (setup_jacobian)
578  {
579  ARKode_mem->ark_lsetup = t_arkode_setup_jacobian<VectorType>;
580 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
581  ARKode_mem->ark_setupNonNull = true;
582 # endif
583  }
584  }
585  else
586  {
587  status =
588  ARKodeSetFixedPoint(arkode_mem, data.maximum_non_linear_iterations);
589  AssertARKode(status);
590  }
591 
592 
593  if (solve_mass_system)
594  {
595  ARKode_mem->ark_msolve = t_arkode_solve_mass<VectorType>;
596 
597  if (setup_mass)
598  {
599  ARKode_mem->ark_msetup = t_arkode_setup_mass<VectorType>;
600 # if DEAL_II_SUNDIALS_VERSION_LT(3, 0, 0)
601  ARKode_mem->ark_MassSetupNonNull = true;
602 # endif
603  }
604  }
605 
606  status = ARKodeSetOrder(arkode_mem, data.maximum_order);
607  AssertARKode(status);
608 
609  if (custom_setup)
611  }
612 
613 # else
614 
615  template <typename VectorType>
616  void
617  ARKode<VectorType>::reset(const double current_time,
618  const double current_time_step,
619  const VectorType &solution)
620  {
621  last_end_time = current_time;
622  if (arkode_mem)
623  ARKStepFree(&arkode_mem);
624 
625  int status;
626  (void)status;
627 
628  // just a view on the memory in solution, all write operations on yy by
629  // ARKODE will automatically be mirrored to solution
630  auto initial_condition_nvector = internal::make_nvector_view(solution);
631 
633  ExcFunctionNotProvided("explicit_function || implicit_function"));
634 
635  arkode_mem = ARKStepCreate(
636  explicit_function ? &t_arkode_explicit_function<VectorType> : nullptr,
637  implicit_function ? &t_arkode_implicit_function<VectorType> : nullptr,
638  current_time,
639  initial_condition_nvector);
640 
641  Assert(arkode_mem != nullptr, ExcInternalError());
642 
644  {
645  const auto abs_tols =
646  internal::make_nvector_view(get_local_tolerances());
647  status =
648  ARKStepSVtolerances(arkode_mem, data.relative_tolerance, abs_tols);
649  AssertARKode(status);
650  }
651  else
652  {
653  status = ARKStepSStolerances(arkode_mem,
656  AssertARKode(status);
657  }
658 
659  // for version 4.0.0 this call must be made before solver settings
660  status = ARKStepSetUserData(arkode_mem, this);
661  AssertARKode(status);
662 
663  setup_system_solver(solution);
664 
665  setup_mass_solver(solution);
666 
667  status = ARKStepSetInitStep(arkode_mem, current_time_step);
668  AssertARKode(status);
669 
670  status = ARKStepSetStopTime(arkode_mem, data.final_time);
671  AssertARKode(status);
672 
673  status = ARKStepSetOrder(arkode_mem, data.maximum_order);
674  AssertARKode(status);
675 
676  if (custom_setup)
678  }
679 
680 
681 
682  template <typename VectorType>
683  void
685  {
686  // no point in setting up a solver when there is no linear system
687  if (!implicit_function)
688  return;
689 
690  int status;
691  (void)status;
692  // Initialize solver
693  // currently only iterative linear solvers are supported
695  {
696  SUNLinearSolver sun_linear_solver;
698  {
699  linear_solver =
700  std::make_unique<internal::LinearSolverWrapper<VectorType>>(
702  sun_linear_solver = *linear_solver;
703  }
704  else
705  {
706  // use default solver from SUNDIALS
707  // TODO give user options
708  auto y_template = internal::make_nvector_view(solution);
709  sun_linear_solver =
710  SUNLinSol_SPGMR(y_template,
711  PREC_NONE,
712  0 /*krylov subvectors, 0 uses default*/);
713  }
714  status = ARKStepSetLinearSolver(arkode_mem, sun_linear_solver, nullptr);
715  AssertARKode(status);
716  status =
717  ARKStepSetJacTimes(arkode_mem,
719  t_arkode_jac_times_setup_function<VectorType> :
720  nullptr,
721  t_arkode_jac_times_vec_function<VectorType>);
722  AssertARKode(status);
724  {
725  status = ARKStepSetPreconditioner(
726  arkode_mem,
728  t_arkode_prec_setup_function<VectorType> :
729  nullptr,
730  t_arkode_prec_solve_function<VectorType>);
731  AssertARKode(status);
732  }
734  {
735  status = ARKStepSetLinear(
737  AssertARKode(status);
738  }
739  }
740  else
741  {
742  auto y_template = internal::make_nvector_view(solution);
743 
744  SUNNonlinearSolver fixed_point_solver =
745  SUNNonlinSol_FixedPoint(y_template,
747 
748  status = ARKStepSetNonlinearSolver(arkode_mem, fixed_point_solver);
749  AssertARKode(status);
750  }
751 
752  status =
753  ARKStepSetMaxNonlinIters(arkode_mem, data.maximum_non_linear_iterations);
754  AssertARKode(status);
755  }
756 
757 
758 
759  template <typename VectorType>
760  void
762  {
763  int status;
764  (void)status;
765 
766  if (mass_times_vector)
767  {
768  SUNLinearSolver sun_mass_linear_solver;
769  if (solve_mass)
770  {
771  mass_solver =
772  std::make_unique<internal::LinearSolverWrapper<VectorType>>(
773  solve_mass);
774  sun_mass_linear_solver = *mass_solver;
775  }
776  else
777  {
778  auto y_template = internal::make_nvector_view(solution);
779  sun_mass_linear_solver =
780  SUNLinSol_SPGMR(y_template,
781  PREC_NONE,
782  0 /*krylov subvectors, 0 uses default*/);
783  }
784  booleantype mass_time_dependent =
785  data.mass_is_time_independent ? SUNFALSE : SUNTRUE;
786  status = ARKStepSetMassLinearSolver(arkode_mem,
787  sun_mass_linear_solver,
788  nullptr,
789  mass_time_dependent);
790  AssertARKode(status);
791  status =
792  ARKStepSetMassTimes(arkode_mem,
794  t_arkode_mass_times_setup_function<VectorType> :
795  nullptr,
796  t_arkode_mass_times_vec_function<VectorType>,
797  this);
798  AssertARKode(status);
799 
801  {
802  status = ARKStepSetMassPreconditioner(
803  arkode_mem,
805  t_arkode_mass_prec_setup_function<VectorType> :
806  nullptr,
807  t_arkode_mass_prec_solve_function<VectorType>);
808  AssertARKode(status);
809  }
810  }
811  }
812 # endif
813 
814 
815 
816  template <typename VectorType>
817  void
819  {
820  reinit_vector = [](VectorType &) {
821  AssertThrow(false, ExcFunctionNotProvided("reinit_vector"));
822  };
823 
824  solver_should_restart = [](const double, VectorType &) -> bool {
825  return false;
826  };
827  }
828 
829 
830 
831  template <typename VectorType>
832  void *
834  {
835  return arkode_mem;
836  }
837 
838  template class ARKode<Vector<double>>;
839  template class ARKode<BlockVector<double>>;
840 
843 
844 # ifdef DEAL_II_WITH_MPI
845 
846 # ifdef DEAL_II_WITH_TRILINOS
849 
850 # endif // DEAL_II_WITH_TRILINOS
851 
852 # ifdef DEAL_II_WITH_PETSC
853 # ifndef PETSC_USE_COMPLEX
854  template class ARKode<PETScWrappers::MPI::Vector>;
856 # endif // PETSC_USE_COMPLEX
857 # endif // DEAL_II_WITH_PETSC
858 
859 # endif // DEAL_II_WITH_MPI
860 
861 } // namespace SUNDIALS
862 
864 
865 #endif
void setup_mass_solver(const VectorType &solution)
Definition: arkode.cc:761
void reset(const double t, const double h, const VectorType &y)
Definition: arkode.cc:505
#define AssertNothrow(cond, exc)
Definition: exceptions.h:1536
const auto & SundialsARKode
Definition: arkode.cc:57
void * arkode_mem
Definition: arkode.h:1309
#define DEAL_II_SUNDIALS_VERSION_LT(major, minor, patch)
Definition: config.h:307
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:1012
double last_end_time
Definition: arkode.h:1321
void * get_arkode_memory() const
Definition: arkode.cc:833
double get_current_time() const
std::function< int(const VectorType &rhs, VectorType &dst)> solve_mass_system
Definition: arkode.h:864
DEAL_II_DEPRECATED std::function< void(VectorType &)> reinit_vector
Definition: arkode.h:623
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
double get_previous_step_size() const
MPI_Comm communicator
Definition: arkode.h:1316
std::function< void(const double t, const VectorType &sol, const unsigned int step_number)> output_step
Definition: arkode.h:1198
void set_next_step_size(const double time_step_size)
std::function< VectorType &()> get_local_tolerances
Definition: arkode.h:1225
unsigned int maximum_non_linear_iterations
Definition: arkode.h:492
std::function< int(const double t)> setup_mass
Definition: arkode.h:843
std::function< bool(const double t, VectorType &sol)> solver_should_restart
Definition: arkode.h:1217
static ::ExceptionBase & ExcMessage(std::string arg1)
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > mass_solver
Definition: arkode.h:1325
ARKode(const AdditionalData &data=AdditionalData(), const MPI_Comm &mpi_comm=MPI_COMM_WORLD)
Definition: arkode.cc:380
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:807
std::function< int(const double t)> mass_times_setup
Definition: arkode.h:921
AdditionalData data
Definition: arkode.h:1304
bool is_at_end() const
#define Assert(cond, exc)
Definition: exceptions.h:1473
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:1071
void setup_system_solver(const VectorType &solution)
Definition: arkode.cc:684
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
std::function< int(const double t, const VectorType &y, VectorType &res)> implicit_function
Definition: arkode.h:662
std::function< int(realtype t, const VectorType &y, const VectorType &fy)> jacobian_times_setup
Definition: arkode.h:991
LinearSolveFunction< VectorType > solve_mass
Definition: arkode.h:1029
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:885
std::function< int(double t, const VectorType &y, const VectorType &fy, int jok, int &jcur, double gamma)> jacobian_preconditioner_setup
Definition: arkode.h:1121
MPI_Comm duplicate_communicator(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:160
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
std::function< int(const double t, const VectorType &y, VectorType &explicit_f)> explicit_function
Definition: arkode.h:643
unsigned int get_step_number() const
int do_evolve_time(VectorType &solution, ::DiscreteTime &time, const bool do_reset)
Definition: arkode.cc:448
std::function< int(const VectorType &v, VectorType &Jv, double t, const VectorType &y, const VectorType &fy)> jacobian_times_vector
Definition: arkode.h:954
std::function< int(double t)> mass_preconditioner_setup
Definition: arkode.h:1178
std::function< void(void *arkode_mem)> custom_setup
Definition: arkode.h:1251
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:755
#define AssertARKode(code)
Definition: arkode.h:62
std::function< int(double t, const VectorType &r, VectorType &z, double tol, int lr)> mass_preconditioner_solve
Definition: arkode.h:1152
long double gamma(const unsigned int n)
void set_functions_to_trigger_an_assert()
Definition: arkode.cc:818
double get_next_time() const
unsigned int solve_ode_incrementally(VectorType &solution, const double intermediate_time, const bool reset_solver=false)
Definition: arkode.cc:429
std::unique_ptr< internal::LinearSolverWrapper< VectorType > > linear_solver
Definition: arkode.h:1324
static ::ExceptionBase & ExcInternalError()