Mass Diffusion Analysis

A mass diffusion analysis:

  • models the transient or steady-state diffusion of one material through another, such as the diffusion of hydrogen through a metal;

  • requires the use of mass diffusion elements; and

  • can be used to model temperature and/or pressure-driven mass diffusion.

This page discusses:

Governing Equations

The governing equations for mass diffusion are an extension of Fick's equations: they allow for nonuniform solubility of the diffusing substance in the base material and for mass diffusion driven by gradients of temperature and pressure. The basic solution variable (used as the degree of freedom at the nodes of the mesh) is the “normalized concentration” (often also referred to as the “activity” of the diffusing material), ϕ=defc/s, where c is the mass concentration of the diffusing material and s is its solubility in the base material. Therefore, when the mesh includes dissimilar materials that share nodes, the normalized concentration is continuous across the interface between the different materials.

For example, a diatomic gas that dissociates during diffusion can be described using Sievert's law: c = s ( p / p 0 ) 1 2 , where p is the partial pressure of the diffusing gas, and the value of p 0 is customarily one atmosphere (for historical consistency). Combining Sievert's law with the definition of normalized concentration given earlier, ϕ = c / s = ( p / p 0 ) 1 2 . Equilibrium requires the partial pressure to be continuous across an interface, so normalized concentration will be continuous as well. If an expression other than Sievert's law defines the relationship between concentration and partial pressure for a diffusing material, solubility should be defined accordingly.

The diffusion problem is defined from the requirement of mass conservation for the diffusing phase:

VdcdtdV+SnJdS=0,

where V is any volume whose surface is S, n is the outward normal to S, J is the flux of concentration of the diffusing phase, and nJ is the concentration flux leaving S.

Diffusion is assumed to be driven by the gradient of a general chemical potential, which gives the behavior

J=-sD[ϕx+κsx(ln(θ-θZ))+κppx],

where D(c,θ,f) is the diffusivity; s(θ,f) is the solubility; κs(c,θ,f) is the “Soret effect” factor, providing diffusion because of temperature gradient; θ is the temperature; θZ is the value of absolute zero on the temperature scale being used; κp(c,θ,f) is the pressure stress factor, providing diffusion driven by the gradient of the equivalent pressure stress, p=def-trace(σ)/3; σ is stress; and f are any predefined field variables.

Whenever D, κ s , or κ p depends on concentration, the problem becomes nonlinear and the system of equations becomes nonsymmetric. In practical cases the dependence on concentration is quite strong, so the nonsymmetric matrix storage and solution scheme is invoked automatically when a mass diffusion analysis is performed (see Defining an Analysis).

Fick's Law

Mass diffusion behavior is often described by Fick's law (Crank, 1956):

J=-Dcx.

Fick's law is offered in Abaqus/Standard as a special case of the general chemical potential relation. To establish the relationship between Fick's law and the general chemical potential, we write Fick's law as

J=-D(sϕx+ϕsx).

In most practical cases s=s(θ), and we can write

J=-sDϕx-Dcssθθx.

The two terms in this equation describe the normalized concentration and temperature-driven diffusion, respectively. The normalized concentration-driven diffusion term is identical to that given in the general relation. The temperature-driven diffusion term in Fick's law is recovered in the general relation if

κs=c(θ-θZ)s2sθ    .

This conversion is done automatically in Abaqus/Standard when you request Fick's law (see Diffusivity).

An extended form of Fick's law can also be chosen by specifying a nonzero value for κp:

J=-D(cx+sκppx).

In this case Abaqus/Standard will still define κs automatically as discussed earlier.

Units

The units of concentration are commonly given as parts per million (P). Similarly, the units of solubility are P. The Soret effect factor is dimensionless. The units of the pressure stress factor are F−1L2, where F is force and L is length; and the units of equivalent pressure stress are FL−2. The diffusivity, D , has units of L2T−1, where T is time. The concentration flux, J , then has units of PLT−1; and the concentration volumetric flux, Q = S n J d S , has units of PL3T−1.

Steady-State Analysis

Steady-state mass diffusion analysis provides the steady-state solution directly: the rate of change of concentration with respect to time is omitted from the governing diffusion equation in steady-state analysis. In nonlinear cases iteration may be necessary to achieve a converged solution.

Since the rate term is removed from the governing equations, the steady-state problem has no intrinsic physically meaningful time scale; nevertheless, you may assign a “time” scale to the analysis step. This time scale is often convenient for output identification and for specifying prescribed normalized concentrations and fluxes with varying magnitudes. Thus, when steady-state analysis is chosen, you specify a “time” increment and a “time” period for the step; Abaqus/Standard then increments through the step accordingly. If a steady-state analysis step is to be followed by a transient analysis step and total time is used in amplitude definitions (Amplitude Curves), the time period should be defined to be negligibly small in the steady-state step. For more details on time scales and time stepping, see Defining an Analysis.

Transient Analysis

Time integration in transient diffusion analysis is done with the backward Euler method (also referred to as the modified Crank-Nicholson operator). This method is unconditionally stable for linear problems.

