Reference documentation for deal.II version 9.5.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// Copyright (C) 2023 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#include <deal.II/base/config.h>
17
19
21
23
25
26
28
75template <typename VectorType = Vector<double>>
77{
78public:
86 {
87 public:
94 {
107 };
108
110 {
124 /*
125 * Use the PETSc SNES solver.
126 */
128 };
129
146 const unsigned int maximum_non_linear_iterations = 200,
147 const double function_tolerance = 1e-8,
148 const double relative_tolerance = 1e-5,
149 const double step_tolerance = 0.0,
150 const unsigned int anderson_subspace_size = 0);
151
158
166
171
179
186 const double relative_tolerance;
187
195
203 };
204
209
215
223 const MPI_Comm & mpi_communicator);
224
229 void
230 select(const typename AdditionalData::SolverType &type);
231
235 void
237
242#ifdef DEAL_II_WITH_TRILINOS
243 void
244 set_data(
247 const Teuchos::RCP<Teuchos::ParameterList> &parameters =
248 Teuchos::rcp(new Teuchos::ParameterList));
249#endif
250
255#ifdef DEAL_II_WITH_SUNDIALS
256 void
259#endif
260
265#ifdef DEAL_II_WITH_PETSC
266 void
268#endif
269
276 void
277 solve(VectorType &initial_guess_and_solution);
278
285 std::function<void(VectorType &)> reinit_vector;
286
300 std::function<void(const VectorType &src, VectorType &dst)> residual;
301
328 std::function<void(const VectorType &current_u)> setup_jacobian;
329
351 std::function<
352 void(const VectorType &rhs, VectorType &dst, const double tolerance)>
354
355private:
360
361private:
366
370#ifdef DEAL_II_WITH_SUNDIALS
372#endif
373
377#ifdef DEAL_II_WITH_TRILINOS
380 Teuchos::RCP<Teuchos::ParameterList> parameters_nox =
381 Teuchos::rcp(new Teuchos::ParameterList);
382#endif
383
387#ifdef DEAL_II_WITH_PETSC
389#endif
390
397 void
398 solve_with_petsc(VectorType &initial_guess_and_solution);
399};
400
401
402
403template <typename VectorType>
404void
406 const AdditionalData &additional_data)
407{
408#ifdef DEAL_II_WITH_SUNDIALS
409 // These if statements pass on the strategy to the other nonlinear solvers
410 if (additional_data.strategy ==
412 additional_data_kinsol.strategy =
414 else if (additional_data.strategy ==
416 additional_data_kinsol.strategy =
418 else if (additional_data.strategy ==
420 additional_data_kinsol.strategy =
422
423 // Setting data points in the KINSOL class from the NonlinearSolverSelector
424 // class
425 additional_data_kinsol.maximum_non_linear_iterations =
426 additional_data.maximum_non_linear_iterations;
427 additional_data_kinsol.function_tolerance =
428 additional_data.function_tolerance;
429 additional_data_kinsol.step_tolerance = additional_data.step_tolerance;
430 additional_data_kinsol.anderson_subspace_size =
431 additional_data.anderson_subspace_size;
432#endif
433
434// Do the same thing we did above but with NOX
435#ifdef DEAL_II_WITH_TRILINOS
436 parameters_nox->set("Nonlinear Solver", "Line Search Based");
437 Teuchos::ParameterList &Line_Search = parameters_nox->sublist("Line Search");
438 Line_Search.set("Method", "Full Step");
439
440 additional_data_nox.max_iter = additional_data.maximum_non_linear_iterations;
441 additional_data_nox.abs_tol = additional_data.function_tolerance;
442 additional_data_nox.rel_tol = additional_data.relative_tolerance;
443#endif
444
445#ifdef DEAL_II_WITH_PETSC
446 additional_data_petsc_snes.options_prefix = "";
447
448 if (additional_data.anderson_subspace_size > 0 &&
449 additional_data.strategy ==
451 {
452 additional_data_petsc_snes.snes_type = "anderson";
453 // TODO additional_data.anderson_subspace_size;
454 }
455 else if (additional_data.strategy ==
457 {
458 additional_data_petsc_snes.snes_type = "newtonls";
459 additional_data_petsc_snes.snes_linesearch_type = "bt";
460 }
461 else if (additional_data.strategy ==
463 {
464 additional_data_petsc_snes.snes_linesearch_type = "newtonls";
465 additional_data_petsc_snes.snes_linesearch_type = "basic";
466 }
467 else if (additional_data.strategy ==
469 additional_data_petsc_snes.snes_type = "nrichardson";
470
471 additional_data_petsc_snes.absolute_tolerance =
472 additional_data.function_tolerance;
473 additional_data_petsc_snes.relative_tolerance =
474 additional_data.relative_tolerance;
475 additional_data_petsc_snes.step_tolerance = additional_data.step_tolerance;
476 additional_data_petsc_snes.maximum_non_linear_iterations =
477 additional_data.maximum_non_linear_iterations;
478 additional_data_petsc_snes.max_n_function_evaluations = -1;
479
480#endif
481}
482
483
484
485template <typename VectorType>
487
488
489
490template <typename VectorType>
492 const AdditionalData &additional_data)
493 : additional_data(additional_data)
494{
496}
497
498
499
500template <typename VectorType>
502 const AdditionalData &additional_data,
503 const MPI_Comm & mpi_communicator)
504 : additional_data(additional_data)
505 , mpi_communicator(mpi_communicator)
506{
508}
509
510
511
512template <typename VectorType>
513void
515 const typename AdditionalData::SolverType &type)
516{
517 additional_data.solver_type = type;
518}
519
520
521
522template <typename VectorType>
524 const SolverType & solver_type,
525 const SolutionStrategy &strategy,
526 const unsigned int maximum_non_linear_iterations,
527 const double function_tolerance,
528 const double relative_tolerance,
529 const double step_tolerance,
530 const unsigned int anderson_subspace_size)
531 : solver_type(solver_type)
532 , strategy(strategy)
533 , maximum_non_linear_iterations(maximum_non_linear_iterations)
534 , function_tolerance(function_tolerance)
535 , relative_tolerance(relative_tolerance)
536 , step_tolerance(step_tolerance)
537 , anderson_subspace_size(anderson_subspace_size)
538{}
539
540
541
542#ifdef DEAL_II_WITH_TRILINOS
543template <typename VectorType>
544void
548 const Teuchos::RCP<Teuchos::ParameterList> &parameters)
549{
551 parameters_nox = parameters;
552}
553#endif
554
555
556
557#ifdef DEAL_II_WITH_SUNDIALS
558template <typename VectorType>
559void
562{
564}
565#endif
566
567
568
569template <typename VectorType>
570void
572 VectorType & /*initial_guess_and_solution*/)
573{
574 AssertThrow(false,
575 ExcMessage("PETSc SNES requires you to use PETSc vectors."));
576}
577
578
579
580#ifdef DEAL_II_WITH_PETSC
581template <>
582void
584 PETScWrappers::MPI::Vector &initial_guess_and_solution)
585{
588
589 nonlinear_solver.residual = residual;
590
591 nonlinear_solver.setup_jacobian = setup_jacobian;
592
593 nonlinear_solver.solve_with_jacobian =
594 [&](const PETScWrappers::MPI::Vector &src,
596 // PETSc does not gives a tolerance, so we have to choose something
597 // reasonable to provide to the user:
598 const double tolerance = 1e-6;
599 solve_with_jacobian(src, dst, tolerance);
600 };
601
602 nonlinear_solver.solve(initial_guess_and_solution);
603}
604#endif
605
606
607
608template <typename VectorType>
609void
611 VectorType &initial_guess_and_solution)
612{
613 // The "auto" solver_type will default to kinsol, however if KINSOL is not
614 // available then we will use NOX.
616 {
617#ifdef DEAL_II_WITH_PETSC
619#endif
620#ifdef DEAL_II_WITH_TRILINOS
622#endif
623#ifdef DEAL_II_WITH_SUNDIALS
625#endif
626
627 // If "auto" is still the solver type we cannot solve the problem
629 AssertThrow(false, ExcMessage("No valid solver type."));
630 }
631
633 {
634#ifdef DEAL_II_WITH_SUNDIALS
637
638 nonlinear_solver.reinit_vector = reinit_vector;
639 nonlinear_solver.residual = residual;
640
641 // We cannot simply set these two functions equal to each other
642 // because they have a different number of inputs.
643 nonlinear_solver.setup_jacobian = [&](const VectorType &current_u,
644 const VectorType & /*current_f*/) {
646 };
647
648 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
649
650 nonlinear_solver.solve(initial_guess_and_solution);
651#else
653 false, ExcMessage("You do not have SUNDIALS configured with deal.II!"));
654#endif
655 }
657 {
658#ifdef DEAL_II_WITH_TRILINOS
659
662
663 // Do the same thing for NOX that we did with KINSOL.
664 nonlinear_solver.residual = residual;
665
666 // setup_jacobian for NOX has the same number of arguments for the same
667 // function in NonlinearSolverSelector.
668 nonlinear_solver.setup_jacobian = setup_jacobian;
669
670 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
671
672 nonlinear_solver.solve(initial_guess_and_solution);
673#else
675 false, ExcMessage("You do not have Trilinos configured with deal.II"));
676#endif
677 }
678 else if (additional_data.solver_type ==
680 {
681 // Forward to internal function specializations, which throws for
682 // non-supported vector types:
683 solve_with_petsc(initial_guess_and_solution);
684 }
685 else
686 {
687 const std::string solvers =
688#ifdef DEAL_II_WITH_SUNDIALS
689 "kinsol\n"
690#endif
691#ifdef DEAL_II_WITH_TRILINOS
692 "NOX\n"
693#endif
694#ifdef DEAL_II_WITH_PETSC
695 "SNES\n"
696#endif
697 ;
698
699 AssertThrow(false,
701 "Invalid nonlinear solver specified. "
702 "The solvers available in your installation are:\n" +
703 solvers));
704 }
705}
706
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:523
TrilinosWrappers::NOXSolver< VectorType >::AdditionalData additional_data_nox
Definition nonlinear.h:379
void set_data(const AdditionalData &additional_data)
Definition nonlinear.h:405
void solve_with_petsc(VectorType &initial_guess_and_solution)
Definition nonlinear.h:571
void set_data(const typename PETScWrappers::NonlinearSolverData &additional_data)
SUNDIALS::KINSOL< VectorType >::AdditionalData additional_data_kinsol
Definition nonlinear.h:371
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition nonlinear.h:300
std::function< void(VectorType &)> reinit_vector
Definition nonlinear.h:285
AdditionalData additional_data
Definition nonlinear.h:359
Teuchos::RCP< Teuchos::ParameterList > parameters_nox
Definition nonlinear.h:380
PETScWrappers::NonlinearSolverData additional_data_petsc_snes
Definition nonlinear.h:388
void solve(VectorType &initial_guess_and_solution)
Definition nonlinear.h:610
void select(const typename AdditionalData::SolverType &type)
Definition nonlinear.h:514
std::function< void(const VectorType &current_u)> setup_jacobian
Definition nonlinear.h:328
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition nonlinear.h:353
std::function< void(const VectorType &src, VectorType &dst)> solve_with_jacobian
Definition petsc_snes.h:432
std::function< void(const VectorType &x, VectorType &res)> residual
Definition petsc_snes.h:374
unsigned int solve(VectorType &x)
std::function< void(const VectorType &x)> setup_jacobian
Definition petsc_snes.h:418
std::function< void(const VectorType &current_u, const VectorType &current_f)> setup_jacobian
Definition kinsol.h:506
unsigned int solve(VectorType &initial_guess_and_solution)
Definition kinsol.cc:187
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition kinsol.h:441
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition kinsol.h:621
std::function< void(VectorType &)> reinit_vector
Definition kinsol.h:426
std::function< void(const VectorType &y, VectorType &x, const double tolerance)> solve_with_jacobian
Definition nox.h:284
std::function< void(const VectorType &u, VectorType &F)> residual
Definition nox.h:205
unsigned int solve(VectorType &solution)
std::function< void(const VectorType &current_u)> setup_jacobian
Definition nox.h:219
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)