|
OOFEM 3.0
|
#include <calmls.h>
Public Member Functions | |
| CylindricalALM (Domain *d, EngngModel *m) | |
| virtual | ~CylindricalALM () |
| ConvergedReason | solve (SparseMtrx &K, FloatArray &R, FloatArray *R0, FloatArray &X, FloatArray &dX, FloatArray &F, const FloatArray &internalForcesEBENorm, double &ReachedLambda, referenceLoadInputModeType rlm, int &nite, TimeStep *) override |
| double | giveCurrentStepLength () override |
| void | setStepLength (double s) override |
| void | initializeFrom (InputRecord &ir) override |
| bool | referenceLoad () const override |
| void | saveContext (DataStream &stream, ContextMode mode) override |
| void | restoreContext (DataStream &stream, ContextMode mode) override |
| void | setDomain (Domain *d) override |
| void | reinitialize () override |
| const char * | giveClassName () const override |
| const char * | giveInputRecordName () const |
| SparseLinearSystemNM * | giveLinearSolver () override |
| Public Member Functions inherited from oofem::SparseNonLinearSystemNM | |
| SparseNonLinearSystemNM (Domain *d, EngngModel *m) | |
| Constructor. | |
| virtual | ~SparseNonLinearSystemNM () |
| Destructor. | |
| virtual void | printState (FILE *outputStream) |
| 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 () |
Protected Types | |
| enum | calm_ControlType { calm_hpc_off = 0 , calm_hpc_on , calml_hpc } |
| CALM mode type; determines the calm step length control. More... | |
| enum | calm_NR_ModeType { calm_modifiedNRM , calm_fullNRM , calm_accelNRM } |
| Controlling mode of Newton-Raphson Method. More... | |
| enum | RootSelectionType { RST_Cos =0 , RST_Dot =1 } |
| typedef std ::set< DofIDItem > | __DofIDSet |
Protected Member Functions | |
| void | convertHPCMap () |
| int | computeDeltaLambda (double &deltaLambda, const FloatArray &dX, const FloatArray &deltaXt, const FloatArray &deltaX_, const FloatArray &R, double RR, double eta, double deltaL, double DeltaLambda0, int neq) |
| void | search (int istep, FloatArray &prod, FloatArray &eta, double amp, double maxeta, double mineta, int &status) |
| bool | checkConvergence (const FloatArray &R, const FloatArray *R0, const FloatArray &F, const FloatArray &X, const FloatArray &ddX, double Lambda, double RR0, double RR, double drProduct, const FloatArray &internalForcesEBENorm, int nite, bool &errorOutOfRange) |
| Evaluates the convergence criteria. | |
| void | do_lineSearch (FloatArray &X, const FloatArray &XInitial, const FloatArray &deltaX_, const FloatArray &deltaXt, const FloatArray &dXm1, FloatArray &dX, FloatArray &ddX, const FloatArray &R, const FloatArray *R0, const FloatArray &F, double &DeltaLambda, double &DeltaLambdam1, double &deltaLambda, double &Lambda, double &ReachedLambda, double RR, double &drProduct, TimeStep *tStep) |
| Perform line search optimization of step length. | |
Protected Attributes | |
| int | nsmax |
| int | maxRestarts |
| double | Psi |
| double | deltaL |
| double | minStepLength |
| double | maxStepLength |
| int | solved |
| int | numberOfRequiredIterations |
| calm_NR_ModeType | calm_NR_Mode |
| calm_NR_ModeType | calm_NR_OldMode |
| int | calm_NR_ModeTick |
| int | calm_MANRMSteps |
| int | minIterations |
| Minimum hard number of iteration.s. | |
| int | calm_hpc_init |
| Variables for HyperPlaneControl. | |
| calm_ControlType | calm_Control |
| FloatArray | calm_HPCWeights |
| IntArray | calm_HPCIndirectDofMask |
| Array containing equation numbers of dofs under indirect control. | |
| IntArray | calm_HPCDmanDofSrcArray |
| Input array containing dofmanagers and corresponding dof numbers under indirect control. | |
| FloatArray | calm_HPCDmanWeightSrcArray |
| Input array of dofman weights (for hpcmode 2). | |
| std ::unique_ptr< SparseLinearSystemNM > | linSolver |
| Linear system solver. | |
| LinSystSolverType | solverType |
| linear system solver ID. | |
| int | lsFlag |
| Line search flag. | |
| double | ls_tolerance |
| Line search tolerance. | |
| double | amplifFactor |
| Line search amplification factor. | |
| double | maxEta |
| Line search parameters (limits). | |
| double | minEta |
| int | nccdg |
| Number of convergence criteria dof groups. | |
| std ::vector< __DofIDSet > | ccDofGroups |
| Convergence criteria dof groups. | |
| FloatArray | rtolf |
| Relative unbalanced force tolerance for each group. | |
| FloatArray | rtold |
| Relative iterative displacement change tolerance for each group. | |
| ParallelContext * | parallel_context |
| Parallel context for computing norms, dot products and such. | |
| RootSelectionType | rootselectiontype |
| Root selection type. | |
| FloatArray | old_dX |
| previous increment of dX, needed by root selection type 1 | |
| 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 } |
Implementation of sparse nonlinear solver with indirect control. It uses Cylindrical Arc-Length Method algorithm.
DESCRIPTION : Perform solution of non-linear problem in the form Kt * deltaR = g where g is defined as g = g(Lambda,X) = Lambda* R - F(X). During iteration process, we want g became zero vector.
This method uses Modified Newton Raphson iteration scheme
If we solve non-linear static we can interpret symbols as follows:
Kt - tangential stiffness dX - increment of displacements g - vector of unbalanced forces (at the end should be zero one) Lambda - Load level Rt - Quasi Total Load vector Rt = R0 + Lambda*R X - total displacement vector F(X) - Nodal representation of (real) internal forces. Psi - control parameter (=0. means displacement control) Bergan_k0 - variable used to compute bergan's parameters of current stiffness. calm_NR_Mode - variable controlling the mode of NRM (ModifiedNR, Full NRM (stiffness update after each iteration), Modified Accelerated NRM (we perform iteration with stiffness matrix updated only after calm_MANRMSteps) calm_NR_OldMode - variable containing the old mode of NRM, which will be restored after calm_NR_ModeTick iterations. calm_NR_ModeTick - see calm_NR_OldMode. calm_MANRMSteps - if calm_NR_Mode == calm_accelNRM, it specifies, that new updated stiffness matrix is assembled after calm_MANRMSteps. calm_Control - variable indicating the ALM control. calm_HPCIndirectDofMask - Mask, telling which dofs are used for HPC. calm_HPCWeights - dofs weights in constrain.
Tasks:
Variable description :
OUTPUT : (after call solveYourselfAt)
|
protected |
|
protected |
CALM mode type; determines the calm step length control.
|
protected |
|
protected |
| oofem::CylindricalALM::CylindricalALM | ( | Domain * | d, |
| EngngModel * | m ) |
References CylindricalALM(), and solve().
Referenced by CylindricalALM().
|
protected |
Evaluates the convergence criteria.
Definition at line 484 of file calmls.C.
References oofem::FloatArray::add(), oofem::FloatArray::at(), CALM_MAX_REL_ERROR_BOUND, calm_SMALL_ERROR_NUM, ccDofGroups, oofem::NumericalMethod::domain, oofem::Element_local, nccdg, OOFEM_LOG_INFO, parallel_context, rtold, rtolf, oofem::FloatArray::subtract(), oofem::FloatArray::sum(), oofem::FloatArray::times(), and oofem::FloatArray::zero().
Referenced by solve().
|
protected |
Definition at line 995 of file calmls.C.
References oofem::FloatArray::at(), calm_Control, calm_hpc_off, calm_hpc_on, calm_HPCIndirectDofMask, calm_HPCWeights, calm_SMALL_NUM, calml_hpc, deltaL, OOFEM_ERROR, OOFEM_WARNING, parallel_context, Psi, rootselectiontype, RST_Cos, oofem::Vec2(), and oofem::Vec4().
Referenced by do_lineSearch(), and solve().
|
protected |
Definition at line 900 of file calmls.C.
References oofem::FloatArray::at(), oofem::IntArray::at(), calm_Control, calm_HPCDmanDofSrcArray, calm_HPCDmanWeightSrcArray, calm_HPCIndirectDofMask, calm_HPCWeights, calml_hpc, oofem::NumericalMethod::domain, OOFEM_WARNING, parallel_context, oofem::FloatArray::resize(), and oofem::IntArray::resize().
Referenced by solve().
|
protected |
Perform line search optimization of step length.
Definition at line 1235 of file calmls.C.
References oofem::FloatArray::add(), amplifFactor, oofem::FloatArray::at(), oofem::FloatArray::clear(), computeDeltaLambda(), deltaL, oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, oofem::FloatArray::giveSize(), oofem::TimeStep::incrementStateCounter(), oofem::InternalRhs, ls_tolerance, maxEta, minEta, OOFEM_LOG_INFO, parallel_context, and search().
Referenced by solve().
|
inlineoverridevirtual |
Reimplemented from oofem::SparseNonLinearSystemNM.
|
inlineoverridevirtual |
|
inline |
Definition at line 257 of file calmls.h.
References _IFT_CylindricalALM_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 977 of file calmls.C.
References oofem::classFactory, oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, linSolver, OOFEM_ERROR, and solverType.
Referenced by initializeFrom(), and solve().
|
overridevirtual |
Reimplemented from oofem::NumericalMethod.
Definition at line 681 of file calmls.C.
References _IFT_CylindricalALM_ccdg, _IFT_CylindricalALM_forcedinitialsteplength, _IFT_CylindricalALM_hpc, _IFT_CylindricalALM_hpcmode, _IFT_CylindricalALM_hpcw, _IFT_CylindricalALM_initialsteplength, _IFT_CylindricalALM_linesearch, _IFT_CylindricalALM_lsearchamp, _IFT_CylindricalALM_lsearchmaxeta, _IFT_CylindricalALM_lsearchtol, _IFT_CylindricalALM_lstype, _IFT_CylindricalALM_manrmsteps, _IFT_CylindricalALM_maxiter, _IFT_CylindricalALM_maxrestarts, _IFT_CylindricalALM_miniterations, _IFT_CylindricalALM_minsteplength, _IFT_CylindricalALM_nccdg, _IFT_CylindricalALM_psi, _IFT_CylindricalALM_reqiterations, _IFT_CylindricalALM_rootselectiontype, _IFT_CylindricalALM_rtold, _IFT_CylindricalALM_rtolf, _IFT_CylindricalALM_rtolv, _IFT_CylindricalALM_steplength, amplifFactor, oofem::IntArray::at(), calm_accelNRM, calm_Control, calm_hpc_init, calm_hpc_off, calm_hpc_on, calm_HPCDmanDofSrcArray, calm_HPCDmanWeightSrcArray, calm_HPCWeights, calm_MANRMSteps, calm_modifiedNRM, calm_NR_Mode, calm_NR_OldMode, calml_hpc, ccDofGroups, deltaL, giveLinearSolver(), oofem::IntArray::giveSize(), IR_GIVE_FIELD, IR_GIVE_OPTIONAL_FIELD, ls_tolerance, lsFlag, maxEta, maxRestarts, maxStepLength, minIterations, minStepLength, nccdg, nsmax, numberOfRequiredIterations, OOFEM_ERROR, Psi, rootselectiontype, RST_Cos, RST_Dot, rtold, rtolf, and solverType.
|
inlineoverridevirtual |
Returns true if reference loads are used (i.e. arc length methods).
Reimplemented from oofem::SparseNonLinearSystemNM.
|
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 250 of file calmls.h.
References calm_hpc_init, and linSolver.
|
overridevirtual |
Reimplemented from oofem::NumericalMethod.
Definition at line 968 of file calmls.C.
References oofem::CIO_IOERR, deltaL, oofem::DataStream::read(), and THROW_CIOERR.
|
overridevirtual |
Reimplemented from oofem::NumericalMethod.
Definition at line 959 of file calmls.C.
References oofem::CIO_IOERR, deltaL, THROW_CIOERR, and oofem::DataStream::write().
|
protected |
Definition at line 1149 of file calmls.C.
References oofem::FloatArray::at(), and oofem::max().
Referenced by do_lineSearch().
|
inlineoverridevirtual |
Reimplemented from oofem::NumericalMethod.
Definition at line 244 of file calmls.h.
References oofem::NumericalMethod::domain, and linSolver.
|
inlineoverridevirtual |
|
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.
Definition at line 108 of file calmls.C.
References oofem::FloatArray::add(), oofem::FloatArray::at(), calm_accelNRM, calm_Control, CALM_DEFAULT_NRM_TICKS, calm_fullNRM, calm_hpc_init, calm_hpc_off, calm_hpc_on, calm_HPCIndirectDofMask, calm_HPCWeights, calm_MANRMSteps, calm_NR_Mode, calm_NR_ModeTick, calm_NR_OldMode, CALM_RESET_STEP_REDUCE, calml_hpc, checkConvergence(), oofem::FloatArray::clear(), computeDeltaLambda(), convertHPCMap(), oofem::CR_CONVERGED, oofem::CR_DIVERGED_ITS, oofem::CR_UNKNOWN, deltaL, do_lineSearch(), oofem::NumericalMethod::domain, oofem::NumericalMethod::engngModel, giveLinearSolver(), oofem::TimeStep::giveNumber(), oofem::FloatArray::giveSize(), oofem::TimeStep::incrementStateCounter(), oofem::TimeStep::incrementSubStepNumber(), oofem::InternalRhs, linSolver, lsFlag, maxRestarts, maxStepLength, minIterations, minStepLength, nccdg, oofem::NonLinearLhs, nsmax, numberOfRequiredIterations, old_dX, OOFEM_ERROR, OOFEM_LOG_INFO, OOFEM_WARNING, parallel_context, Psi, oofem::FloatArray::resize(), rootselectiontype, RST_Cos, oofem::sgn(), solved, oofem::FloatArray::subtract(), oofem::FloatArray::times(), oofem::Vec2(), and oofem::FloatArray::zero().
Referenced by CylindricalALM().
|
protected |
Line search amplification factor.
Definition at line 198 of file calmls.h.
Referenced by do_lineSearch(), and initializeFrom().
|
protected |
Definition at line 179 of file calmls.h.
Referenced by computeDeltaLambda(), convertHPCMap(), initializeFrom(), and solve().
|
protected |
Variables for HyperPlaneControl.
Definition at line 178 of file calmls.h.
Referenced by initializeFrom(), reinitialize(), and solve().
|
protected |
Input array containing dofmanagers and corresponding dof numbers under indirect control.
Definition at line 184 of file calmls.h.
Referenced by convertHPCMap(), and initializeFrom().
|
protected |
Input array of dofman weights (for hpcmode 2).
Definition at line 186 of file calmls.h.
Referenced by convertHPCMap(), and initializeFrom().
|
protected |
Array containing equation numbers of dofs under indirect control.
Definition at line 182 of file calmls.h.
Referenced by computeDeltaLambda(), convertHPCMap(), and solve().
|
protected |
Definition at line 180 of file calmls.h.
Referenced by computeDeltaLambda(), convertHPCMap(), initializeFrom(), and solve().
|
protected |
Definition at line 172 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Definition at line 170 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
|
protected |
Definition at line 170 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Convergence criteria dof groups.
Definition at line 206 of file calmls.h.
Referenced by checkConvergence(), and initializeFrom().
|
protected |
Definition at line 168 of file calmls.h.
Referenced by computeDeltaLambda(), do_lineSearch(), giveCurrentStepLength(), initializeFrom(), restoreContext(), saveContext(), setStepLength(), and solve().
|
protected |
Linear system solver.
Definition at line 189 of file calmls.h.
Referenced by giveLinearSolver(), reinitialize(), setDomain(), and solve().
|
protected |
Line search tolerance.
Definition at line 196 of file calmls.h.
Referenced by do_lineSearch(), and initializeFrom().
|
protected |
Line search flag.
Definition at line 194 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Line search parameters (limits).
Definition at line 200 of file calmls.h.
Referenced by do_lineSearch(), and initializeFrom().
|
protected |
Definition at line 166 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Definition at line 168 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Definition at line 200 of file calmls.h.
Referenced by do_lineSearch().
|
protected |
Minimum hard number of iteration.s.
Definition at line 175 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Definition at line 168 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Number of convergence criteria dof groups.
Definition at line 204 of file calmls.h.
Referenced by checkConvergence(), initializeFrom(), and solve().
|
protected |
Definition at line 165 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
Definition at line 169 of file calmls.h.
Referenced by initializeFrom(), and solve().
|
protected |
|
protected |
Parallel context for computing norms, dot products and such.
Definition at line 212 of file calmls.h.
Referenced by checkConvergence(), computeDeltaLambda(), convertHPCMap(), do_lineSearch(), and solve().
|
protected |
Definition at line 167 of file calmls.h.
Referenced by computeDeltaLambda(), initializeFrom(), and solve().
|
protected |
Root selection type.
Definition at line 225 of file calmls.h.
Referenced by computeDeltaLambda(), initializeFrom(), and solve().
|
protected |
Relative iterative displacement change tolerance for each group.
Definition at line 210 of file calmls.h.
Referenced by checkConvergence(), and initializeFrom().
|
protected |
Relative unbalanced force tolerance for each group.
Definition at line 208 of file calmls.h.
Referenced by checkConvergence(), and initializeFrom().
|
protected |
|
protected |
linear system solver ID.
Definition at line 191 of file calmls.h.
Referenced by giveLinearSolver(), and initializeFrom().