Automatic or fixed time incrementation can be used for transient analysis. The automatic time incrementation scheme is generally preferred because the response is usually simple diffusion: the rate of change of normalized concentration varies widely during the step and requires different time increments to maintain accuracy in the time integration.

Spurious Oscillations due to Small Time Increments

In transient mass diffusion analysis with second-order elements there is a relationship between the minimum usable time step and the element size. A simple guideline is

Δt>16DΔ,2

where Δt is the time increment, D is the diffusivity, and Δ is a typical element dimension (such as the length of a side of an element). If time increments smaller than this value are used, spurious oscillations can appear in the solution. Abaqus/Standard provides no check on the initial time increment defined for a mass diffusion analysis; you must ensure that the given value does not violate the above criterion.

In transient analysis using first-order elements the solubility terms are lumped, which eliminates such oscillations but can lead to locally inaccurate solutions for small time increments. If smaller time increments are required, a finer mesh should be used in regions where the normalized concentration changes occur.

Generally there is no upper limit on the time increment because the integration procedure is unconditionally stable unless nonlinearities cause numerical problems.

Automatic Incrementation

The automatic time incrementation scheme for mass diffusion problems is based on the user-specified maximum normalized concentration change allowed at any node during an increment, ΔCmax.

Fixed Time Incrementation

If you choose fixed time incrementation, fixed time increments equal to the size of the user-specified initial time increment, Δt0, will be used.

Ending a Transient Analysis

Transient mass diffusion analysis can be terminated by completing a specified time period, or it can be continued until steady-state conditions are reached. By default, the analysis will end when the given time period has been completed. Alternatively, you can specify that the analysis will end when steady state is reached or the time period ends, whichever comes first. Steady state is defined as the point in time when all normalized concentrations change at less than a user-defined rate.

Initial Conditions

An initial normalized concentration of the diffusing material at specific nodes that belong to mass diffusion elements can be defined (Initial Conditions). For an analysis in which mass diffusion is driven by gradients of temperature and/or pressure (Diffusivity), the initial temperature and pressure stress fields in a model can also be defined.

Boundary Conditions

Boundary conditions can be applied to nodal degree of freedom 11 in any mass diffusion element to prescribe values of normalized concentration (Boundary Conditions). Such values can be specified as functions of time.

Any boundary condition changes to be applied during a mass diffusion step should be given in the respective step using appropriate amplitude definitions to specify their “time” variations (Amplitude Curves). If boundary conditions are specified for the step without amplitude references, they are assumed to change either linearly with “time” during the step or instantly at the start of the step, according to the user-specified or default time variation associated with the step (see Defining an Analysis).

Loads

Concentration fluxes are the only loads that can be applied in a mass diffusion analysis step.

Modifying or Removing Concentration Fluxes

Concentrated or distributed concentration fluxes can be added, modified, or removed as described in About Loads.

Specifying Time-Dependent Concentration Fluxes

The magnitude of a concentrated or a distributed concentration flux can be controlled by referring to an amplitude curve (see Amplitude Curves). If different magnitude variations are needed for different fluxes, the flux definitions can be repeated, with each referring to its own amplitude curve.

Defining Nonuniform Distributed Concentration Fluxes in a User Subroutine

To define nonuniform distributed concentration fluxes, the variation of the flux magnitude throughout a step can be defined in user subroutine DFLUX. If a reference flux magnitude is specified directly, it will be ignored. As a result, any amplitude reference in the flux definition is also ignored.

Predefined Fields

Predefined temperatures, equivalent pressure stresses, and field variables can be specified in a mass diffusion analysis.

Prescribing Temperatures

Temperatures are applied to nodes in temperature-driven mass diffusion analyses by defining a temperature field; absolute zero on the temperature scale used is defined as described in Specifying the Value of Absolute Zero. Alternatively, the temperature field can be obtained from a previous heat transfer analysis. Time-dependent temperature variations are possible with either approach.

A simple interface is provided that uses the Abaqus/Standard results file from the heat transfer analysis to define the temperature field at different times in the mass diffusion analysis. Abaqus/Standard assumes that the nodes in the mass diffusion analysis have the same numbers as the nodes in the previous heat transfer analysis. Values in the results file are ignored at nodes that exist in the heat transfer analysis but not in the mass diffusion analysis, and the temperatures at nodes that did not exist in the heat transfer analysis will not be set by reading the results file.

For specific details on prescribing temperatures, see Predefined Temperature.

Prescribing Equivalent Pressure Stresses

Equivalent pressure stress values can be given at nodes by specifying them directly as a predefined field in the mass diffusion analysis or indirectly by reading the equivalent pressure stresses from the results file of a previous stress/displacement, fully coupled temperature-displacement, or fully coupled thermal-electrical-structural analysis. Regardless of the manner in which they are specified, pressures should be entered according to the Abaqus convention that equivalent pressure stresses are positive when they are compressive.

