Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20:02+00:00
\(\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
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
214
223
228 void
229 select(const typename AdditionalData::SolverType &type);
230
234 void
236
241#ifdef DEAL_II_TRILINOS_WITH_NOX
242 void
243 set_data(
246 const Teuchos::RCP<Teuchos::ParameterList> &parameters =
247 Teuchos::rcp(new Teuchos::ParameterList));
248#endif
249
254#ifdef DEAL_II_WITH_SUNDIALS
255 void
258#endif
259
264#ifdef DEAL_II_WITH_PETSC
265 void
267#endif
268
275 void
276 solve(VectorType &initial_guess_and_solution);
277
284 std::function<void(VectorType &)> reinit_vector;
285
299 std::function<void(const VectorType &src, VectorType &dst)> residual;
300
327 std::function<void(const VectorType &current_u)> setup_jacobian;
328
350 std::function<
351 void(const VectorType &rhs, VectorType &dst, const double tolerance)>
353
354private:
359
360private:
365
369#ifdef DEAL_II_WITH_SUNDIALS
371#endif
372
376#ifdef DEAL_II_TRILINOS_WITH_NOX
379 Teuchos::RCP<Teuchos::ParameterList> parameters_nox =
380 Teuchos::rcp(new Teuchos::ParameterList);
381#endif
382
386#ifdef DEAL_II_WITH_PETSC
388#endif
389
396 void
397 solve_with_petsc(VectorType &initial_guess_and_solution);
398};
399
400
401
402template <typename VectorType>
403void
405 const AdditionalData &additional_data)
406{
407#ifdef DEAL_II_WITH_SUNDIALS
408 // These if statements pass on the strategy to the other nonlinear solvers
409 if (additional_data.strategy ==
411 additional_data_kinsol.strategy =
413 else if (additional_data.strategy ==
415 additional_data_kinsol.strategy =
417 else if (additional_data.strategy ==
419 additional_data_kinsol.strategy =
421
422 // Setting data points in the KINSOL class from the NonlinearSolverSelector
423 // class
424 additional_data_kinsol.maximum_non_linear_iterations =
425 additional_data.maximum_non_linear_iterations;
426 additional_data_kinsol.function_tolerance =
427 additional_data.function_tolerance;
428 additional_data_kinsol.step_tolerance = additional_data.step_tolerance;
429 additional_data_kinsol.anderson_subspace_size =
430 additional_data.anderson_subspace_size;
431#endif
432
433// Do the same thing we did above but with NOX
434#ifdef DEAL_II_TRILINOS_WITH_NOX
435 parameters_nox->set("Nonlinear Solver", "Line Search Based");
436 Teuchos::ParameterList &Line_Search = parameters_nox->sublist("Line Search");
437 Line_Search.set("Method", "Full Step");
438
439 additional_data_nox.max_iter = additional_data.maximum_non_linear_iterations;
440 additional_data_nox.abs_tol = additional_data.function_tolerance;
441 additional_data_nox.rel_tol = additional_data.relative_tolerance;
442#endif
443
444#ifdef DEAL_II_WITH_PETSC
445 additional_data_petsc_snes.options_prefix = "";
446
447 if (additional_data.anderson_subspace_size > 0 &&
448 additional_data.strategy ==
450 {
451 additional_data_petsc_snes.snes_type = "anderson";
452 // TODO additional_data.anderson_subspace_size;
453 }
454 else if (additional_data.strategy ==
456 {
457 additional_data_petsc_snes.snes_type = "newtonls";
458 additional_data_petsc_snes.snes_linesearch_type = "bt";
459 }
460 else if (additional_data.strategy ==
462 {
463 additional_data_petsc_snes.snes_linesearch_type = "newtonls";
464 additional_data_petsc_snes.snes_linesearch_type = "basic";
465 }
466 else if (additional_data.strategy ==
468 additional_data_petsc_snes.snes_type = "nrichardson";
469
470 additional_data_petsc_snes.absolute_tolerance =
471 additional_data.function_tolerance;
472 additional_data_petsc_snes.relative_tolerance =
473 additional_data.relative_tolerance;
474 additional_data_petsc_snes.step_tolerance = additional_data.step_tolerance;
475 additional_data_petsc_snes.maximum_non_linear_iterations =
476 additional_data.maximum_non_linear_iterations;
477 additional_data_petsc_snes.max_n_function_evaluations = -1;
478
479#endif
480}
481
482
483
484template <typename VectorType>
486
487
488
489template <typename VectorType>
491 const AdditionalData &additional_data)
492 : additional_data(additional_data)
493{
495}
496
497
498
499template <typename VectorType>
501 const AdditionalData &additional_data,
502 const MPI_Comm &mpi_communicator)
503 : additional_data(additional_data)
504 , mpi_communicator(mpi_communicator)
505{
507}
508
509
510
511template <typename VectorType>
512void
514 const typename AdditionalData::SolverType &type)
515{
516 additional_data.solver_type = type;
517}
518
519
520
521template <typename VectorType>
523 const SolverType &solver_type,
524 const SolutionStrategy &strategy,
525 const unsigned int maximum_non_linear_iterations,
526 const double function_tolerance,
527 const double relative_tolerance,
528 const double step_tolerance,
529 const unsigned int anderson_subspace_size)
530 : solver_type(solver_type)
531 , strategy(strategy)
532 , maximum_non_linear_iterations(maximum_non_linear_iterations)
533 , function_tolerance(function_tolerance)
534 , relative_tolerance(relative_tolerance)
535 , step_tolerance(step_tolerance)
536 , anderson_subspace_size(anderson_subspace_size)
537{}
538
539
540
541#ifdef DEAL_II_TRILINOS_WITH_NOX
542template <typename VectorType>
543void
547 const Teuchos::RCP<Teuchos::ParameterList> &parameters)
548{
550 parameters_nox = parameters;
551}
552#endif
553
554
555
556#ifdef DEAL_II_WITH_SUNDIALS
557template <typename VectorType>
558void
564#endif
565
566
567
568template <typename VectorType>
569void
571 VectorType & /*initial_guess_and_solution*/)
572{
573 AssertThrow(false,
574 ExcMessage("PETSc SNES requires you to use PETSc vectors."));
575}
576
577
578
579#ifdef DEAL_II_WITH_PETSC
580template <>
581void
583 PETScWrappers::MPI::Vector &initial_guess_and_solution)
584{
585 PETScWrappers::NonlinearSolver<PETScWrappers::MPI::Vector> nonlinear_solver(
587
588 nonlinear_solver.residual = residual;
589
590 nonlinear_solver.setup_jacobian = setup_jacobian;
591
592 nonlinear_solver.solve_with_jacobian =
593 [&](const PETScWrappers::MPI::Vector &src,
595 // PETSc does not gives a tolerance, so we have to choose something
596 // reasonable to provide to the user:
597 const double tolerance = 1e-6;
598 solve_with_jacobian(src, dst, tolerance);
599 };
600
601 nonlinear_solver.solve(initial_guess_and_solution);
602}
603#endif
604
605
606
607template <typename VectorType>
608void
610 VectorType &initial_guess_and_solution)
611{
612 // The "auto" solver_type will default to kinsol, however if KINSOL is not
613 // available then we will use NOX.
615 {
616#ifdef DEAL_II_WITH_PETSC
618#endif
619#ifdef DEAL_II_TRILINOS_WITH_NOX
621#endif
622#ifdef DEAL_II_WITH_SUNDIALS
624#endif
625
626 // If "auto" is still the solver type we cannot solve the problem
628 AssertThrow(false, ExcMessage("No valid solver type."));
629 }
630
632 {
633#ifdef DEAL_II_WITH_SUNDIALS
636
637 nonlinear_solver.reinit_vector = reinit_vector;
638 nonlinear_solver.residual = residual;
639
640 // We cannot simply set these two functions equal to each other
641 // because they have a different number of inputs.
642 nonlinear_solver.setup_jacobian = [&](const VectorType &current_u,
643 const VectorType & /*current_f*/) {
645 };
646
647 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
648
649 nonlinear_solver.solve(initial_guess_and_solution);
650#else
652 false, ExcMessage("You do not have SUNDIALS configured with deal.II!"));
653#endif
654 }
656 {
657#ifdef DEAL_II_TRILINOS_WITH_NOX
658
661
662 // Do the same thing for NOX that we did with KINSOL.
663 nonlinear_solver.residual = residual;
664
665 // setup_jacobian for NOX has the same number of arguments for the same
666 // function in NonlinearSolverSelector.
667 nonlinear_solver.setup_jacobian = setup_jacobian;
668
669 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
670
671 nonlinear_solver.solve(initial_guess_and_solution);
672#else
674 false, ExcMessage("You do not have Trilinos configured with deal.II"));
675#endif
676 }
677 else if (additional_data.solver_type ==
679 {
680 // Forward to internal function specializations, which throws for
681 // non-supported vector types:
682 solve_with_petsc(initial_guess_and_solution);
683 }
684 else
685 {
686 const std::string solvers =
687#ifdef DEAL_II_WITH_SUNDIALS
688 "kinsol\n"
689#endif
690#ifdef DEAL_II_TRILINOS_WITH_NOX
691 "NOX\n"
692#endif
693#ifdef DEAL_II_WITH_PETSC
694 "SNES\n"
695#endif
696 ;
697
698 AssertThrow(false,
700 "Invalid nonlinear solver specified. "
701 "The solvers available in your installation are:\n" +
702 solvers));
703 }
704}
705
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:522
TrilinosWrappers::NOXSolver< VectorType >::AdditionalData additional_data_nox
Definition nonlinear.h:378
void set_data(const AdditionalData &additional_data)
Definition nonlinear.h:404
void solve_with_petsc(VectorType &initial_guess_and_solution)
Definition nonlinear.h:570
void set_data(const typename PETScWrappers::NonlinearSolverData &additional_data)
SUNDIALS::KINSOL< VectorType >::AdditionalData additional_data_kinsol
Definition nonlinear.h:370
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition nonlinear.h:299
std::function< void(VectorType &)> reinit_vector
Definition nonlinear.h:284
AdditionalData additional_data
Definition nonlinear.h:358
Teuchos::RCP< Teuchos::ParameterList > parameters_nox
Definition nonlinear.h:379
PETScWrappers::NonlinearSolverData additional_data_petsc_snes
Definition nonlinear.h:387
void solve(VectorType &initial_guess_and_solution)
Definition nonlinear.h:609
void select(const typename AdditionalData::SolverType &type)
Definition nonlinear.h:513
std::function< void(const VectorType &current_u)> setup_jacobian
Definition nonlinear.h:327
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition nonlinear.h:352
std::function< void(const VectorType &current_u, const VectorType &current_f)> setup_jacobian
Definition kinsol.h:562
unsigned int solve(VectorType &initial_guess_and_solution)
Definition kinsol.cc:223
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition kinsol.h:497
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition kinsol.h:605
std::function< void(VectorType &)> reinit_vector
Definition kinsol.h:482
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_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)