deal.II version GIT relicensing-2659-g040196caa3 2025-02-18 14:20:01+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#ifndef dealii_nonlinear_h
16#define dealii_nonlinear_h
17
18#include <deal.II/base/config.h>
19
22
24#include <deal.II/lac/vector.h>
25
27
29
30
32
79template <typename VectorType = Vector<double>>
81{
82public:
90 {
91 public:
112
114 {
128 /*
129 * Use the PETSc SNES solver.
130 */
132 };
133
150 const unsigned int maximum_non_linear_iterations = 200,
151 const double function_tolerance = 1e-8,
152 const double relative_tolerance = 1e-5,
153 const double step_tolerance = 0.0,
154 const unsigned int anderson_subspace_size = 0);
155
162
170
175
183
190 const double relative_tolerance;
191
199
207 };
208
213
222
231
236 void
237 select(const typename AdditionalData::SolverType &type);
238
242 void
244
249#ifdef DEAL_II_TRILINOS_WITH_NOX
250 void
251 set_data(
254 const Teuchos::RCP<Teuchos::ParameterList> &parameters =
255 Teuchos::rcp(new Teuchos::ParameterList));
256#endif
257
262#ifdef DEAL_II_WITH_SUNDIALS
263 void
266#endif
267
272#ifdef DEAL_II_WITH_PETSC
273 void
275#endif
276
283 void
285
292 std::function<void(VectorType &)> reinit_vector;
293
307 std::function<void(const VectorType &src, VectorType &dst)> residual;
308
335 std::function<void(const VectorType &current_u)> setup_jacobian;
336
358 std::function<
359 void(const VectorType &rhs, VectorType &dst, const double tolerance)>
361
362private:
367
368private:
373
377#ifdef DEAL_II_WITH_SUNDIALS
379#endif
380
384#ifdef DEAL_II_TRILINOS_WITH_NOX
387 Teuchos::RCP<Teuchos::ParameterList> parameters_nox =
388 Teuchos::rcp(new Teuchos::ParameterList);
389#endif
390
394#ifdef DEAL_II_WITH_PETSC
396#endif
397
404 void
406};
407
408
409
410template <typename VectorType>
411void
413 const AdditionalData &additional_data)
414{
415#ifdef DEAL_II_WITH_SUNDIALS
416 // These if statements pass on the strategy to the other nonlinear solvers
417 if (additional_data.strategy ==
419 additional_data_kinsol.strategy =
421 else if (additional_data.strategy ==
423 additional_data_kinsol.strategy =
425 else if (additional_data.strategy ==
427 additional_data_kinsol.strategy =
429
430 // Setting data points in the KINSOL class from the NonlinearSolverSelector
431 // class
432 additional_data_kinsol.maximum_non_linear_iterations =
433 additional_data.maximum_non_linear_iterations;
434 additional_data_kinsol.function_tolerance =
435 additional_data.function_tolerance;
436 additional_data_kinsol.step_tolerance = additional_data.step_tolerance;
437 additional_data_kinsol.anderson_subspace_size =
438 additional_data.anderson_subspace_size;
439#endif
440
441// Do the same thing we did above but with NOX
442#ifdef DEAL_II_TRILINOS_WITH_NOX
443 parameters_nox->set("Nonlinear Solver", "Line Search Based");
444 Teuchos::ParameterList &Line_Search = parameters_nox->sublist("Line Search");
445 Line_Search.set("Method", "Full Step");
446
447 additional_data_nox.max_iter = additional_data.maximum_non_linear_iterations;
448 additional_data_nox.abs_tol = additional_data.function_tolerance;
449 additional_data_nox.rel_tol = additional_data.relative_tolerance;
450#endif
451
452#ifdef DEAL_II_WITH_PETSC
453 additional_data_petsc_snes.options_prefix = "";
454
455 if (additional_data.anderson_subspace_size > 0 &&
456 additional_data.strategy ==
458 {
459 additional_data_petsc_snes.snes_type = "anderson";
460 // TODO additional_data.anderson_subspace_size;
461 }
462 else if (additional_data.strategy ==
464 {
465 additional_data_petsc_snes.snes_type = "newtonls";
466 additional_data_petsc_snes.snes_linesearch_type = "bt";
467 }
468 else if (additional_data.strategy ==
470 {
471 additional_data_petsc_snes.snes_linesearch_type = "newtonls";
472 additional_data_petsc_snes.snes_linesearch_type = "basic";
473 }
474 else if (additional_data.strategy ==
476 additional_data_petsc_snes.snes_type = "nrichardson";
477
478 additional_data_petsc_snes.absolute_tolerance =
479 additional_data.function_tolerance;
480 additional_data_petsc_snes.relative_tolerance =
481 additional_data.relative_tolerance;
482 additional_data_petsc_snes.step_tolerance = additional_data.step_tolerance;
483 additional_data_petsc_snes.maximum_non_linear_iterations =
484 additional_data.maximum_non_linear_iterations;
485 additional_data_petsc_snes.max_n_function_evaluations = -1;
486
487#endif
488}
489
490
491
492template <typename VectorType>
496
497
498
499template <typename VectorType>
501 const AdditionalData &additional_data)
502 : additional_data(additional_data)
503 , mpi_communicator(MPI_COMM_SELF)
504{
506}
507
508
509
510template <typename VectorType>
512 const AdditionalData &additional_data,
513 const MPI_Comm &mpi_communicator)
514 : additional_data(additional_data)
515 , mpi_communicator(mpi_communicator)
516{
518}
519
520
521
522template <typename VectorType>
523void
525 const typename AdditionalData::SolverType &type)
526{
527 additional_data.solver_type = type;
528}
529
530
531
532template <typename VectorType>
534 const SolverType &solver_type,
535 const SolutionStrategy &strategy,
536 const unsigned int maximum_non_linear_iterations,
537 const double function_tolerance,
538 const double relative_tolerance,
539 const double step_tolerance,
540 const unsigned int anderson_subspace_size)
541 : solver_type(solver_type)
542 , strategy(strategy)
543 , maximum_non_linear_iterations(maximum_non_linear_iterations)
544 , function_tolerance(function_tolerance)
545 , relative_tolerance(relative_tolerance)
546 , step_tolerance(step_tolerance)
547 , anderson_subspace_size(anderson_subspace_size)
548{}
549
550
551
552#ifdef DEAL_II_TRILINOS_WITH_NOX
553template <typename VectorType>
554void
558 const Teuchos::RCP<Teuchos::ParameterList> &parameters)
559{
561 parameters_nox = parameters;
562}
563#endif
564
565
566
567#ifdef DEAL_II_WITH_SUNDIALS
568template <typename VectorType>
569void
575#endif
576
577
578
579template <typename VectorType>
580void
582 VectorType & /*initial_guess_and_solution*/)
583{
584 AssertThrow(false,
585 ExcMessage("PETSc SNES requires you to use PETSc vectors."));
586}
587
588
589
590#ifdef DEAL_II_WITH_PETSC
591template <>
592void
595{
598
599 nonlinear_solver.residual = residual;
600
601 nonlinear_solver.setup_jacobian = setup_jacobian;
602
603 nonlinear_solver.solve_with_jacobian =
604 [&](const PETScWrappers::MPI::Vector &src,
606 // PETSc does not gives a tolerance, so we have to choose something
607 // reasonable to provide to the user:
608 const double tolerance = 1e-6;
609 solve_with_jacobian(src, dst, tolerance);
610 };
611
613}
614#endif
615
616
617
618template <typename VectorType>
619void
621 VectorType &initial_guess_and_solution)
622{
623 // The "auto" solver_type will default to kinsol, however if KINSOL is not
624 // available then we will use NOX.
626 {
627#ifdef DEAL_II_WITH_PETSC
629#endif
630#ifdef DEAL_II_TRILINOS_WITH_NOX
632#endif
633#ifdef DEAL_II_WITH_SUNDIALS
635#endif
636
637 // If "auto" is still the solver type we cannot solve the problem
639 AssertThrow(false, ExcMessage("No valid solver type."));
640 }
641
643 {
644#ifdef DEAL_II_WITH_SUNDIALS
647
648 nonlinear_solver.reinit_vector = reinit_vector;
649 nonlinear_solver.residual = residual;
650
651 // We cannot simply set these two functions equal to each other
652 // because they have a different number of inputs.
653 nonlinear_solver.setup_jacobian = [&](const VectorType &current_u,
654 const VectorType & /*current_f*/) {
656 };
657
658 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
659
661#else
663 false, ExcMessage("You do not have SUNDIALS configured with deal.II!"));
664#endif
665 }
667 {
668#ifdef DEAL_II_TRILINOS_WITH_NOX
669
672
673 // Do the same thing for NOX that we did with KINSOL.
674 nonlinear_solver.residual = residual;
675
676 // setup_jacobian for NOX has the same number of arguments for the same
677 // function in NonlinearSolverSelector.
678 nonlinear_solver.setup_jacobian = setup_jacobian;
679
680 nonlinear_solver.solve_with_jacobian = solve_with_jacobian;
681
683#else
685 false, ExcMessage("You do not have Trilinos configured with deal.II"));
686#endif
687 }
688 else if (additional_data.solver_type ==
690 {
691 // Forward to internal function specializations, which throws for
692 // non-supported vector types:
694 }
695 else
696 {
697 const std::string solvers = ""
698#ifdef DEAL_II_WITH_SUNDIALS
699 "kinsol\n"
700#endif
701#ifdef DEAL_II_TRILINOS_WITH_NOX
702 "NOX\n"
703#endif
704#ifdef DEAL_II_WITH_PETSC
705 "SNES\n"
706#endif
707 ;
708
709 AssertThrow(false,
711 "Invalid nonlinear solver specified. "
712 "The solvers available in your installation are:\n" +
713 solvers));
714 }
715}
716
718
719#endif
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:533
TrilinosWrappers::NOXSolver< VectorType >::AdditionalData additional_data_nox
Definition nonlinear.h:386
void set_data(const AdditionalData &additional_data)
Definition nonlinear.h:412
void solve_with_petsc(VectorType &initial_guess_and_solution)
Definition nonlinear.h:581
void set_data(const typename PETScWrappers::NonlinearSolverData &additional_data)
SUNDIALS::KINSOL< VectorType >::AdditionalData additional_data_kinsol
Definition nonlinear.h:378
std::function< void(const VectorType &src, VectorType &dst)> residual
Definition nonlinear.h:307
std::function< void(VectorType &)> reinit_vector
Definition nonlinear.h:292
AdditionalData additional_data
Definition nonlinear.h:366
Teuchos::RCP< Teuchos::ParameterList > parameters_nox
Definition nonlinear.h:387
PETScWrappers::NonlinearSolverData additional_data_petsc_snes
Definition nonlinear.h:395
void solve(VectorType &initial_guess_and_solution)
Definition nonlinear.h:620
void select(const typename AdditionalData::SolverType &type)
Definition nonlinear.h:524
std::function< void(const VectorType &current_u)> setup_jacobian
Definition nonlinear.h:335
std::function< void(const VectorType &rhs, VectorType &dst, const double tolerance)> solve_with_jacobian
Definition nonlinear.h:360
friend class Tensor
Definition tensor.h:865
#define DEAL_II_DEPRECATED
Definition config.h:228
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:518
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:519
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)