A simple interface is provided that uses the Abaqus/Standard results file from a mechanical analysis to define the equivalent pressure stresses at different times in the mass diffusion analysis. Abaqus/Standard assumes that the nodes in the mass diffusion analysis have the same numbers as the nodes in the previous mechanical analysis. Values in the results file are ignored at nodes that exist in the mechanical analysis but not in the mass diffusion analysis, and the pressures at nodes that did not exist in the mechanical analysis will not be set by reading the results file.

For specific details on prescribing equivalent pressure stresses, see Predefined Pressure Stress.

Specifying Predefined Field Variables

You can specify values of predefined field variables during a mass diffusion analysis. These values affect only field-variable-dependent material properties, if any. See Predefined Field Variables.

Material Options

Both diffusivity (Diffusivity) and solubility (Solubility) must be defined in a mass diffusion analysis. Optionally, a Soret effect factor and a pressure stress factor can be defined to introduce mass diffusion caused by temperature and pressure gradients, respectively. The use of Fick's law also introduces temperature-driven mass diffusion since a Soret effect factor is calculated automatically.

Elements

Mass diffusion analysis can be performed using only the two-dimensional, three-dimensional, and axisymmetric solid elements that are included in the Abaqus/Standard heat transfer/mass diffusion element library.

Output

In addition to the standard output identifiers available in Abaqus/Standard (Abaqus/Standard Output Variable Identifiers), the following variables have special meaning in mass diffusion analyses:

Element integration point variables:

CONC

Mass concentration.

ISOL

Amount of solute at the integration point, calculated as the product of the mass concentration and the integration point volume.

MFL

Magnitude and components of the concentration flux vector (excluding the terms due to pressure and temperature gradients).

MFLM

Magnitude of the concentration flux vector.

MFLn

Component n of the concentration flux vector (n = 1, 2, 3).

TEMP

Magnitude of the applied temperature field.

Whole element variables:

ESOL

Amount of solute in the element, calculated as the sum of ISOL over all the element integration points.

NFLUX

Fluxes at the nodes of the element caused by mass diffusion in the element.

FLUXS

Distributed mass flux applied to an element.

Whole or partial model variables:

SOL

Amount of solute in the model or specified element set, calculated as the sum of ESOL over all the elements in the model or set.

Nodal variables:

CFL

All concentrated flux values.

CFLn

Concentrated flux value n at a node (n = 11).

NNC

All normalized concentration values at a node.

NNCn

Normalized concentration degree of freedom n at a node (n = 11).

RFL

All reaction flux values (conjugate to normalized concentration).

RFLn

Reaction flux value n at a node (n = 11) (conjugate to normalized concentration).

Input File Template

The following template is representative of a three-step mass diffusion analysis. The first step establishes an initial steady-state concentration distribution of a diffusing material. In the second step equivalent pressure stresses are read from a fully coupled temperature-displacement analysis and the transient mass diffusion response is obtained for the case of mechanical loading of the body. In the final step a temperature field is read from a fully coupled temperature-displacement analysis and the transient mass diffusion response is calculated for the case of heating and cooling the body in which diffusion occurs. An example problem that follows this template is Thermomechanical diffusion of hydrogen in a bending beam.

HEADINGMATERIAL,NAME=mat1
SOLUBILITY
Data lines to define solubility
DIFFUSIVITY
Data lines to define diffusivity
KAPPA,TYPE=TEMP
Data lines to define diffusion driven by temperature gradients
KAPPA,TYPE=PRESS
Data lines to define diffusion driven by gradients of equivalent pressure stress
INITIAL CONDITIONS,TYPE=TEMPERATURE
Data lines to define an initial temperature field
INITIAL CONDITIONS,TYPE=CONCENTRATION
Data lines to define initial nodal values of normalized concentration
INITIAL CONDITIONS,TYPE=PRESSURE STRESS
Data lines to define initial nodal values of equivalent pressure stress
AMPLITUDE,NAME=name
Data lines to define amplitude variations
**
STEP
Step 1 - steady-state solution
MASS DIFFUSION,STEADY STATE
Data line to define incrementation
BOUNDARY
Data lines to prescribe nodal values of normalized concentration
EL FILE
Data lines to define element integration output to the results file
NODE FILE
Data lines to define nodal output to the results file
END STEP
**
STEP
Step 2 - transient analysis driven by pressure stress gradients
MASS DIFFUSION,DCMAX=dcmax,END=SS
Data line to define incrementation
BOUNDARY
Data lines to prescribe nodal values of normalized concentration
PRESSURE STRESS,FILE=name
EL FILE
Data lines to define element integration output to the results file
NODE FILE
Data lines to define nodal output to the results file
END STEP
**
STEP
Step 3 - transient analysis driven by temperature gradients
MASS DIFFUSION,DCMAX=dcmax,END=SS
Data line to define incrementation
BOUNDARY
Data lines to prescribe nodal values of normalized concentration
TEMPERATURE,FILE=name
EL FILE
Data lines to define element integration output to the results file
NODE FILE
Data lines to define nodal output to the results file
END STEP

References

  1. Crank J. The Mathematics of Diffusion, Clarendon Press , Oxford , 1956.