xref: /petsc/src/snes/impls/al/alimpl.h (revision 447e0e2062395afab88d46916d90376f9cc66fab)
1*97276fddSZach Atkins /*
2*97276fddSZach Atkins    Private context for a Newton arc length continuation
3*97276fddSZach Atkins    method for solving systems of nonlinear equations
4*97276fddSZach Atkins */
5*97276fddSZach Atkins 
6*97276fddSZach Atkins #pragma once
7*97276fddSZach Atkins #include <petsc/private/snesimpl.h>
8*97276fddSZach Atkins #include <petscsnes.h>
9*97276fddSZach Atkins 
10*97276fddSZach Atkins typedef struct {
11*97276fddSZach Atkins   PetscInt                   max_continuation_steps; /* maximum number of continuation steps */
12*97276fddSZach Atkins   PetscReal                  step_size;              /* radius of quadratic constraint surface */
13*97276fddSZach Atkins   PetscReal                  psisq;                  /* load factor regularization parameter */
14*97276fddSZach Atkins   PetscReal                  lambda_update;          /* accumulated update to lambda over the current increment */
15*97276fddSZach Atkins   PetscReal                  lambda;                 /* load parameter */
16*97276fddSZach Atkins   PetscReal                  lambda_max;             /* maximum value of the load parameter */
17*97276fddSZach Atkins   PetscReal                  lambda_min;             /* minimum value of the load parameter */
18*97276fddSZach Atkins   PetscBool                  scale_rhs;              /* should the RHS vector be scaled by the load parameter? */
19*97276fddSZach Atkins   SNESNewtonALCorrectionType correction_type;        /* type of correction scheme to use */
20*97276fddSZach Atkins   PetscBool                  copied_rhs;             /* has the right-hand side vector been copied? */
21*97276fddSZach Atkins 
22*97276fddSZach Atkins   Vec vec_rhs_orig; /* original right-hand side vector, used if `scale_rhs == PETSC_TRUE` */
23*97276fddSZach Atkins 
24*97276fddSZach Atkins   SNESFunctionFn *computealfunction; /* user-provided function to compute the tangent load vector  */
25*97276fddSZach Atkins   void           *alctx;             /* user-provided context for the tangent load vector computation */
26*97276fddSZach Atkins } SNES_NEWTONAL;
27*97276fddSZach Atkins 
28*97276fddSZach Atkins PETSC_INTERN const char NewtonALExactCitation[];
29*97276fddSZach Atkins PETSC_INTERN PetscBool  NewtonALExactCitationSet;
30*97276fddSZach Atkins PETSC_INTERN const char NewtonALNormalCitation[];
31*97276fddSZach Atkins PETSC_INTERN PetscBool  NewtonALNormalCitationSet;
32