www.mooseframework.org
NonlinearSystemBase.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #ifndef NONLINEARSYSTEMBASE_H
16 #define NONLINEARSYSTEMBASE_H
17 
18 #include "SystemBase.h"
19 #include "KernelWarehouse.h"
20 #include "ConstraintWarehouse.h"
21 #include "MooseObjectWarehouse.h"
22 
23 #include "libmesh/transient_system.h"
24 #include "libmesh/nonlinear_implicit_system.h"
25 
26 // Forward declarations
27 class FEProblemBase;
29 class JacobianBlock;
30 class TimeIntegrator;
31 class Predictor;
32 class ElementDamper;
33 class NodalDamper;
34 class GeneralDamper;
36 class IntegratedBC;
37 class NodalBC;
38 class PresetNodalBC;
39 class DGKernel;
40 class InterfaceKernel;
41 class ScalarKernel;
42 class DiracKernel;
43 class NodalKernel;
44 class Split;
45 
46 // libMesh forward declarations
47 namespace libMesh
48 {
49 template <typename T>
50 class NumericVector;
51 template <typename T>
52 class SparseMatrix;
53 }
54 
61 {
62 public:
63  NonlinearSystemBase(FEProblemBase & problem, System & sys, const std::string & name);
64  virtual ~NonlinearSystemBase();
65 
66  virtual void init() override;
67 
71  void turnOffJacobian();
72 
73  virtual void addExtraVectors() override;
74  virtual void solve() override = 0;
75  virtual void restoreSolutions() override;
76 
80  virtual void stopSolve() = 0;
81 
82  virtual NonlinearSolver<Number> * nonlinearSolver() = 0;
83 
84  virtual unsigned int getCurrentNonlinearIterationNumber() = 0;
85 
90  virtual bool computingInitialResidual() { return _computing_initial_residual; }
91 
92  // Setup Functions ////
93  virtual void initialSetup();
94  virtual void timestepSetup();
95 
96  virtual void setupFiniteDifferencedPreconditioner() = 0;
97  void setupFieldDecomposition();
98 
99  bool haveFiniteDifferencedPreconditioner() { return _use_finite_differenced_preconditioner; }
100  bool haveFieldSplitPreconditioner() { return _use_field_split_preconditioner; }
101 
106  virtual bool converged() = 0;
107 
114  void
115  addTimeIntegrator(const std::string & type, const std::string & name, InputParameters parameters);
116 
123  virtual void
124  addKernel(const std::string & kernel_name, const std::string & name, InputParameters parameters);
125 
126  virtual void addEigenKernels(std::shared_ptr<KernelBase> /*kernel*/, THREAD_ID /*tid*/){};
127 
134  virtual void addNodalKernel(const std::string & kernel_name,
135  const std::string & name,
136  InputParameters parameters);
137 
144  void addScalarKernel(const std::string & kernel_name,
145  const std::string & name,
146  InputParameters parameters);
147 
154  void addBoundaryCondition(const std::string & bc_name,
155  const std::string & name,
156  InputParameters parameters);
157 
164  void
165  addConstraint(const std::string & c_name, const std::string & name, InputParameters parameters);
166 
173  void addDiracKernel(const std::string & kernel_name,
174  const std::string & name,
175  InputParameters parameters);
176 
183  void
184  addDGKernel(std::string dg_kernel_name, const std::string & name, InputParameters parameters);
185 
192  void addInterfaceKernel(std::string interface_kernel_name,
193  const std::string & name,
194  InputParameters parameters);
195 
202  void
203  addDamper(const std::string & damper_name, const std::string & name, InputParameters parameters);
204 
211  void
212  addSplit(const std::string & split_name, const std::string & name, InputParameters parameters);
213 
218  std::shared_ptr<Split> getSplit(const std::string & name);
219 
220  void zeroVectorForResidual(const std::string & vector_name);
221 
222  void setInitialSolution();
223 
227  void setConstraintSlaveValues(NumericVector<Number> & solution, bool displaced);
228 
236  void constraintResiduals(NumericVector<Number> & residual, bool displaced);
237 
243  void computeResidual(NumericVector<Number> & residual, Moose::KernelType type = Moose::KT_ALL);
244 
248  void
249  findImplicitGeometricCouplingEntries(GeometricSearchData & geom_search_data,
250  std::map<dof_id_type, std::vector<dof_id_type>> & graph);
251 
257  void addImplicitGeometricCouplingEntries(SparseMatrix<Number> & jacobian,
258  GeometricSearchData & geom_search_data);
259 
266  void constraintJacobians(SparseMatrix<Number> & jacobian, bool displaced);
267 
269  void setVariableGlobalDoFs(const std::string & var_name);
270  const std::vector<dof_id_type> & getVariableGlobalDoFs() { return _var_all_dof_indices; }
271 
276  void computeJacobian(SparseMatrix<Number> & jacobian,
277  Moose::KernelType kernel_type = Moose::KT_ALL);
278 
287  void computeJacobianBlocks(std::vector<JacobianBlock *> & blocks);
288 
295  Real computeDamping(const NumericVector<Number> & solution, const NumericVector<Number> & update);
296 
300  void computeTimeDerivatives();
301 
305  void onTimestepBegin();
306 
312  virtual void subdomainSetup(SubdomainID subdomain, THREAD_ID tid);
313 
314  virtual void setSolution(const NumericVector<Number> & soln);
315 
319  void updateActive(THREAD_ID tid);
320 
327  virtual void setSolutionUDot(const NumericVector<Number> & udot);
328 
329  virtual NumericVector<Number> & solutionUDot() override;
330  virtual NumericVector<Number> & residualVector(Moose::KernelType type) override;
331  virtual bool hasResidualVector(Moose::KernelType type) const override;
332 
333  virtual const NumericVector<Number> *& currentSolution() override { return _current_solution; }
334 
335  virtual void serializeSolution();
336  virtual NumericVector<Number> & serializedSolution() override;
337 
338  virtual NumericVector<Number> & residualCopy() override;
339  virtual NumericVector<Number> & residualGhosted() override;
340 
341  virtual NumericVector<Number> & RHS() = 0;
342 
343  virtual void augmentSparsity(SparsityPattern::Graph & sparsity,
344  std::vector<dof_id_type> & n_nz,
345  std::vector<dof_id_type> & n_oz) override;
346 
351  void setPreconditioner(std::shared_ptr<MoosePreconditioner> pc);
352 
358  {
359  _use_finite_differenced_preconditioner = use;
360  }
361 
367  void setDecomposition(const std::vector<std::string> & decomposition);
368 
372  void useFieldSplitPreconditioner(bool use = true) { _use_field_split_preconditioner = use; }
373 
383  {
384  _add_implicit_geometric_coupling_entries_to_jacobian = add;
385  }
386 
392  void assembleConstraintsSeparately(bool separately = true)
393  {
394  _assemble_constraints_separately = separately;
395  }
396 
400  void setupDampers();
406  void reinitIncrementAtQpsForDampers(THREAD_ID tid, const std::set<MooseVariable *> & damped_vars);
407 
413  void reinitIncrementAtNodeForDampers(THREAD_ID tid,
414  const std::set<MooseVariable *> & damped_vars);
415 
418  void checkKernelCoverage(const std::set<SubdomainID> & mesh_subdomains) const;
419  bool containsTimeKernel();
421 
425  unsigned int nNonlinearIterations() { return _n_iters; }
426 
430  unsigned int nLinearIterations() { return _n_linear_iters; }
431 
435  unsigned int nResidualEvaluations() { return _n_residual_evaluations; }
436 
440  Real finalNonlinearResidual() { return _final_residual; }
441 
446  Real nonlinearNorm() { return _last_nl_rnorm; }
447 
452  void printAllVariableNorms(bool state) { _print_all_var_norms = state; }
453 
454  void debuggingResiduals(bool state) { _debugging_residuals = state; }
455 
457 
458  void setPredictor(std::shared_ptr<Predictor> predictor);
459  Predictor * getPredictor() { return _predictor.get(); }
460 
461  TimeIntegrator * getTimeIntegrator() { return _time_integrator.get(); }
462 
463  void setPCSide(MooseEnum pcs);
464 
465  Moose::PCSideType getPCSide() { return _pc_side; }
466 
467  void setMooseKSPNormType(MooseEnum kspnorm);
468 
470 
475  bool needMaterialOnSide(BoundaryID bnd_id, THREAD_ID tid) const;
476 
481  bool needMaterialOnSide(SubdomainID subdomain_id, THREAD_ID tid) const;
482 
486  bool doingDG() const;
487 
489 
492  const KernelWarehouse & getKernelWarehouse() { return _kernels; }
493  const KernelWarehouse & getTimeKernelWarehouse() { return _time_kernels; }
494  const KernelWarehouse & getNonTimeKernelWarehouse() { return _non_time_kernels; }
495  const KernelWarehouse & getEigenKernelWarehouse() { return _eigen_kernels; }
496  const KernelWarehouse & getNonEigenKernelWarehouse() { return _non_eigen_kernels; }
497  const MooseObjectWarehouse<DGKernel> & getDGKernelWarehouse() { return _dg_kernels; }
499  {
500  return _interface_kernels;
501  }
502  const MooseObjectWarehouse<DiracKernel> & getDiracKernelWarehouse() { return _dirac_kernels; }
503  const MooseObjectWarehouse<NodalKernel> & getNodalKernelWarehouse(THREAD_ID tid);
504  const MooseObjectWarehouse<IntegratedBC> & getIntegratedBCWarehouse() { return _integrated_bcs; }
506  {
507  return _element_dampers;
508  }
509  const MooseObjectWarehouse<NodalDamper> & getNodalDamperWarehouse() { return _nodal_dampers; }
510  const ConstraintWarehouse & getConstraintWarehouse() { return _constraints; };
512 
516  bool hasSaveIn() const { return _has_save_in || _has_nodalbc_save_in; }
517 
521  bool hasDiagSaveIn() const { return _has_diag_save_in || _has_nodalbc_diag_save_in; }
522 
523  virtual NumericVector<Number> & solution() override { return *_sys.solution; }
524 
525  virtual System & system() override { return _sys; }
526  virtual const System & system() const override { return _sys; }
527 
528  virtual NumericVector<Number> * solutionPreviousNewton() override
529  {
530  return _solution_previous_nl;
531  }
532 
533  virtual void setPreviousNewtonSolution(const NumericVector<Number> & soln);
534 
535 public:
537  System & _sys;
538  // FIXME: make these protected and create getters/setters
544  std::vector<unsigned int> _current_l_its;
545  unsigned int _current_nl_its;
547 
548 protected:
553  void computeResidualInternal(Moose::KernelType type = Moose::KT_ALL);
554 
559  void computeNodalBCs(NumericVector<Number> & residual, Moose::KernelType type = Moose::KT_ALL);
560 
561  void computeJacobianInternal(SparseMatrix<Number> & jacobian, Moose::KernelType kernel_type);
562 
563  void computeDiracContributions(SparseMatrix<Number> * jacobian = NULL);
564 
565  void computeScalarKernelsJacobians(SparseMatrix<Number> & jacobian);
566 
570  void enforceNodalConstraintsResidual(NumericVector<Number> & residual);
571  void enforceNodalConstraintsJacobian(SparseMatrix<Number> & jacobian);
572 
574  const NumericVector<Number> * _current_solution;
576  NumericVector<Number> * _residual_ghosted;
577 
579  NumericVector<Number> & _serialized_solution;
580 
582  NumericVector<Number> * _solution_previous_nl;
583 
585  NumericVector<Number> & _residual_copy;
586 
588  std::shared_ptr<TimeIntegrator> _time_integrator;
590  NumericVector<Number> * _u_dot;
592  Number _du_dot_du;
594  NumericVector<Number> * _Re_time;
596  NumericVector<Number> * _Re_non_time;
597 
610 
612 
619 
622 
625 
628 
631 
634 
636  MooseObjectWarehouseBase<Split> _splits; // use base b/c there are no setup methods
637 
640 
641 protected:
643  NumericVector<Number> * _increment_vec;
645  std::shared_ptr<MoosePreconditioner> _preconditioner;
650 
653 #ifdef LIBMESH_HAVE_PETSC
654  MatFDColoring _fdcoloring;
655 #endif
656  bool _have_decomposition;
659  std::string _decomposition_split;
662 
665 
668 
671 
678 
680  bool _doing_dg;
681 
683  std::vector<std::string> _vecs_to_zero_for_residual;
684 
685  unsigned int _n_iters;
686  unsigned int _n_linear_iters;
687 
690 
692 
694  std::shared_ptr<Predictor> _predictor;
695 
697 
699 
702 
705 
708 
711 
712  void getNodeDofs(dof_id_type node_id, std::vector<dof_id_type> & dofs);
713 
714  std::vector<dof_id_type> _var_all_dof_indices;
715 };
716 
717 #endif /* NONLINEARSYSTEMBASE_H */
virtual const NumericVector< Number > *& currentSolution() override
The solution vector that is currently being operated on.
const ConstraintWarehouse & getConstraintWarehouse()
NumericVector< Number > * _Re_time
residual vector for time contributions
std::vector< dof_id_type > _var_all_dof_indices
Helper class for holding the preconditioning blocks to fill.
unsigned int nLinearIterations()
Return the number of linear iterations.
Base class for deriving general dampers.
Definition: GeneralDamper.h:34
Real nonlinearNorm()
Return the last nonlinear norm.
Base class for split-based preconditioners.
Definition: Split.h:30
unsigned int nNonlinearIterations()
Return the number of non-linear iterations.
const MooseObjectWarehouse< DiracKernel > & getDiracKernelWarehouse()
TimeIntegrator * getTimeIntegrator()
bool _debugging_residuals
true if debugging residuals
NumericVector< Number > * _Re_non_time
residual vector for non-time contributions
bool _assemble_constraints_separately
Whether or not to assemble the residual and Jacobian after the application of each constraint...
const std::vector< dof_id_type > & getVariableGlobalDoFs()
const MooseObjectWarehouse< IntegratedBC > & getIntegratedBCWarehouse()
Moose::PCSideType getPCSide()
void addImplicitGeometricCouplingEntriesToJacobian(bool add=true)
If called with true this will add entries into the jacobian to link together degrees of freedom that ...
Base class for predictors.
Definition: Predictor.h:39
subdomain_id_type SubdomainID
Definition: MooseTypes.h:77
const KernelWarehouse & getNonEigenKernelWarehouse()
bool _has_nodalbc_diag_save_in
If there is a nodal BC having diag_save_in.
const KernelWarehouse & getEigenKernelWarehouse()
MooseObjectWarehouse< NodalBC > _nodal_bcs
KernelWarehouse _eigen_kernels
void debuggingResiduals(bool state)
MooseObjectWarehouseBase< Split > _splits
Decomposition splits.
NumericVector< Number > & _serialized_solution
Serialized version of the solution vector.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const MooseObjectWarehouse< DGKernel > & getDGKernelWarehouse()
bool _has_nodalbc_save_in
If there is a nodal BC having save_in.
MooseObjectWarehouse< NodalDamper > _nodal_dampers
Nodal Dampers for each thread.
void useFieldSplitPreconditioner(bool use=true)
If called with true this system will use a field split preconditioner matrix.
static PetscErrorCode Vec Mat Mat pc
Moose::MooseKSPNormType _ksp_norm
KSP norm type.
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
Base class for a system (of equations)
Definition: SystemBase.h:91
std::shared_ptr< TimeIntegrator > _time_integrator
Time integrator.
Base class for deriving any boundary condition that works at nodes.
Definition: NodalBC.h:38
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
bool _has_save_in
If there is any Kernel or IntegratedBC having save_in.
MooseObjectWarehouse< ScalarKernel > _time_scalar_kernels
Base class for MOOSE preconditioners.
MooseObjectWarehouse< NodalKernel > _nodal_kernels
NodalKernels for each thread.
MooseKSPNormType
Norm type for converge test.
Definition: MooseTypes.h:233
bool _need_residual_ghosted
Whether or not a ghosted copy of the residual needs to be made.
Nonlinear system to be solved.
virtual NumericVector< Number > & solution() override
virtual const System & system() const override
NumericVector< Number > & _residual_copy
Copy of the residual vector.
bool _doing_dg
true if DG is active (optimization reasons)
const MooseObjectWarehouse< ElementDamper > & getElementDamperWarehouse()
void useFiniteDifferencedPreconditioner(bool use=true)
If called with true this system will use a finite differenced form of the Jacobian as the preconditio...
FEProblemBase & _fe_problem
MooseObjectWarehouse< ScalarKernel > _non_time_scalar_kernels
virtual bool computingInitialResidual()
Returns true if this system is currently computing the initial residual for a solve.
std::vector< unsigned int > _current_l_its
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
Definition: DGKernel.h:47
std::shared_ptr< MoosePreconditioner > _preconditioner
Preconditioner.
MooseObjectWarehouse< InterfaceKernel > _interface_kernels
std::vector< std::string > _vecs_to_zero_for_residual
vectors that will be zeroed before a residual computation
An inteface for the _console for outputting to the Console object.
const std::vector< numeric_index_type > & n_nz
bool _use_finite_differenced_preconditioner
Whether or not to use a finite differenced preconditioner.
KernelWarehouse _non_eigen_kernels
MooseObjectWarehouse< ScalarKernel > _scalar_kernels
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
Base class for deriving nodal dampers.
Definition: NodalDamper.h:37
Base class for deriving any boundary condition of a integrated type.
Definition: IntegratedBC.h:33
unsigned int nResidualEvaluations()
Return the total number of residual evaluations done so far in this calculation.
bool _need_serialized_solution
Whether or not a copy of the residual needs to be made.
InterfaceKernel is responsible for interfacing physics across subdomains.
KernelType
Definition: MooseTypes.h:162
MooseObjectWarehouse< DGKernel > _dg_kernels
bool hasSaveIn() const
Weather or not the nonlinear system has save-ins.
ConstraintWarehouse _constraints
Constraints storage object.
Base class for deriving element dampers.
Definition: ElementDamper.h:37
bool _add_implicit_geometric_coupling_entries_to_jacobian
Whether or not to add implicit geometric couplings to the Jacobian for FDP.
KernelWarehouse _kernels
bool _compute_initial_residual_before_preset_bcs
MooseObjectWarehouse< DiracKernel > _dirac_kernels
Dirac Kernel storage for each thread.
void printAllVariableNorms(bool state)
Force the printing of all variable norms after each solve.
MatType type
unsigned int _n_residual_evaluations
Total number of residual evaluations that have been performed.
bool _has_diag_save_in
If there is any Kernel or IntegratedBC having diag_save_in.
Base class for time integrators.
bool hasDiagSaveIn() const
Weather or not the nonlinear system has diagonal Jacobian save-ins.
const MooseObjectWarehouse< NodalDamper > & getNodalDamperWarehouse()
bool _need_residual_copy
Whether or not a copy of the residual needs to be made.
NumericVector< Number > * _solution_previous_nl
Solution vector of the previous nonlinear iterate.
MooseObjectWarehouse< ElementDamper > _element_dampers
Element Dampers for each thread.
virtual System & system() override
Get the reference to the libMesh system.
Real finalNonlinearResidual()
Return the final nonlinear residual.
Moose::PCSideType _pc_side
Preconditioning side.
Holds kernels and provides some services.
bool haveFiniteDifferencedPreconditioner()
const NumericVector< Number > * _current_solution
solution vector from nonlinear solver
PCSideType
Preconditioning side.
Definition: MooseTypes.h:222
virtual NumericVector< Number > * solutionPreviousNewton() override
NumericVector< Number > * _residual_ghosted
ghosted form of the residual
KernelWarehouse _time_kernels
Moose::MooseKSPNormType getMooseKSPNormType()
NumericVector< Number > * _increment_vec
increment vector
const KernelWarehouse & getTimeKernelWarehouse()
bool _use_field_split_preconditioner
Whether or not to use a FieldSplitPreconditioner matrix based on the decomposition.
MooseObjectWarehouse< IntegratedBC > _integrated_bcs
void assembleConstraintsSeparately(bool separately=true)
Indicates whether to assemble residual and Jacobian after each constraint application.
std::string _decomposition_split
Name of the top-level split of the decomposition.
const KernelWarehouse & getNonTimeKernelWarehouse()
Base class for creating new types of boundary conditions.
Definition: NodalKernel.h:50
std::shared_ptr< Predictor > _predictor
If predictor is active, this is non-NULL.
unsigned int _num_residual_evaluations
MooseObjectWarehouse< GeneralDamper > _general_dampers
General Dampers.
const MooseObjectWarehouse< InterfaceKernel > & getInterfaceKernelWarehouse()
A DiracKernel is used when you need to add contributions to the residual by means of multiplying some...
Definition: DiracKernel.h:50
const std::vector< numeric_index_type > & n_oz
Warehouse for storing constraints.
virtual void addEigenKernels(std::shared_ptr< KernelBase >, THREAD_ID)
boundary_id_type BoundaryID
Definition: MooseTypes.h:75
unsigned int THREAD_ID
Definition: MooseTypes.h:79
KernelWarehouse _non_time_kernels
const KernelWarehouse & getKernelWarehouse()
Access functions to Warehouses from outside NonlinearSystemBase.
NumericVector< Number > * _u_dot
solution vector for u^dot
MooseObjectWarehouse< PresetNodalBC > _preset_nodal_bcs