assembleFEMatrices

Assemble finite element matrices

Syntax

FEM = assembleFEMatrices(model)

FEM = assembleFEMatrices(model,bcmethod)

Description

FEM = assembleFEMatrices(model) returns a structural array containing finite element matrices. Model attributes, such as coefficients, material properties, boundary conditions, and so on, must not depend on time or solution.

FEM = assembleFEMatrices(model,bcmethod) assembles finite element matrices and imposes boundary conditions using the method specified by bcmethod.

Examples

Finite Element Matrices for 2-D Problem

Create a PDE model for the Poisson equation on an L-shaped membrane with zero Dirichlet boundary conditions.

model = createpde(1);
geometryFromEdges(model,@lshapeg);
specifyCoefficients(model,'m',0,'d',0,'c',1,'a',0,'f',1);
applyBoundaryCondition(model,'edge',1:model.Geometry.NumEdges,'u',0);

Generate a mesh and obtain the default finite element matrices for the problem and mesh.

generateMesh(model,'Hmax',0.2);
FEM = assembleFEMatrices(model)

FEM = struct with fields:
    K: [401x401 double]
    A: [401x401 double]
    F: [401x1 double]
    Q: [401x401 double]
    G: [401x1 double]
    H: [80x401 double]
    R: [80x1 double]
    M: [401x401 double]

Finite Element Matrices with nullspace Method

Create a PDE model for the Poisson equation on an L-shaped membrane with zero Dirichlet boundary conditions.

model = createpde(1);
geometryFromEdges(model,@lshapeg);
specifyCoefficients(model,'m',0,'d',0,'c',1,'a',0,'f',1);
applyBoundaryCondition(model,'edge',1:model.Geometry.NumEdges,'u',0);

Generate a mesh and obtain the nullspace finite element matrices for the problem and mesh.

generateMesh(model,'Hmax',0.2);
FEM = assembleFEMatrices(model,'nullspace')

FEM = struct with fields:
    Kc: [321x321 double]
    Fc: [321x1 double]
     B: [401x321 double]
    ud: [401x1 double]
     M: [321x321 double]

Obtain the solution to the PDE.

u = FEM.B*(FEM.Kc\FEM.Fc) + FEM.ud;

Compare this result to the solution given by solvepde. The two solutions are identical.

u1 = solvepde(model);
norm(u - u1.NodalSolution)

ans = 0

Finite Element Matrices for Thermal Model

Assemble the intermediate finite element matrices for a thermal problem.

Create a transient thermal model and include the geometry of the built-in function squareg.

thermalmodel = createpde('thermal','transient');
geometryFromEdges(thermalmodel,@squareg);

Plot the geometry with the edge labels.

pdegplot(thermalmodel,'EdgeLabels','on')
xlim([-1.1 1.1])
ylim([-1.1 1.1])

Specify the thermal conductivity, mass density, and specific heat of the material.

thermalProperties(thermalmodel,'ThermalConductivity',0.2, ...
                               'MassDensity',2.7*10^(-6), ...
                               'SpecificHeat',920); 

Set the boundary and initial conditions.

thermalBC(thermalmodel,'Edge',1:4,'Temperature',100);
thermalIC(thermalmodel,0,'Face',1);

Generate a mesh and obtain the default finite element matrices.

generateMesh(thermalmodel);
FEM = assembleFEMatrices(thermalmodel)

FEM = struct with fields:
    K: [1541x1541 double]
    A: [1541x1541 double]
    F: [1541x1 double]
    Q: [1541x1541 double]
    G: [1541x1 double]
    H: [144x1541 double]
    R: [144x1 double]
    M: [1541x1541 double]

Input Arguments

model — Model object

PDEModel object | ThermalModel object | StructuralModel object

Model object, specified as a PDEModel object, ThermalModel object, or StructuralModel object.

Example: model = createpde(1)

Example: thermalmodel = createpde('thermal','steadystate')

Example: structuralmodel = createpde('structural','static-solid')

bcmethod — Method for including boundary conditions

'none' (default) | 'nullspace' | 'stiff-spring'

Method for including boundary conditions, specified as 'none', 'nullspace', or 'stiff-spring'. For more information, see Algorithms.

Example: FEM = assembleFEMatrices(model,'nullspace')

Data Types: char | string

Output Arguments

FEM — Finite element matrices

structural array

Finite element matrices, returned as a structural array. The fields in the structural array depend on bcmethod as follows:

• If the value is 'none', then the fields are K, A, F, Q, G, H, R, and M.

• If the value is 'nullspace', then the fields are Kc, Fc, B, ud, and M.

• If the value is 'stiff-spring', then the fields are Ks, Fs, and M.

For more information, see Algorithms.

Tips

• The mass matrix M is nonzero when the model is time-dependent. By using this matrix, you can solve a model with Rayleigh damping. For an example, see Dynamics of Damped Cantilever Beam.

• For a thermal model, the m and a coefficients are zeros. The thermal conductivity maps to the c coefficient. The product of the mass density and the specific heat maps to the d coefficient. The internal heat source maps to the f coefficient.

• For a structural model, the a coefficient is zero. The Young's modulus and Poisson's ratio map to the c coefficient. The mass density maps to the m coefficient. The body loads map to the f coefficient. When you specify the damping model by using the Rayleigh damping parameters Alpha and Beta, the discretized damping matrix C is computed by using the mass matrix M and the stiffness matrix K as C = Alpha*M+Beta*K.

Algorithms

The full finite element matrices and vectors are as follows:

• K is the stiffness matrix, the integral of the c coefficient against the basis functions.

• M is the mass matrix, the integral of the m or d coefficient against the basis functions.

• A is the integral of the a coefficient against the basis functions.

• F is the integral of the f coefficient against the basis functions.

• Q is the integral of the q boundary condition against the basis functions.

• G is the integral of the g boundary condition against the basis functions.

• The H and R matrices come directly from the Dirichlet conditions and the mesh.

Given these matrices, the 'nullspace' technique generates the combined finite element matrices [Kc,Fc,B,ud] as follows. The combined stiffness matrix is for the reduced linear system Kc = K + M + Q. The corresponding combined load vector is Fc = F + G. The B matrix spans the null space of the columns of H (the Dirichlet condition matrix representing hu = r). The R vector represents the Dirichlet conditions in Hu = R. The ud vector represents the boundary condition solutions for the Dirichlet conditions.

From the 'nullspace' matrices, you can compute the solution u as

u = B*(Kc\Fc) + ud.

Note

Internally, for time-independent problems, solvepde uses the 'nullspace' technique, and calculates solutions using u = B*(Kc\Fc) + ud.

The 'stiff-spring' technique returns a matrix Ks and a vector Fs that together represent a different type of combined finite element matrices. The approximate solution u is u = Ks\Fs.

Compared to the 'nullspace' technique, the 'stiff-spring' technique generates matrices more quickly, but generally gives less accurate solutions.