Reference documentation for deal.II version 9.6.0
\(\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\}}\)
Loading...
Searching...
No Matches
nonlinear.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#include <deal.II/base/config.h>
16
18
20
22
24
25
27
74template <typename VectorType = Vector<double>>
76{
77public:
85 {
86 public:
107
109 {
123 /*
124 * Use the PETSc SNES solver.
125 */
127 };
128
143 AdditionalData(const SolverType &solver_type = automatic,
144 const SolutionStrategy &strategy = linesearch,
145 const unsigned int maximum_non_linear_iterations = 200,
146 const double function_tolerance = 1e-8,
147 const double relative_tolerance = 1e-5,
148 const double step_tolerance = 0.0,
149 const unsigned int anderson_subspace_size = 0);
150
157
165
170
178
185 const double relative_tolerance;
186
194
202 };
203
208
217
226
231 void
232 select(const typename AdditionalData::SolverType &type);
233
237 void
239
244#ifdef DEAL_II_TRILINOS_WITH_NOX
245 void
246 set_data(
249 const Teuchos::RCP<Teuchos::ParameterList> &parameters =
250 Teuchos::rcp(new Teuchos::ParameterList));
251#endif
252
257#ifdef DEAL_II_WITH_SUNDIALS
258 void
261#endif
262
267#ifdef DEAL_II_WITH_PETSC
268 void
270#endif
271
278 void
279 solve(VectorType &initial_guess_and_solution);
280
287 std::function<void(VectorType &)> reinit_vector;
288
302 std::function<void(const VectorType &src, VectorType &dst)> residual;
303
330 std::function<void(const VectorType &current_u)> setup_jacobian;
331
353 std::function<
354 void(const VectorType &rhs, VectorType &dst, const double tolerance)>
356
357private:
362
363private:
368
372#ifdef DEAL_II_WITH_SUNDIALS
374#endif
375
379#ifdef DEAL_II_TRILINOS_WITH_NOX
382 Teuchos::RCP<Teuchos::ParameterList> parameters_nox =
383 Teuchos::rcp(new Teuchos::ParameterList);
384#endif
385
389#ifdef DEAL_II_WITH_PETSC
391#endif
392
399 void
400 solve_with_petsc(VectorType &initial_guess_and_solution);
401};
402
403
404
405template <typename VectorType>
406void
408 const AdditionalData &additional_data)
409{
410#ifdef DEAL_II_WITH_SUNDIALS
411 // These if statements pass on the strategy to the other nonlinear solvers
412 if (additional_data.strategy ==
414 additional_data_kinsol.strategy =
416 else if (additional_data.strategy ==
418 additional_data_kinsol.strategy =
420 else if (additional_data.strategy ==
422 additional_data_kinsol.strategy =
424
425 // Setting data points in the KINSOL class from the NonlinearSolverSelector
426 // class
427 additional_data_kinsol.maximum_non_linear_iterations =
428 additional_data.maximum_non_linear_iterations;
429 additional_data_kinsol.function_tolerance =
430 additional_data.function_tolerance;
431 additional_data_kinsol.step_tolerance = additional_data.step_tolerance;
432 additional_data_kinsol.anderson_subspace_size =
433 additional_data.anderson_subspace_size;
434#endif
435
436// Do the same thing we did above but with NOX
437#ifdef DEAL_II_TRILINOS_WITH_NOX
438 parameters_nox->set("Nonlinear Solver", "Line Search Based");
439 Teuchos::ParameterList &Line_Search = parameters_nox->sublist("Line Search");
440 Line_Search.set("Method", "Full Step");
441
442 additional_data_nox.max_iter = additional_data.maximum_non_linear_iterations;
443 additional_data_nox.abs_tol = additional_data.function_tolerance;
444 additional_data_nox.rel_tol = additional_data.relative_tolerance;
445#endif
446
447#ifdef DEAL_II_WITH_PETSC
448 additional_data_petsc_snes.options_prefix = "";
449
450 if (additional_data.anderson_subspace_size > 0 &&
451 additional_data.strategy ==
453 {
454 additional_data_petsc_snes.snes_type = "anderson";
455 // TODO additional_data.anderson_subspace_size;
456 }
457 else if (additional_data.strategy ==
459 {
460 additional_data_petsc_snes.snes_type = "newtonls";
461 additional_data_petsc_snes.snes_linesearch_type = "bt";
462 }
463 else if (additional_data.strategy ==
465 {
466 additional_data_petsc_snes.snes_linesearch_type = "newtonls";
467 additional_data_petsc_snes.snes_linesearch_type = "basic";
468 }
469 else if (additional_data.strategy ==
471 additional_data_petsc_snes.snes_type = "nrichardson";
472
473 additional_data_petsc_snes.absolute_tolerance =
474 additional_data.function_tolerance;
475 additional_data_petsc_snes.relative_tolerance =
476 additional_data.relative_tolerance;
477 additional_data_petsc_snes.step_tolerance = additional_data.step_tolerance;
478 additional_data_petsc_snes.maximum_non_linear_iterations =
479 additional_data.maximum_non_linear_iterations;
480 additional_data_petsc_snes.max_n_function_evaluations = -1;
481
482#endif
483}
484
485
486
487template <typename VectorType>
489 : mpi_communicator(MPI_COMM_SELF)
490{}
491
492
493
494template <typename VectorType>
496 const AdditionalData &additional_data)
497 : additional_data(additional_data)
498 , mpi_communicator(MPI_COMM_SELF)
499{
501}
502
503
504
505template <typename VectorType>
507 const AdditionalData &additional_data,
508 const MPI_Comm &mpi_communicator)
509 : additional_data(additional_data)
510 , mpi_communicator(mpi_communicator)
511{
513}
514
515
516
517template <typename VectorType>
518void
520 const typename AdditionalData::SolverType &type)
521{
522 additional_data.solver_type = type;
523}
524
525
526
527template <typename VectorType>
529 const SolverType &solver_type,
530 const SolutionStrategy &strategy,
531 const unsigned int maximum_non_linear_iterations,
532 const double function_tolerance,
533 const double relative_tolerance,
534 const double step_tolerance,
535 const unsigned int anderson_subspace_size)
536 : solver_type(solver_type)
537 , strategy(strategy)
538 , maximum_non_linear_iterations(maximum_non_linear_iterations)
539 , function_tolerance(function_tolerance)
540 , relative_tolerance(relative_tolerance)
541 , step_tolerance(step_tolerance)
542 , anderson_subspace_size(anderson_subspace_size)
543{}
544
545
546
547#ifdef DEAL_II_TRILINOS_WITH_NOX
548template <typename VectorType>
549void
553 const Teuchos::RCP<Teuchos::ParameterList> &parameters)
554{
556 parameters_nox = parameters;
557}
558#endif
559
560
561
562#ifdef DEAL_II_WITH_SUNDIALS
563template <typename VectorType>
564void
570#endif
571
572
573
574template <typename VectorType>
575void
577 VectorType & /*initial_guess_and_solution*/)
578{
579 AssertThrow(false,
580 ExcMessage("PETSc SNES requires you to use PETSc vectors."));
581}
582
583
584
585#ifdef DEAL_II_WITH_PETSC
586template <>
587void
589 PETScWrappers::MPI::Vector &initial_guess_and_solution)
590{
593
594 nonlinear_solver.residual = residual;
595
596 nonlinear_solver.setup_jacobian = setup_jacobian;
597
598 nonlinear_solver.solve_with_jacobian =
599 [&](const PETScWrappers::MPI::Vector &src,
601 // PETSc does not gives a tolerance, so we have to choose something
602 // reasonable to provide to the user:
603 const double tolerance = 1e-6;
604 solve_with_jacobian(src, dst, tolerance);
605 };
606
607 nonlinear_solver.solve(initial_guess_and_solution);
608}
609#endif
610
611
612
613template <typename VectorType>
614void
616 VectorType &initial_guess_and_solution)
617{
618 // The "auto" solver_type will default to kinsol, however if KINSOL is not
619 // available then we will use NOX.
621 {
622#ifdef DEAL_II_WITH_PETSC
624#endif
625#ifdef DEAL_II_TRILINOS_WITH_NOX
627#endif
628#ifdef DEAL_II_WITH_SUNDIALS
630#endif
631
632 // If "auto" is still the solver type we cannot solve the problem
634 AssertThrow(false, ExcMessage("No valid solver type."));
635 }
636
638 {
639#ifdef DEAL_II_WITH_SUNDIALS
642
643 nonlinear_solver.reinit_vector = reinit_vector;
644 nonlinear_solver.residual = residual;
645
646 // We cannot simply set these two functions equal to each other
647 // because they have a different number of inputs.
648 nonlinear_solver.setup_jacobian = [&](const VectorType &current_u,
649 const VectorType & /*current_f*/) {
651 };
652
653 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
654
655 nonlinear_solver.solve(initial_guess_and_solution);
656#else
658 false, ExcMessage("You do not have SUNDIALS configured with deal.II!"));
659#endif
660 }
662 {
663#ifdef DEAL_II_TRILINOS_WITH_NOX
664
667
668 // Do the same thing for NOX that we did with KINSOL.
669 nonlinear_solver.residual = residual;
670
671 // setup_jacobian for NOX has the same number of arguments for the same
672 // function in NonlinearSolverSelector.
673 nonlinear_solver.setup_jacobian = setup_jacobian;
674
675 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
676
677 nonlinear_solver.solve(initial_guess_and_solution);
678#else
680 false, ExcMessage("You do not have Trilinos configured with deal.II"));
681#endif
682 }
683 else if (additional_data.solver_type ==
685 {
686 // Forward to internal function specializations, which throws for
687 // non-supported vector types:
688 solve_with_petsc(initial_guess_and_solution);
689 }
690 else
691 {
692 const std::string solvers =
693#ifdef DEAL_II_WITH_SUNDIALS
694 "kinsol\n"
695#endif
696#ifdef DEAL_II_TRILINOS_WITH_NOX
697 "NOX\n"
698#endif
699#ifdef DEAL_II_WITH_PETSC
700 "SNES\n"
701#endif
702 ;
703
704 AssertThrow(false,
706 "Invalid nonlinear solver specified. "
707 "The solvers available in your installation are:\n" +
708 solvers));
709 }
710}
711
AdditionalData(const SolverType &solver_type=automatic, const SolutionStrategy &strategy=linesearch, const unsigned int maximum_non_linear_iterations=200, const double function_tolerance=1e-8, const double relative_tolerance=1e-5, const double step_tolerance=0.0, const unsigned int anderson_subspace_size=0)
Definition nonlinear.h:528
TrilinosWrappers::NOXSolver< VectorType >::AdditionalData additional_data_nox
Definition nonlinear.h:381
void set_data(const AdditionalData &additional_data)
Definition nonlinear.h:407
void solve_with_petsc(VectorType &initial_guess_and_solution)
Definition nonlinear.h:576
void set_data(const typename PETScWrappers::NonlinearSolverData &additional_data)
SUNDIALS::KINSOL< VectorType >::AdditionalData additional_data_kinsol
Definition nonlinear.h:373
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition nonlinear.h:302
std::function< void(VectorType &)> reinit_vector
Definition nonlinear.h:287
AdditionalData additional_data
Definition nonlinear.h:361
Teuchos::RCP< Teuchos::ParameterList > parameters_nox
Definition nonlinear.h:382
PETScWrappers::NonlinearSolverData additional_data_petsc_snes
Definition nonlinear.h:390
void solve(VectorType &initial_guess_and_solution)
Definition nonlinear.h:615
void select(const typename AdditionalData::SolverType &type)
Definition nonlinear.h:519
std::function< void(const VectorType &current_u)> setup_jacobian
Definition nonlinear.h:330
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition nonlinear.h:355
std::function< void(const VectorType &src, VectorType &dst)> solve_with_jacobian
Definition petsc_snes.h:445
std::function< void(const VectorType &x, VectorType &res)> residual
Definition petsc_snes.h:387
unsigned int solve(VectorType &x)
std::function< void(const VectorType &x)> setup_jacobian
Definition petsc_snes.h:431
std::function< void(const VectorType &current_u, const VectorType &current_f)> setup_jacobian
Definition kinsol.h:561
unsigned int solve(VectorType &initial_guess_and_solution)
Definition kinsol.cc:223
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition kinsol.h:496
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition kinsol.h:604
std::function< void(VectorType &)> reinit_vector
Definition kinsol.h:481
std::function< void(const VectorType &y, VectorType &x, const double tolerance)> solve_with_jacobian
Definition nox.h:286
std::function< void(const VectorType &u, VectorType &F)> residual
Definition nox.h:204
unsigned int solve(VectorType &solution)
std::function< void(const VectorType &current_u)> setup_jacobian
Definition nox.h:218
#define DEAL_II_DEPRECATED
Definition config.h:207
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)