18.1. Reaction Diffusion Solver

The phase-field model solves the phase-field evolution using a reaction-diffusion solver implemented in Sierra/SM. The reaction-diffusion command block is as follows:

BEGIN REACTION DIFFUSION rxndiffname
  BLOCK = <string list>block_names
  REMOVE BLOCK = <string list>remove_block_names
  INCLUDE ALL BLOCKS

  ASSEMBLY = <string list>assembly_names
  REMOVE ASSEMBLY = <string list>remove_assembly_names

  ACTIVE PERIODS = <string list> active_period_names
  INACTIVE PERIODS = <string list> inactive_period_names

  PHASE FIELD BOUND CONSTRAINTS

  BEGIN PRESCRIBED FIELD
    SURFACE = <string>surf_name
    NODE SET|NODESET = <string>node_set_name
    FUNCTION = <string>func_name
  END
  BEGIN PRESCRIBED FLUX
    SURFACE = <string>surf_name
    NODE SET|NODESET = <string>node_set_name
    FUNCTION = <string>func_name
  END

  GRADIENT CONFIGURATION = MODEL|CURRENT(MODEL)
  USE FINITE ELEMENT MODEL = <string>model_name
      [Model Coordinates Are <string> nodal_variable_name]

  SOLVE AT INITIALIZATION = OFF|ON(OFF)
  SOLVE STEP INCREMENT = <integer>solv_step_incr(1)
  SOLVE TIME INCREMENT = <real>solv_time_incr

  EQUATION SYSTEM = LINEAR|NONLINEAR(LINEAR)
  USE LINEAR SOLVER = FALSE|TRUE(TRUE)

  BEGIN SOLVER OPTIONS
    ACCEPTABLE RESIDUAL = <real>accept_resid(1.0e-14)
    ACCEPTABLE RELATIVE RESIDUAL = <real>accept_rel_resid(1.0e-12)
    TARGET RESIDUAL = <real>target_resid(1.0e-15)
    TARGET RELATIVE RESIDUAL = <real>target_rel_resid(1.0e-13)
    MAXIMUM ITERATIONS = <integer>max_iter(1000)
    LINE SEARCH BACKTRACK [<real>min_step(1.0e-8) <real>reduction(0.1)
        <real>init_step(1.0)
        BISECTION|CUTBACK|SECANT|TRUST_REGION|CRITICAL_POINT(BISECTION)
        <real>dispinc(0.5) <real>rotinc(0.0) ]
    LINEAR SOLVER = <string>linear_solver_name(DEFAULT_FETI_SOLVER)
    PRECONDITIONER = PROBE|IDENTITY(PROBE)
  END

END

The commands regarding BLOCK, ASSEMBLY, and ACTIVE PERIODS are common to Sierra/SM, and allow the user to specify reaction-diffusion solution only on the blocks or time periods of interest.

The SOLVE AT INITIALIZATION command provides a way to solve the phase field before any mechanical loading is applied. This is particularly useful when applying a phase boundary condition (e.g. prescribed field, prescribed flux), as the phase field will not be in equilibrium prior to the mechanics solve in the first time-step.

Explicit only

The commands SOLVE STEP INCREMENT and SOLVE TIME INCREMENT offer a way of decreasing the frequency of reaction diffusion solves, particularly for the use-case involving implicit reaction diffusion solves with explicit time incrementation of the solid mechanics solution. Individually, the two commands define the period (in steps, or in time) between successive reaction-diffusion solves. It is an error to specify both commands simultaneously.

The command PHASE FIELD BOUND CONSTRAINTS offers a way of ensuring satisfaction of the phase field bound constraints. Specifically, this command restricts the reaction-diffusion solution to the interval \(\phi \in [0,1]\) and also ensures that the field evolves monotonically \(\Delta \phi \leq 0\). Implementation is done using a bound-constrained conjugate gradient method [[1]] with modifications around the active set selection.

The commands BEGIN PRESCRIBED FIELD and BEGIN PRESCRIBED FLUX offer a means to specify Dirichlet or Neumann boundary conditions on the phase field. The prescribed field is commonly used to instantiate the model with an initial damage field, such as a phase field pre-crack. Often this provides for better numerical performance than a meshed sharp crack with no initial phase field. The prescribed flux has not been tested, as it lacks a clear physical motivation. It is recommended to use the SOLVE AT INITIALIZATION command with these options to ensure that the phase field is in equilibrium before performing any mechanical solution.

The command GRADIENT CONFIGURATION provides the option to compute the phase-field gradient in model or current coordinates. Model coordinates are recommended.

The command USE FINITE ELEMENT MODEL provides the user the option to specify the finite element model for the phase field. By default, the same finite element model as the Adagio region is used. Using a different finite element mesh for the reaction-diffusion solve has not been tested and is not recommended.

The commands EQUATION SYSTEM and USE LINEAR SYSTEM are used to determine whether to solve the reaction-diffusion equation as a linear or nonlinear partial differential equation. For the current Phase Field FeFp implementation, the phase-field evolution equation is a linear partial differential equation when the degradation function parameter \(\gamma = 0\) (corresponding to \(\psi_c = \frac{3 G_c}{16 \ell}\)), so the linear reaction-diffusion solve may be used: LINEAR, TRUE; this uses the default FETI solver to solve the linear system. For other cases (\(\gamma \ne 0\)), it is the phase-field evolution equation is a nonlinear partial differential equation, so it is recommended to use the nonlinear reaction-diffusion solver: NONLINEAR, FALSE; this uses a preconditioned conjugate gradient method to iteratively solve the system.

The BEGIN SOLVER OPTIONS block offers options to customize the preconditioned conjugate gradient solve (EQUATION SYSTEM = NONLINEAR, USE LINEAR SYSTEM = FALSE), with commands that reflect a limited set of the commands in the Sierra/SM BEGIN CG block. Two nodal preconditioners are available: identity and probe. The probe preconditioner defines its entries as the quotient of the phase system residual divided by the nodal phase stiffness; the phase stiffness is computed as the derivative of the phase force-internal at that node with respect to phase, estimated using a forward finite difference (\(\Delta \phi = -1.0\cdot10^{-10}\)).