www.mooseframework.org
NonlinearSystemBase.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "SystemBase.h"
13 #include "ConstraintWarehouse.h"
14 #include "MooseObjectWarehouse.h"
16 #include "PerfGraphInterface.h"
17 #include "ComputeMortarFunctor.h"
18 #include "MooseHashing.h"
19 
20 #include "libmesh/transient_system.h"
21 #include "libmesh/nonlinear_implicit_system.h"
22 
23 // Forward declarations
24 class FEProblemBase;
26 class JacobianBlock;
27 class TimeIntegrator;
28 class Predictor;
29 class ElementDamper;
30 class NodalDamper;
31 class GeneralDamper;
33 class IntegratedBCBase;
34 class NodalBCBase;
35 class PresetNodalBC;
36 template <ComputeStage>
37 class ADPresetNodalBC;
38 class DGKernelBase;
40 class ScalarKernel;
41 class DiracKernel;
42 class NodalKernel;
43 class Split;
44 class KernelBase;
45 class BoundaryCondition;
46 
47 // libMesh forward declarations
48 namespace libMesh
49 {
50 template <typename T>
51 class NumericVector;
52 template <typename T>
53 class SparseMatrix;
54 } // namespace libMesh
55 
63  public PerfGraphInterface
64 {
65 public:
66  NonlinearSystemBase(FEProblemBase & problem, System & sys, const std::string & name);
67  virtual ~NonlinearSystemBase();
68 
69  virtual void init() override;
70 
74  void turnOffJacobian();
75 
76  virtual void addExtraVectors() override;
77  virtual void solve() override = 0;
78  virtual void restoreSolutions() override;
79 
83  virtual void stopSolve() = 0;
84 
85  virtual NonlinearSolver<Number> * nonlinearSolver() = 0;
86 
87  virtual unsigned int getCurrentNonlinearIterationNumber() = 0;
88 
94 
95  // Setup Functions ////
96  virtual void initialSetup();
97  virtual void timestepSetup();
98 
99  virtual void setupFiniteDifferencedPreconditioner() = 0;
101 
103  {
105  }
107 
112  virtual bool converged() = 0;
113 
120  void addTimeIntegrator(const std::string & type,
121  const std::string & name,
122  InputParameters parameters) override;
124 
129  void addDotVectors();
130 
137  virtual void
138  addKernel(const std::string & kernel_name, const std::string & name, InputParameters parameters);
139 
146  virtual void addNodalKernel(const std::string & kernel_name,
147  const std::string & name,
148  InputParameters parameters);
149 
156  void addScalarKernel(const std::string & kernel_name,
157  const std::string & name,
158  InputParameters parameters);
159 
166  void addBoundaryCondition(const std::string & bc_name,
167  const std::string & name,
168  InputParameters parameters);
169 
176  void
177  addConstraint(const std::string & c_name, const std::string & name, InputParameters parameters);
178 
185  void addDiracKernel(const std::string & kernel_name,
186  const std::string & name,
187  InputParameters parameters);
188 
195  void
196  addDGKernel(std::string dg_kernel_name, const std::string & name, InputParameters parameters);
197 
204  void addInterfaceKernel(std::string interface_kernel_name,
205  const std::string & name,
206  InputParameters parameters);
207 
214  void
215  addDamper(const std::string & damper_name, const std::string & name, InputParameters parameters);
216 
223  void
224  addSplit(const std::string & split_name, const std::string & name, InputParameters parameters);
225 
230  std::shared_ptr<Split> getSplit(const std::string & name);
231 
232  void zeroVectorForResidual(const std::string & vector_name);
233 
234  void setInitialSolution();
235 
239  void setConstraintSlaveValues(NumericVector<Number> & solution, bool displaced);
240 
248  void constraintResiduals(NumericVector<Number> & residual, bool displaced);
249 
255  void computeResidualTag(NumericVector<Number> & residual, TagID tag_id);
256 
260  void computeResidualTags(const std::set<TagID> & tags);
261 
265  void computeResidual(NumericVector<Number> & residual, TagID tag_id);
266 
270  void
272  std::map<dof_id_type, std::vector<dof_id_type>> & graph);
273 
280 
287  void constraintJacobians(bool displaced);
288 
290  void setVariableGlobalDoFs(const std::string & var_name);
291  const std::vector<dof_id_type> & getVariableGlobalDoFs() { return _var_all_dof_indices; }
292 
296  void computeJacobianTags(const std::set<TagID> & tags);
297 
301  void computeJacobian(SparseMatrix<Number> & jacobian, const std::set<TagID> & tags);
302 
306  void computeJacobian(SparseMatrix<Number> & jacobian);
307 
316  void computeJacobianBlocks(std::vector<JacobianBlock *> & blocks);
317 
318  void computeJacobianBlocks(std::vector<JacobianBlock *> & blocks, const std::set<TagID> & tags);
319 
326  Real computeDamping(const NumericVector<Number> & solution, const NumericVector<Number> & update);
327 
331  void computeTimeDerivatives();
332 
336  void onTimestepBegin();
337 
343  virtual void subdomainSetup(SubdomainID subdomain, THREAD_ID tid);
344 
345  virtual void setSolution(const NumericVector<Number> & soln);
346 
350  void updateActive(THREAD_ID tid);
351 
358  virtual void setSolutionUDot(const NumericVector<Number> & udot);
359 
366  virtual void setSolutionUDotDot(const NumericVector<Number> & udotdot);
367 
368  NumericVector<Number> * solutionUDot() override { return _u_dot; }
369  NumericVector<Number> * solutionUDotDot() override { return _u_dotdot; }
370  NumericVector<Number> * solutionUDotOld() override { return _u_dot_old; }
371  NumericVector<Number> * solutionUDotDotOld() override { return _u_dotdot_old; }
372  const NumericVector<Number> * solutionUDot() const override { return _u_dot; }
373  const NumericVector<Number> * solutionUDotDot() const override { return _u_dotdot; }
374  const NumericVector<Number> * solutionUDotOld() const override { return _u_dot_old; }
375  const NumericVector<Number> * solutionUDotDotOld() const override { return _u_dotdot_old; }
376 
380  NumericVector<Number> & getResidualTimeVector();
381 
385  NumericVector<Number> & getResidualNonTimeVector();
386 
390  NumericVector<Number> & residualVector(TagID tag);
391 
392  const NumericVector<Number> * const & currentSolution() const override
393  {
394  return _current_solution;
395  }
396 
397  virtual void serializeSolution();
398  virtual NumericVector<Number> & serializedSolution() override;
399 
400  virtual NumericVector<Number> & residualCopy() override;
401  virtual NumericVector<Number> & residualGhosted() override;
402 
403  virtual NumericVector<Number> & RHS() = 0;
404 
405  virtual void augmentSparsity(SparsityPattern::Graph & sparsity,
406  std::vector<dof_id_type> & n_nz,
407  std::vector<dof_id_type> & n_oz) override;
408 
413  void setPreconditioner(std::shared_ptr<MoosePreconditioner> pc);
414  MoosePreconditioner const * getPreconditioner() const;
415 
421  {
423  }
424 
430  void setDecomposition(const std::vector<std::string> & decomposition);
431 
436 
446  {
448  }
449 
455  void assembleConstraintsSeparately(bool separately = true)
456  {
458  }
459 
463  void setupDampers();
469  void reinitIncrementAtQpsForDampers(THREAD_ID tid, const std::set<MooseVariable *> & damped_vars);
470 
477  const std::set<MooseVariable *> & damped_vars);
478 
481  void checkKernelCoverage(const std::set<SubdomainID> & mesh_subdomains) const;
482  bool containsTimeKernel();
484 
488  unsigned int nNonlinearIterations() const { return _n_iters; }
489 
493  unsigned int nLinearIterations() const { return _n_linear_iters; }
494 
498  unsigned int nResidualEvaluations() const { return _n_residual_evaluations; }
499 
503  Real finalNonlinearResidual() const { return _final_residual; }
504 
509  Real nonlinearNorm() const { return _last_nl_rnorm; }
510 
515  void printAllVariableNorms(bool state) { _print_all_var_norms = state; }
516 
517  void debuggingResiduals(bool state) { _debugging_residuals = state; }
518 
520 
521  void setPredictor(std::shared_ptr<Predictor> predictor);
522  Predictor * getPredictor() { return _predictor.get(); }
523 
524  void setPCSide(MooseEnum pcs);
525 
527 
528  void setMooseKSPNormType(MooseEnum kspnorm);
529 
531 
536  bool needBoundaryMaterialOnSide(BoundaryID bnd_id, THREAD_ID tid) const;
537 
542  bool needSubdomainMaterialOnSide(SubdomainID subdomain_id, THREAD_ID tid) const;
543 
547  bool doingDG() const;
548 
550 
555  {
556  return _ad_jacobian_kernels;
557  }
560  {
561  return _interface_kernels;
562  }
566  {
567  return _element_dampers;
568  }
570  {
571  return _nodal_dampers;
572  }
575 
579  bool hasSaveIn() const { return _has_save_in || _has_nodalbc_save_in; }
580 
585 
586  NumericVector<Number> & solution() override { return *_sys.solution; }
587  const NumericVector<Number> & solution() const override { return *_sys.solution; }
588 
589  virtual System & system() override { return _sys; }
590  virtual const System & system() const override { return _sys; }
591 
592  NumericVector<Number> * solutionPreviousNewton() override { return _solution_previous_nl; }
593  const NumericVector<Number> * solutionPreviousNewton() const override
594  {
595  return _solution_previous_nl;
596  }
597 
598  virtual void setSolutionUDotOld(const NumericVector<Number> & u_dot_old);
599 
600  virtual void setSolutionUDotDotOld(const NumericVector<Number> & u_dotdot_old);
601 
602  virtual void setPreviousNewtonSolution(const NumericVector<Number> & soln);
603 
604  virtual TagID timeVectorTag() override { return _Re_time_tag; }
605 
606  virtual TagID nonTimeVectorTag() override { return _Re_non_time_tag; }
607 
608  virtual TagID residualVectorTag() override { return _Re_tag; }
609 
610  virtual TagID systemMatrixTag() override { return _Ke_system_tag; }
611 
612 public:
614  System & _sys;
615  // FIXME: make these protected and create getters/setters
621  std::vector<unsigned int> _current_l_its;
622  unsigned int _current_nl_its;
624 
625 protected:
630  void computeResidualInternal(const std::set<TagID> & tags);
631 
636  void computeNodalBCs(NumericVector<Number> & residual);
637 
641  void computeNodalBCs(NumericVector<Number> & residual, const std::set<TagID> & tags);
642 
646  void computeNodalBCs(const std::set<TagID> & tags);
647 
651  void computeJacobianInternal(const std::set<TagID> & tags);
652 
653  void computeDiracContributions(const std::set<TagID> & tags, bool is_jacobian);
654 
655  void computeScalarKernelsJacobians(const std::set<TagID> & tags);
656 
660  void enforceNodalConstraintsResidual(NumericVector<Number> & residual);
662 
666  void mortarResidualConstraints(bool displaced);
667 
671  void mortarJacobianConstraints(bool displaced);
672 
673 protected:
675  const NumericVector<Number> * _current_solution;
677  NumericVector<Number> * _residual_ghosted;
678 
680  NumericVector<Number> & _serialized_solution;
681 
683  NumericVector<Number> * _solution_previous_nl;
684 
686  NumericVector<Number> & _residual_copy;
687 
689  NumericVector<Number> * _u_dot;
691  NumericVector<Number> * _u_dotdot;
692 
694  NumericVector<Number> * _u_dot_old;
696  NumericVector<Number> * _u_dotdot_old;
697 
699  Number _du_dot_du;
702 
705 
707  std::set<TagID> _nl_vector_tags;
708 
710  std::set<TagID> _nl_matrix_tags;
711 
713  NumericVector<Number> * _Re_time;
714 
718  NumericVector<Number> * _Re_non_time;
719 
722 
725 
728 
736 
738 
746 
749 
752 
755 
758 
761 
763  MooseObjectWarehouseBase<Split> _splits; // use base b/c there are no setup methods
764 
767 
769  NumericVector<Number> * _increment_vec;
771  std::shared_ptr<MoosePreconditioner> _preconditioner;
776 
779 #ifdef LIBMESH_HAVE_PETSC
780  MatFDColoring _fdcoloring;
781 #endif
782  bool _have_decomposition;
785  std::string _decomposition_split;
788 
791 
794 
797 
804 
806  bool _doing_dg;
807 
809  std::vector<std::string> _vecs_to_zero_for_residual;
810 
811  unsigned int _n_iters;
812  unsigned int _n_linear_iters;
813 
816 
818 
820  std::shared_ptr<Predictor> _predictor;
821 
823 
825 
828 
831 
834 
837 
838  void getNodeDofs(dof_id_type node_id, std::vector<dof_id_type> & dofs);
839 
840  std::vector<dof_id_type> _var_all_dof_indices;
841 
854 
855 private:
857  std::unordered_map<std::pair<BoundaryID, BoundaryID>,
860 
862  std::unordered_map<std::pair<BoundaryID, BoundaryID>,
865 
867  std::unordered_map<std::pair<BoundaryID, BoundaryID>,
870 
872  std::unordered_map<std::pair<BoundaryID, BoundaryID>,
875 };
virtual void setSolutionUDotDotOld(const NumericVector< Number > &u_dotdot_old)
const ConstraintWarehouse & getConstraintWarehouse() const
NumericVector< Number > & getResidualTimeVector()
Return a numeric vector that is associated with the time tag.
NumericVector< Number > * _Re_time
residual vector for time contributions
void computeJacobianBlocks(std::vector< JacobianBlock *> &blocks)
Computes several Jacobian blocks simultaneously, summing their contributions into smaller preconditio...
std::vector< dof_id_type > _var_all_dof_indices
unsigned int nLinearIterations() const
Return the number of linear iterations.
Helper class for holding the preconditioning blocks to fill.
TagID _Re_time_tag
Tag for time contribution residual.
virtual void setSolutionUDotDot(const NumericVector< Number > &udotdot)
Set transient term used by residual and Jacobian evaluation.
void reinitIncrementAtNodeForDampers(THREAD_ID tid, const std::set< MooseVariable *> &damped_vars)
Compute the incremental change in variables at nodes for dampers.
Base class for deriving general dampers.
Definition: GeneralDamper.h:26
std::unordered_map< std::pair< BoundaryID, BoundaryID >, ComputeMortarFunctor< ComputeStage::RESIDUAL > > _undisplaced_mortar_residual_functors
Functors for computing residuals from undisplaced mortar constraints.
virtual bool converged()=0
Returns the convergence state.
MoosePreconditioner const * getPreconditioner() const
void setupDampers()
Setup damping stuff (called before we actually start)
void zeroVectorForResidual(const std::string &vector_name)
Base class for split-based preconditioners.
Definition: Split.h:29
virtual void timestepSetup()
void setConstraintSlaveValues(NumericVector< Number > &solution, bool displaced)
Sets the value of constrained variables in the solution vector.
bool _debugging_residuals
true if debugging residuals
NumericVector< Number > * _Re_non_time
residual vector for non-time contributions
unsigned int TagID
Definition: MooseTypes.h:162
virtual NonlinearSolver< Number > * nonlinearSolver()=0
Real computeDamping(const NumericVector< Number > &solution, const NumericVector< Number > &update)
Compute damping.
virtual void setPreviousNewtonSolution(const NumericVector< Number > &soln)
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< ElementDamper > & getElementDamperWarehouse() const
MooseObjectTagWarehouse< DGKernelBase > _dg_kernels
virtual void init() override
Initialize the system.
void reinitIncrementAtQpsForDampers(THREAD_ID tid, const std::set< MooseVariable *> &damped_vars)
Compute the incremental change in variables at QPs for dampers.
Moose::PCSideType getPCSide()
MooseObjectTagWarehouse< KernelBase > & getADJacobianKernelWarehouse()
void addImplicitGeometricCouplingEntriesToJacobian(bool add=true)
If called with true this will add entries into the jacobian to link together degrees of freedom that ...
virtual void setSolution(const NumericVector< Number > &soln)
bool haveFieldSplitPreconditioner() const
Base class for predictors.
Definition: Predictor.h:33
void addDiracKernel(const std::string &kernel_name, const std::string &name, InputParameters parameters)
Adds a Dirac kernel.
bool _has_nodalbc_diag_save_in
If there is a nodal BC having diag_save_in.
NumericVector< Number > * _u_dotdot_old
old solution vector for u^dotdot
void checkKernelCoverage(const std::set< SubdomainID > &mesh_subdomains) const
void addDotVectors()
Add u_dot, u_dotdot, u_dot_old and u_dotdot_old vectors if requested by the time integrator.
void getNodeDofs(dof_id_type node_id, std::vector< dof_id_type > &dofs)
MooseObjectTagWarehouse< KernelBase > & getKernelWarehouse()
Access functions to Warehouses from outside NonlinearSystemBase.
std::set< TagID > _nl_vector_tags
Vector tags to temporarily store all tags associated with the current system.
MooseObjectTagWarehouse< NodalBCBase > _nodal_bcs
void addBoundaryCondition(const std::string &bc_name, const std::string &name, InputParameters parameters)
Adds a boundary condition.
void setDecomposition(const std::vector< std::string > &decomposition)
If called with a single string, it is used as the name of a the top-level decomposition split...
bool needSubdomainMaterialOnSide(SubdomainID subdomain_id, THREAD_ID tid) const
Indicates whether this system needs material properties on internal sides.
std::unordered_map< std::pair< BoundaryID, BoundaryID >, ComputeMortarFunctor< ComputeStage::JACOBIAN > > _undisplaced_mortar_jacobian_functors
Functors for computing jacobians from undisplaced mortar constraints.
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...
bool _has_nodalbc_save_in
If there is a nodal BC having save_in.
MooseObjectWarehouse< NodalDamper > _nodal_dampers
Nodal Dampers for each thread.
void computeTimeDerivatives()
Computes the time derivative vector.
void useFieldSplitPreconditioner(bool use=true)
If called with true this system will use a field split preconditioner matrix.
MooseObjectTagWarehouse< IntegratedBCBase > & getIntegratedBCWarehouse()
virtual void setSolutionUDotOld(const NumericVector< Number > &u_dot_old)
bool hasDiagSaveIn() const
Weather or not the nonlinear system has diagonal Jacobian save-ins.
static PetscErrorCode Vec Mat Mat pc
void findImplicitGeometricCouplingEntries(GeometricSearchData &geom_search_data, std::map< dof_id_type, std::vector< dof_id_type >> &graph)
Finds the implicit sparsity graph between geometrically related dofs.
void setPCSide(MooseEnum pcs)
Moose::MooseKSPNormType _ksp_norm
KSP norm type.
const NumericVector< Number > * solutionUDotOld() const override
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
Base class for a system (of equations)
Definition: SystemBase.h:92
void computeJacobian(SparseMatrix< Number > &jacobian, const std::set< TagID > &tags)
Associate jacobian to systemMatrixTag, and then form a matrix for all the tags.
virtual void addExtraVectors() override
Method called during initialSetup to add extra system vector if they are required by the simulation...
Real finalNonlinearResidual() const
Return the final nonlinear residual.
void computeResidualTags(const std::set< TagID > &tags)
Form multiple tag-associated residual vectors for all the given tags.
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
bool _has_save_in
If there is any Kernel or IntegratedBC having save_in.
TagID _Ke_non_time_tag
Tag for non-time contribution Jacobian.
TagID _Ke_system_tag
Tag for system contribution Jacobian.
void addConstraint(const std::string &c_name, const std::string &name, InputParameters parameters)
Adds a Constraint.
virtual TagID systemMatrixTag() override
Return the Matrix Tag ID for System.
bool _have_decomposition
Whether or not the system can be decomposed into splits.
unsigned int PerfID
Definition: MooseTypes.h:163
Serves as a base class for DGKernel and ADDGKernel.
Definition: DGKernelBase.h:46
void constraintResiduals(NumericVector< Number > &residual, bool displaced)
Add residual contributions from Constraints.
Base class for MOOSE preconditioners.
MooseObjectTagWarehouse< KernelBase > _ad_jacobian_kernels
Base class for automatic differentiation nodal BCs that (pre)set the solution vector entries...
virtual void update()
Update the system (doing libMesh magic)
Definition: SystemBase.C:1020
MooseObjectTagWarehouse< DiracKernel > & getDiracKernelWarehouse()
void computeResidual(NumericVector< Number > &residual, TagID tag_id)
Form a residual vector for a given tag.
MooseKSPNormType
Norm type for converge test.
Definition: MooseTypes.h:577
MooseObjectTagWarehouse< InterfaceKernelBase > & getInterfaceKernelWarehouse()
bool _need_residual_ghosted
Whether or not a ghosted copy of the residual needs to be made.
Nonlinear system to be solved.
virtual const std::string & name() const
Definition: SystemBase.C:1088
void onTimestepBegin()
Called at the beginning of the time step.
virtual TagID residualVectorTag() override
NumericVector< Number > * solutionPreviousNewton() override
TagID _Re_non_time_tag
Tag for non-time contribution residual.
std::set< TagID > _nl_matrix_tags
Matrix tags to temporarily store all tags associated with the current system.
NumericVector< Number > & _residual_copy
Copy of the residual vector.
void setPredictor(std::shared_ptr< Predictor > predictor)
virtual void setupFiniteDifferencedPreconditioner()=0
bool _doing_dg
true if DG is active (optimization reasons)
void useFiniteDifferencedPreconditioner(bool use=true)
If called with true this system will use a finite differenced form of the Jacobian as the preconditio...
void addDamper(const std::string &damper_name, const std::string &name, InputParameters parameters)
Adds a damper.
bool needBoundaryMaterialOnSide(BoundaryID bnd_id, THREAD_ID tid) const
Indicated whether this system needs material properties on boundaries.
FEProblemBase & _fe_problem
virtual void restoreSolutions() override
Restore current solutions (call after your solve failed)
virtual void stopSolve()=0
Quit the current solve as soon as possible.
const NumericVector< Number > * solutionUDot() const override
virtual bool computingInitialResidual()
Returns true if this system is currently computing the initial residual for a solve.
std::vector< unsigned int > _current_l_its
PerfID _compute_residual_tags_timer
Timers.
This is the common base class for the two main kernel types implemented in MOOSE, EigenKernel and Ker...
Definition: KernelBase.h:44
NumericVector< Number > * solutionUDotDot() override
void computeDiracContributions(const std::set< TagID > &tags, bool is_jacobian)
void updateActive(THREAD_ID tid)
Update active objects of Warehouses owned by NonlinearSystemBase.
const NumericVector< Number > * solutionUDotDotOld() const override
void setVariableGlobalDoFs(const std::string &var_name)
set all the global dof indices for a nonlinear variable
virtual TagID nonTimeVectorTag() override
virtual const System & system() const override
boundary_id_type BoundaryID
virtual void addTimeIntegrator(const std::string &, const std::string &, InputParameters)
Definition: SystemBase.h:721
virtual void serializeSolution()
void computeJacobianTags(const std::set< TagID > &tags)
Computes multiple (tag associated) Jacobian matricese.
std::shared_ptr< MoosePreconditioner > _preconditioner
Preconditioner.
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
std::shared_ptr< Split > getSplit(const std::string &name)
Retrieves a split by name.
bool haveFiniteDifferencedPreconditioner() const
MooseObjectTagWarehouse< ScalarKernel > _scalar_kernels
void mortarJacobianConstraints(bool displaced)
Do mortar constraint jacobian computation.
void addImplicitGeometricCouplingEntries(GeometricSearchData &geom_search_data)
Adds entries to the Jacobian in the correct positions for couplings coming from dofs being coupled th...
bool _use_finite_differenced_preconditioner
Whether or not to use a finite differenced preconditioner.
MooseObjectTagWarehouse< KernelBase > _kernels
virtual void initialSetup()
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
Base class for deriving nodal dampers.
Definition: NodalDamper.h:32
Base class for nodal BCs that (pre)set the solution vector entries.
Definition: PresetNodalBC.h:22
NumericVector< Number > * solutionUDot() override
void addSplit(const std::string &split_name, const std::string &name, InputParameters parameters)
Adds a split.
subdomain_id_type SubdomainID
virtual TagID timeVectorTag() override
Ideally, we should not need this API.
bool _need_serialized_solution
Whether or not a copy of the residual needs to be made.
Base class for deriving any boundary condition that works at nodes.
Definition: NodalBCBase.h:32
NonlinearSystemBase(FEProblemBase &problem, System &sys, const std::string &name)
void computeResidualInternal(const std::set< TagID > &tags)
Compute the residual for a given tag.
Interface for objects that needs transient capabilities.
virtual void solve() override=0
Solve the system (using libMesh magic)
const NumericVector< Number > & solution() const override
ConstraintWarehouse _constraints
Constraints storage object.
void turnOffJacobian()
Turn off the Jacobian (must be called before equation system initialization)
MooseObjectTagWarehouse< NodalKernel > _nodal_kernels
NodalKernels for each thread.
virtual void addNodalKernel(const std::string &kernel_name, const std::string &name, InputParameters parameters)
Adds a NodalKernel.
NumericVector< Number > & solution() override
const NumericVector< Number > *const & currentSolution() const override
The solution vector that is currently being operated on.
Base class for deriving element dampers.
Definition: ElementDamper.h:33
NumericVector< Number > * solutionUDotOld() override
bool _add_implicit_geometric_coupling_entries_to_jacobian
Whether or not to add implicit geometric couplings to the Jacobian for FDP.
const NumericVector< Number > * solutionPreviousNewton() const override
const NumericVector< Number > * solutionUDotDot() const override
bool _compute_initial_residual_before_preset_bcs
Base class for creating new types of boundary conditions.
unsigned int nNonlinearIterations() const
Return the number of non-linear iterations.
void printAllVariableNorms(bool state)
Force the printing of all variable norms after each solve.
void addDGKernel(std::string dg_kernel_name, const std::string &name, InputParameters parameters)
Adds a DG kernel.
virtual NumericVector< Number > & RHS()=0
MatType type
unsigned int _n_residual_evaluations
Total number of residual evaluations that have been performed.
virtual void addKernel(const std::string &kernel_name, const std::string &name, InputParameters parameters)
Adds a kernel.
NumericVector< Number > * _u_dot_old
old solution vector for u^dot
void computeJacobianInternal(const std::set< TagID > &tags)
Form multiple matrices for all the tags.
bool _has_diag_save_in
If there is any Kernel or IntegratedBC having diag_save_in.
MooseObjectWarehouse< ADPresetNodalBC< RESIDUAL > > _ad_preset_nodal_bcs
Base class for time integrators.
void addScalarKernel(const std::string &kernel_name, const std::string &name, InputParameters parameters)
Adds a scalar kernel.
void constraintJacobians(bool displaced)
Add jacobian contributions from Constraints.
bool hasSaveIn() const
Weather or not the nonlinear system has save-ins.
bool _need_residual_copy
Whether or not a copy of the residual needs to be made.
virtual unsigned int getCurrentNonlinearIterationNumber()=0
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.
NumericVector< Number > * _u_dotdot
solution vector for u^dotdot
virtual void setSolutionUDot(const NumericVector< Number > &udot)
Set transient term used by residual and Jacobian evaluation.
Moose::PCSideType _pc_side
Preconditioning side.
Real nonlinearNorm() const
Return the last nonlinear norm.
void computeScalarKernelsJacobians(const std::set< TagID > &tags)
const NumericVector< Number > * _current_solution
solution vector from nonlinear solver
PCSideType
Preconditioning side.
Definition: MooseTypes.h:566
Base class for deriving any boundary condition of a integrated type.
NumericVector< Number > * _residual_ghosted
ghosted form of the residual
TagID _Re_tag
Used for the residual vector from PETSc.
void computeNodalBCs(NumericVector< Number > &residual)
Enforces nodal boundary conditions.
Moose::MooseKSPNormType getMooseKSPNormType()
NumericVector< Number > * _increment_vec
increment vector
InterfaceKernelBase is the base class for all InterfaceKernel type classes.
bool doingDG() const
Getter for _doing_dg.
virtual NumericVector< Number > & serializedSolution() override
Returns a reference to a serialized version of the solution vector for this subproblem.
void computeResidualTag(NumericVector< Number > &residual, TagID tag_id)
Computes residual for a given tag.
void addInterfaceKernel(std::string interface_kernel_name, const std::string &name, InputParameters parameters)
Adds an interface kernel.
std::unordered_map< std::pair< BoundaryID, BoundaryID >, ComputeMortarFunctor< ComputeStage::RESIDUAL > > _displaced_mortar_residual_functors
Functors for computing residuals from displaced mortar constraints.
bool _use_field_split_preconditioner
Whether or not to use a FieldSplitPreconditioner matrix based on the decomposition.
NumericVector< Number > * solutionUDotDotOld() override
void addTimeIntegrator(const std::string &type, const std::string &name, InputParameters parameters) override
Add a time integrator.
void assembleConstraintsSeparately(bool separately=true)
Indicates whether to assemble residual and Jacobian after each constraint application.
const MooseObjectWarehouse< NodalDamper > & getNodalDamperWarehouse() const
std::unordered_map< std::pair< BoundaryID, BoundaryID >, ComputeMortarFunctor< ComputeStage::JACOBIAN > > _displaced_mortar_jacobian_functors
Functors for computing jacobians from displaced mortar constraints.
unsigned int nResidualEvaluations() const
Return the total number of residual evaluations done so far in this calculation.
std::string _decomposition_split
Name of the top-level split of the decomposition.
MooseObjectTagWarehouse< InterfaceKernelBase > _interface_kernels
MooseObjectTagWarehouse< DGKernelBase > & getDGKernelWarehouse()
NumericVector< Number > & residualVector(TagID tag)
Return a residual vector that is associated with the residual tag.
NumericVector< Number > & getResidualNonTimeVector()
Return a numeric vector that is associated with the nontime tag.
Base class for creating new types of boundary conditions.
Definition: NodalKernel.h:48
std::shared_ptr< Predictor > _predictor
If predictor is active, this is non-NULL.
void enforceNodalConstraintsResidual(NumericVector< Number > &residual)
Enforce nodal constraints.
unsigned int _num_residual_evaluations
MooseObjectWarehouse< GeneralDamper > _general_dampers
General Dampers.
virtual void subdomainSetup(SubdomainID subdomain, THREAD_ID tid)
Called from assembling when we hit a new subdomain.
MooseObjectTagWarehouse< IntegratedBCBase > _integrated_bcs
A DiracKernel is used when you need to add contributions to the residual by means of multiplying some...
Definition: DiracKernel.h:45
virtual NumericVector< Number > & residualCopy() override
const std::vector< numeric_index_type > & n_oz
virtual NumericVector< Number > & residualGhosted() override
MooseObjectTagWarehouse< DiracKernel > _dirac_kernels
Dirac Kernel storage for each thread.
void mortarResidualConstraints(bool displaced)
Do mortar constraint residual computation.
void setMooseKSPNormType(MooseEnum kspnorm)
Warehouse for storing constraints.
unsigned int THREAD_ID
Definition: MooseTypes.h:161
NumericVector< Number > * _u_dot
solution vector for u^dot
virtual void augmentSparsity(SparsityPattern::Graph &sparsity, std::vector< dof_id_type > &n_nz, std::vector< dof_id_type > &n_oz) override
Will modify the sparsity pattern to add logical geometric connections.
MooseObjectWarehouse< PresetNodalBC > _preset_nodal_bcs
void setPreconditioner(std::shared_ptr< MoosePreconditioner > pc)
Sets a preconditioner.