|
OOFEM 3.0
|
#include <nrsolver.h>
Public Member Functions | |
| NRSolver (Domain *d, EngngModel *m) | |
| virtual | ~NRSolver () |
| ConvergedReason | solve (SparseMtrx &k, FloatArray &R, FloatArray *R0, FloatArray &X, FloatArray &dX, FloatArray &F, const FloatArray &internalForcesEBENorm, double &l, referenceLoadInputModeType rlm, int &nite, TimeStep *) override |
| void | printState (FILE *outputStream) override |
| void | initializeFrom (InputRecord &ir) override |
| const char * | giveClassName () const override |
| virtual const char * | giveInputRecordName () const |
| void | setDomain (Domain *d) override |
| void | reinitialize () override |
| SparseLinearSystemNM * | giveLinearSolver () override |
| Public Member Functions inherited from oofem::SparseNonLinearSystemNM | |
| SparseNonLinearSystemNM (Domain *d, EngngModel *m) | |
| Constructor. | |
| virtual | ~SparseNonLinearSystemNM () |
| Destructor. | |
| virtual double | giveCurrentStepLength () |
| virtual void | setStepLength (double s) |
| virtual bool | referenceLoad () const |
| void | initializeFrom (InputRecord &ir) override |
| virtual void | convertPertMap () |
| virtual void | applyPerturbation (FloatArray *displacement) |
| std::string | errorInfo (const char *func) const |
| Error printing helper. | |
| Public Member Functions inherited from oofem::NumericalMethod | |
| NumericalMethod (Domain *d, EngngModel *m) | |
| virtual | ~NumericalMethod () |
| Destructor. | |
| EngngModel * | giveEngngModel () |
| virtual void | saveContext (DataStream &stream, ContextMode mode) |
| virtual void | restoreContext (DataStream &stream, ContextMode mode) |
Protected Types | |
| enum | nrsolver_ModeType { nrsolverModifiedNRM , nrsolverFullNRM , nrsolverAccelNRM } |
Protected Member Functions | |
| LineSearchNM * | giveLineSearchSolver () |
| Constructs and returns a line search solver. | |
| void | initPrescribedEqs () |
| Initiates prescribed equations. | |
| void | applyConstraintsToStiffness (SparseMtrx &k) |
| void | applyConstraintsToLoadIncrement (int nite, const SparseMtrx &k, FloatArray &R, referenceLoadInputModeType rlm, TimeStep *tStep) |
| bool | checkConvergence (FloatArray &RT, FloatArray &F, FloatArray &rhs, FloatArray &ddX, FloatArray &X, double RRT, const FloatArray &internalForcesEBENorm, int nite, bool &errorOutOfRange) |
Protected Attributes | |
| int | nsmax |
| int | minIterations |
| double | minStepLength |
| nrsolver_ModeType | NR_Mode |
| nrsolver_ModeType | NR_OldMode |
| int | NR_ModeTick |
| int | MANRMSteps |
| std ::unique_ptr< SparseLinearSystemNM > | linSolver |
| linear system solver | |
| LinSystSolverType | solverType |
| linear system solver ID | |
| SparseMtrx::SparseMtrxVersionType | smConstraintVersion |
| sparse matrix version, used to control constrains application to stiffness | |
| int | numberOfPrescribedDofs |
| number of prescribed displacement | |
| bool | prescribedDofsFlag |
| IntArray | prescribedDofs |
| Array of pairs identifying prescribed dofs (node, dof). | |
| FloatArray | prescribedDofsValues |
| Array of prescribed values. | |
| int | prescribedDisplacementTF |
| Load Time Function of prescribed values. | |
| IntArray | prescribedEqs |
| Array of prescribed equations. | |
| bool | prescribedEqsInitFlag |
| Flag indicating that prescribedEqs were initialized. | |
| FloatArray | lastReactions |
| Computed reactions. They are stored in order to print them in printState method. | |
| bool | lsFlag |
| Flag indicating whether to use line-search. | |
| std ::unique_ptr< LineSearchNM > | linesearchSolver |
| Line search solver. | |
| bool | mCalcStiffBeforeRes |
| Flag indicating if the stiffness should be evaluated before the residual in the first iteration. | |
| bool | constrainedNRFlag |
| Flag indicating whether to use constrained Newton. | |
| double | constrainedNRalpha |
| Scale factor for dX, dX_new = alpha * dX. | |
| int | constrainedNRminiter |
| Minimum number of iterations before constraint is activated. | |
| FloatArray | rtolf |
| Relative unbalanced force tolerance for each group. | |
| FloatArray | rtold |
| Relative iterative displacement change tolerance for each group. | |
| FloatArray | forceErrVec |
| FloatArray | forceErrVecOld |
| bool | solutionDependentExternalForcesFlag |
| Solution dependent external forces - updating then each NR iteration. | |
| std ::map< int, double > | dg_forceScale |
| Optional user supplied scale of forces used in convergence check. | |
| double | maxIncAllowed |
| Protected Attributes inherited from oofem::SparseNonLinearSystemNM | |
| double | deltaL =0.0 |
| Load level. | |
| double | randPertAmplitude =0.0 |
| Amplitude of a random perturbation applied on the solution before the iteration process. | |
| int | randSeed =0 |
| bool | pert_init_needed |
| IntArray | igp_PertDmanDofSrcArray |
| FloatArray | igp_PertWeightArray |
| IntArray | igp_Map |
| FloatArray | igp_Weight |
| Protected Attributes inherited from oofem::NumericalMethod | |
| Domain * | domain |
| Pointer to domain. | |
| EngngModel * | engngModel |
| Pointer to engineering model. | |
Additional Inherited Members | |
| Public Types inherited from oofem::SparseNonLinearSystemNM | |
| enum | referenceLoadInputModeType { rlm_total =0 , rlm_incremental =1 } |
This class implements Newton-Raphson Method, derived from abstract NumericalMethod class for solving non-linear problems. Except the traditional load control, it also provides direct displacement control without requiring BC applied, so that the equation renumbering is not required when combining arc-length and Newton-Raphson solvers within a simulation.
The direct displacement control is achieved by adding a large number alpha to the corresponding diagonal member of K and replacing the right-hand side by alpha*prescribed_value. If alpha is very much larger than other stiffness coefficients then this alteration effectively replaces the corresponding equation by alpha*unknown_value = alpha*prescribed_value that is, the required condition, but the whole system remains symmetric and minimal changes are necessary in the computational sequence. The above artifice has been introduced by Payne and Irons.
Definition at line 97 of file nrsolver.h.
|
protected |
| Enumerator | |
|---|---|
| nrsolverModifiedNRM | |
| nrsolverFullNRM | |
| nrsolverAccelNRM | |
Definition at line 100 of file nrsolver.h.
| oofem::NRSolver::NRSolver | ( | Domain * | d, |
| EngngModel * | m ) |
References initializeFrom(), NRSolver(), printState(), and solve().
Referenced by NRSolver(), and oofem::StaggeredSolver::StaggeredSolver().
|
virtual |
Definition at line 98 of file nrsolver.C.
|
protected |
Definition at line 512 of file nrsolver.C.
References oofem::FloatArray::at(), oofem::SparseMtrx::at(), oofem::PetscSparseMtrx::createVecGlobal(), oofem::diag(), oofem::NumericalMethod::engngModel, oofem::TimeStep::giveTargetTime(), oofem::TimeStep::giveTimeIncrement(), oofem::TimeStep::isTheFirstStep(), numberOfPrescribedDofs, prescribedDisplacementTF, prescribedDofsValues, prescribedEqs, oofem::SparseNonLinearSystemNM::rlm_total, solverType, and oofem::ST_Petsc.
Referenced by solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 468 of file nrsolver.C.
References oofem::SparseMtrx::at(), oofem::PetscSparseMtrx::createVecGlobal(), oofem::diag(), oofem::PetscSparseMtrx::giveMtrx(), oofem::SparseMtrx::giveVersion(), numberOfPrescribedDofs, prescribedEqs, and smConstraintVersion.
Referenced by solve(), and oofem::StaggeredSolver::solve().
|
protected |
Determines whether or not the solution has reached convergence.
Definition at line 594 of file nrsolver.C.
References oofem::__DofIDItemToString(), oofem::ParallelContext::accumulate(), oofem::FloatArray::at(), oofem::IntArray::at(), constrainedNRFlag, dg_forceScale, oofem::NumericalMethod::domain, oofem::Element_local, oofem::NumericalMethod::engngModel, forceErrVec, forceErrVecOld, oofem::FloatArray::giveSize(), oofem::ParallelContext::isLocal(), oofem::ParallelContext::localNorm(), oofem::macroScale, nrsolver_ERROR_NORM_SMALL_NUM, NRSOLVER_MAX_REL_ERROR_BOUND, OOFEM_LOG_INFO, oofem::FloatArray::resize(), rtold, rtolf, oofem::FloatArray::zero(), and oofem::IntArray::zero().
Referenced by oofem::DynamicRelaxationSolver::solve(), solve(), and oofem::StaggeredSolver::solve().
|
inlineoverridevirtual |
Reimplemented from oofem::SparseNonLinearSystemNM.
Reimplemented in oofem::StaggeredSolver.
Definition at line 179 of file nrsolver.h.
|
inlinevirtual |
Reimplemented in oofem::DynamicRelaxationSolver, and oofem::StaggeredSolver.
Definition at line 180 of file nrsolver.h.
References _IFT_NRSolver_Name.
|
overridevirtual |
Constructs (if necessary) and returns a linear solver. Public method because some problems require it for sensitivity analysis, etc. even for nonlinear problems (e.g. tangent relations in multiscale simulations).
Reimplemented from oofem::SparseNonLinearSystemNM.
Definition at line 406 of file nrsolver.C.
References oofem::classFactory, oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, linSolver, OOFEM_ERROR, and solverType.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Constructs and returns a line search solver.
Definition at line 424 of file nrsolver.C.
References oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, and linesearchSolver.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
overridevirtual |
Reimplemented from oofem::NumericalMethod.
Reimplemented in oofem::StaggeredSolver.
Definition at line 103 of file nrsolver.C.
References _IFT_NRSolver_calcstiffbeforeres, _IFT_NRSolver_constrainedNRminiter, _IFT_NRSolver_ddfunc, _IFT_NRSolver_ddm, _IFT_NRSolver_ddv, _IFT_NRSolver_forceScale, _IFT_NRSolver_forceScaleDofs, _IFT_NRSolver_linesearch, _IFT_NRSolver_lstype, _IFT_NRSolver_manrmsteps, _IFT_NRSolver_maxinc, _IFT_NRSolver_maxiter, _IFT_NRSolver_miniterations, _IFT_NRSolver_minsteplength, _IFT_NRSolver_rtold, _IFT_NRSolver_rtolf, _IFT_NRSolver_rtolv, _IFT_NRSolver_solutionDependentExternalForces, constrainedNRFlag, constrainedNRminiter, dg_forceScale, giveLinearSolver(), giveLineSearchSolver(), oofem::IntArray::giveSize(), oofem::InputRecord::hasField(), IR_GIVE_FIELD, IR_GIVE_OPTIONAL_FIELD, lsFlag, MANRMSteps, maxIncAllowed, mCalcStiffBeforeRes, minIterations, minStepLength, NR_Mode, NR_OldMode, nrsolverAccelNRM, nrsolverModifiedNRM, nsmax, numberOfPrescribedDofs, OOFEM_ERROR, prescribedDisplacementTF, prescribedDofs, prescribedDofsFlag, prescribedDofsValues, rtold, rtolf, solutionDependentExternalForcesFlag, solverType, and oofem::Vec1().
Referenced by NRSolver().
|
protected |
Initiates prescribed equations.
Definition at line 434 of file nrsolver.C.
References oofem::IntArray::at(), oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, oofem::ParallelContext::isLocal(), numberOfPrescribedDofs, prescribedDofs, prescribedEqs, and prescribedEqsInitFlag.
Referenced by solve(), and oofem::StaggeredSolver::solve().
|
overridevirtual |
Prints status message of receiver to output stream. Prints the message corresponding to last solve.
| outputStream | Stream to print state to. |
Reimplemented from oofem::SparseNonLinearSystemNM.
Definition at line 572 of file nrsolver.C.
References lastReactions, numberOfPrescribedDofs, and prescribedDofs.
Referenced by NRSolver().
|
inlineoverridevirtual |
Reinitializes the receiver. This is used, when topology of problem has changed (for example after adaptive refinement or load transfer in parallel applications). This is necessary for numerical methods, that cache some data between solution steps and that may depend on domain or problem topology. The caching of data by receiver is intended only for speeding up the calculation, but numerical method must be always able to generate this data again. This method clears receiver cached data dependent on topology, when it changes.
Reimplemented from oofem::NumericalMethod.
Definition at line 191 of file nrsolver.h.
References linSolver.
|
inlineoverridevirtual |
Reimplemented from oofem::NumericalMethod.
Definition at line 182 of file nrsolver.h.
References oofem::NumericalMethod::domain, linesearchSolver, and linSolver.
|
overridevirtual |
Solves the given sparse linear system of equations \( s R + R_0 - F(X) = 0 \). Total load vector is not passed, it is defined as \( s R + R_0 \), where s is scale factor.
| K | Coefficient matrix ( \(\displaystyle K = \frac{\partial F}{\partial X} \); stiffness matrix). |
| R | Reference incremental RHS (incremental load). |
| R0 | Initial RHS (initial load). |
| X | Total solution (total displacement). |
| dX | Increment of solution (incremental displacements). |
| F | InternalRhs (real internal forces). |
| internalForcesEBENorm | Norm of internal nodal forces (evaluated on element by element basis) (split into each DOF id). |
| s | RHS scale factor (load level). |
| rlm | Reference load mode. |
| nite | Number of iterations needed. |
| tStep | Time step to solve for. |
Implements oofem::SparseNonLinearSystemNM.
Reimplemented in oofem::StaggeredSolver.
Definition at line 205 of file nrsolver.C.
References oofem::FloatArray::add(), applyConstraintsToLoadIncrement(), applyConstraintsToStiffness(), oofem::FloatArray::at(), oofem::FloatArray::beDifferenceOf(), checkConvergence(), constrainedNRalpha, constrainedNRFlag, constrainedNRminiter, oofem::CR_CONVERGED, oofem::CR_DIVERGED_ITS, oofem::CR_DIVERGED_TOL, oofem::CR_UNKNOWN, oofem::SparseNonLinearSystemNM::deltaL, oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, oofem::ExternalRhs, forceErrVec, giveLinearSolver(), giveLineSearchSolver(), oofem::FloatArray::giveSize(), oofem::TimeStep::incrementStateCounter(), oofem::TimeStep::incrementSubStepNumber(), initPrescribedEqs(), oofem::InternalRhs, lastReactions, linSolver, oofem::ParallelContext::localNorm(), lsFlag, oofem::macroScale, MANRMSteps, maxIncAllowed, mCalcStiffBeforeRes, minIterations, oofem::NonLinearLhs, NR_Mode, nrsolverAccelNRM, nrsolverFullNRM, nsmax, numberOfPrescribedDofs, OOFEM_LOG_INFO, prescribedDofs, prescribedDofsFlag, prescribedEqs, prescribedEqsInitFlag, oofem::FloatArray::resize(), rtold, rtolf, solutionDependentExternalForcesFlag, oofem::FloatArray::times(), and oofem::FloatArray::zero().
Referenced by NRSolver().
|
protected |
Scale factor for dX, dX_new = alpha * dX.
Definition at line 145 of file nrsolver.h.
Referenced by solve(), and oofem::StaggeredSolver::solve().
|
protected |
Flag indicating whether to use constrained Newton.
Definition at line 143 of file nrsolver.h.
Referenced by checkConvergence(), oofem::StaggeredSolver::checkConvergenceDofIdArray(), initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Minimum number of iterations before constraint is activated.
Definition at line 147 of file nrsolver.h.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Optional user supplied scale of forces used in convergence check.
Definition at line 164 of file nrsolver.h.
Referenced by checkConvergence(), and initializeFrom().
|
protected |
Definition at line 155 of file nrsolver.h.
Referenced by checkConvergence(), oofem::StaggeredSolver::checkConvergenceDofIdArray(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 156 of file nrsolver.h.
Referenced by checkConvergence(), and oofem::StaggeredSolver::checkConvergenceDofIdArray().
|
protected |
Computed reactions. They are stored in order to print them in printState method.
Definition at line 135 of file nrsolver.h.
Referenced by printState(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Line search solver.
Definition at line 139 of file nrsolver.h.
Referenced by giveLineSearchSolver(), and setDomain().
|
protected |
linear system solver
Definition at line 109 of file nrsolver.h.
Referenced by giveLinearSolver(), reinitialize(), setDomain(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Flag indicating whether to use line-search.
Definition at line 137 of file nrsolver.h.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 106 of file nrsolver.h.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 166 of file nrsolver.h.
Referenced by initializeFrom(), and solve().
|
protected |
Flag indicating if the stiffness should be evaluated before the residual in the first iteration.
Definition at line 141 of file nrsolver.h.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 102 of file nrsolver.h.
Referenced by initializeFrom(), oofem::DynamicRelaxationSolver::solve(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 103 of file nrsolver.h.
Referenced by initializeFrom().
|
protected |
Definition at line 104 of file nrsolver.h.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Definition at line 105 of file nrsolver.h.
|
protected |
Definition at line 104 of file nrsolver.h.
Referenced by initializeFrom().
|
protected |
Definition at line 102 of file nrsolver.h.
Referenced by initializeFrom(), oofem::DynamicRelaxationSolver::solve(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
number of prescribed displacement
Definition at line 115 of file nrsolver.h.
Referenced by applyConstraintsToLoadIncrement(), applyConstraintsToStiffness(), initializeFrom(), initPrescribedEqs(), printState(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Load Time Function of prescribed values.
Definition at line 129 of file nrsolver.h.
Referenced by applyConstraintsToLoadIncrement(), and initializeFrom().
|
protected |
Array of pairs identifying prescribed dofs (node, dof).
Definition at line 125 of file nrsolver.h.
Referenced by initializeFrom(), initPrescribedEqs(), printState(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Flag indicating that some dofs are controlled under displacement control. In parallel mode, numberOfPrescribedDofs is local (related to specific partition) so its nonzero value does not mean that there are no prescribed dofs on other partitions.
Definition at line 122 of file nrsolver.h.
Referenced by initializeFrom(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Array of prescribed values.
Definition at line 127 of file nrsolver.h.
Referenced by applyConstraintsToLoadIncrement(), and initializeFrom().
|
protected |
Array of prescribed equations.
Definition at line 131 of file nrsolver.h.
Referenced by applyConstraintsToLoadIncrement(), applyConstraintsToStiffness(), initPrescribedEqs(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Flag indicating that prescribedEqs were initialized.
Definition at line 133 of file nrsolver.h.
Referenced by initPrescribedEqs(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Relative iterative displacement change tolerance for each group.
Definition at line 152 of file nrsolver.h.
Referenced by checkConvergence(), oofem::StaggeredSolver::checkConvergenceDofIdArray(), initializeFrom(), oofem::DynamicRelaxationSolver::solve(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
Relative unbalanced force tolerance for each group.
Definition at line 150 of file nrsolver.h.
Referenced by checkConvergence(), oofem::StaggeredSolver::checkConvergenceDofIdArray(), initializeFrom(), oofem::DynamicRelaxationSolver::solve(), solve(), and oofem::StaggeredSolver::solve().
|
protected |
sparse matrix version, used to control constrains application to stiffness
Definition at line 113 of file nrsolver.h.
Referenced by applyConstraintsToStiffness().
|
protected |
Solution dependent external forces - updating then each NR iteration.
Definition at line 159 of file nrsolver.h.
Referenced by initializeFrom(), and solve().
|
protected |
linear system solver ID
Definition at line 111 of file nrsolver.h.
Referenced by applyConstraintsToLoadIncrement(), giveLinearSolver(), and initializeFrom().