Documentation

This is machine translation

Translated by Microsoft
Mouse over text to see original. Click the button below to return to the English verison of the page.

assempde

Assemble finite element matrices and solve elliptic PDE

assempde is not recommended. Use solvepde instead.

Syntax

  • [K,M,F,Q,G,H,R] = assempde(___)
    example
  • [K,M,F,Q,G,H,R] = assempde(___,[],sdl)
  • u = assempde(K,M,F,Q,G,H,R)
    example
  • [Ks,Fs] = assempde(K,M,F,Q,G,H,R)
  • [Kc,Fc,B,ud] = assempde(K,M,F,Q,G,H,R)

Description

example

u = assempde(model,c,a,f) solves the PDE

(cu)+au=f

with geometry, boundary conditions, and finite element mesh in model, and coefficients c, a, and f. If the PDE is a system of equations (model.PDESystemSize > 1), then assempde solves the system of equations

(cu)+au=f

example

u = assempde(b,p,e,t,c,a,f) solves the PDE with boundary conditions b, and finite element mesh (p,e,t).

example

[Kc,Fc,B,ud] = assempde(___), for any of the previous input syntaxes, assembles finite element matrices using the reduced linear system form, which eliminates any Dirichlet boundary conditions from the system of linear equations. You can calculate the solution u at node points by the command u = B*(Kc\Fc) + ud. See Reduced Linear System.

example

[Ks,Fs] = assempde(___) assembles finite element matrices that represent any Dirichlet boundary conditions using a stiff-spring approximation. You can calculate the solution u at node points by the command u = Ks\Fs. See Stiff-Spring Approximation.

example

[K,M,F,Q,G,H,R] = assempde(___) assembles finite element matrices that represent the PDE problem. This syntax returns all the matrices involved in converting the problem to finite element form. See Algorithms.

[K,M,F,Q,G,H,R] = assempde(___,[],sdl) restricts the finite element matrices to those that include the subdomain specified by the subdomain labels in sdl. The empty argument is required in this syntax for historic and compatibility reasons.

example

u = assempde(K,M,F,Q,G,H,R) returns the solution u based on the full collection of finite element matrices.

[Ks,Fs] = assempde(K,M,F,Q,G,H,R) returns finite element matrices that approximate Dirichlet boundary conditions using the stiff-spring approximation. See Algorithms.

[Kc,Fc,B,ud] = assempde(K,M,F,Q,G,H,R) returns finite element matrices that eliminate any Dirichlet boundary conditions from the system of linear equations. See Algorithms.

Examples

collapse all

Solve an elliptic PDE on an L-shaped region.

Create a scalar PDE model. Incorporate the geometry of an L-shaped region.

model = createpde;
geometryFromEdges(model,@lshapeg);

Apply zero Dirichlet boundary conditions to all edges.

applyBoundaryCondition(model,'Edge',1:model.Geometry.NumEdges,'u',0);

Generate a finite element mesh.

generateMesh(model);

Solve the PDE $- \nabla \cdot \left( c \nabla u \right) + au = f$ with parameters c = 1, a = 0, and f = 5.

c = 1;
a = 0;
f = 5;
u = assempde(model,c,a,f);

Plot the solution.

pdeplot(model,'XYData',u)

Solve a 3-D elliptic PDE using a PDE model.

Create a PDE model container, import a 3-D geometry description, and view the geometry.

model = createpde;
importGeometry(model,'Block.stl');
pdegplot(model,'FaceLabels','on','FaceAlpha',0.5)

Set zero Dirichlet conditions on faces 1 through 4 (the edges). Set Neumann conditions with g = -1 on face 6 and g = 1 on face 5.

applyBoundaryCondition(model,'Face',1:4,'u',0);
applyBoundaryCondition(model,'Face',6,'g',-1);
applyBoundaryCondition(model,'Face',5,'g',1);

Set coefficients c = 1, a = 0, and f = 0.1.

c = 1;
a = 0;
f = 0.1;

Create a mesh and solve the problem.

generateMesh(model);
u = assempde(model,c,a,f);

Plot the solution on the surface.

pdeplot3D(model,'ColorMapData',u)

Solve a 2-D PDE using the older syntax for mesh.

Create a circle geometry.

g = @circleg;

Set zero Dirichlet boundary conditions.

b = @circleb1;

Create a mesh for the geometry.

[p,e,t] = initmesh(g);

Solve the PDE $- \nabla \cdot \left( c \nabla u \right) + au = f$ with parameters c = 1, a = 0, and f = sin(x).

c = 1;
a = 0;
f = 'sin(x)';
u = assempde(b,p,e,t,c,a,f);

Plot the solution.

pdeplot(p,e,t,'XYData',u)

Obtain the finite-element matrices that represent the problem using a reduced linear algebra representation of Dirichlet boundary conditions.

Create a scalar PDE model. Import a simple 3-D geometry.

model = createpde;
importGeometry(model,'Block.stl');

Set zero Dirichlet boundary conditions on all the geometry faces.

applyBoundaryCondition(model,'dirichlet','Face',1:model.Geometry.NumFaces,'u',0);

Generate a mesh for the geometry.

generateMesh(model);

Obtain finite element matrices K, F, B, and ud that represent the equation $- \nabla \cdot \left( c \nabla u \right) + au = f$ with parameters $c = 1$, $a = 0$, and $f = \log \left(1 + x + \frac{y}{1 + z} \right)$.

c = 1;
a = 0;
f = 'log(1+x+y./(1+z))';
[K,F,B,ud] = assempde(model,c,a,f);

You can obtain the solution u of the PDE at mesh nodes by executing the command

u = B*(K\F) + ud;

Generally, this solution is slightly more accurate than the stiff-spring solution, as calculated in the next example.

Obtain the stiff-spring approximation of finite element matrices.

Create a scalar PDE model. Import a simple 3-D geometry.

model = createpde;
importGeometry(model,'Block.stl');

Set zero Dirichlet boundary conditions on all the geometry faces.

applyBoundaryCondition(model,'Face',1:model.Geometry.NumFaces,'u',0);

Generate a mesh for the geometry.

generateMesh(model);

Obtain finite element matrices Ks and Fs that represent the equation $- \nabla \cdot \left( c \nabla u \right) + au = f$ with parameters $c = 1$, $a = 0$, and $f = \log \left(1 + x + \frac{y}{1 + z} \right)$.

c = 1;
a = 0;
f = 'log(1+x+y./(1+z))';
[Ks,Fs] = assempde(model,c,a,f);

You can obtain the solution u of the PDE at mesh nodes by executing the command

u = Ks\Fs;

Generally, this solution is slightly less accurate than the reduced linear algebra solution, as calculated in the previous example.

Obtain the full collection of finite element matrices for an elliptic problem.

Import geometry and set up an elliptic problem with Dirichlet boundary conditions. The Torus.stl geometry has only one face, so you need set only one boundary condition.

model = createpde();
importGeometry(model,'Torus.stl');
applyBoundaryCondition(model,'face',1,'u',0);
c = 1;
a = 0;
f = 1;
generateMesh(model);

Create the finite element matrices that represent this problem.

[K,M,F,Q,G,H,R] = assempde(model,c,a,f);

Most of the resulting matrices are quite sparse. G, M, Q, and R are all zero sparse matrices.

howsparse = @(x)nnz(x)/numel(x);
disp(['Maximum fraction of nonzero entries in K or H is ',...
       num2str(max(howsparse(K),howsparse(H)))])
Maximum fraction of nonzero entries in K or H is 0.0008067

To find the solution to the PDE, call assempde again.

u = assempde(K,M,F,Q,G,H,R);

Input Arguments

collapse all

PDE model, specified as a PDEModel object.

Example: model = createpde(1)

PDE coefficient, specified as a scalar or matrix, as a character array, or as a coefficient function. c represents the c coefficient in the scalar PDE

(cu)+au=f

or in the system of PDEs

(cu)+au=f

You can specifyc in various ways, detailed in c Coefficient for Systems. See also Specify Scalar PDE Coefficients in Character Form, Specify 2-D Scalar Coefficients in Function Form, and Specify 3-D PDE Coefficients in Function Form.

Example: 'cosh(x+y.^2)'

Data Types: double | char | function_handle
Complex Number Support: Yes

PDE coefficient, specified as a scalar or matrix, as a character array, or as a coefficient function. a represents the a coefficient in the scalar PDE

(cu)+au=f

or in the system of PDEs

(cu)+au=f

You can specifya in various ways, detailed in a or d Coefficient for Systems. See also Specify Scalar PDE Coefficients in Character Form, Specify 2-D Scalar Coefficients in Function Form, and Specify 3-D PDE Coefficients in Function Form.

Example: 2*eye(3)

Data Types: double | char | function_handle
Complex Number Support: Yes

PDE coefficient, specified as a scalar or matrix, as a character array, or as a coefficient function. f represents the f coefficient in the scalar PDE

(cu)+au=f

or in the system of PDEs

(cu)+au=f

You can specifyf in various ways, detailed in f Coefficient for Systems. See also Specify Scalar PDE Coefficients in Character Form, Specify 2-D Scalar Coefficients in Function Form, and Specify 3-D PDE Coefficients in Function Form.

Example: char('sin(x)';'cos(y)';'tan(z)')

Data Types: double | char | function_handle
Complex Number Support: Yes

Boundary conditions, specified as a boundary matrix or boundary file. Pass a boundary file as a function handle or as a file name.

Example: b = 'circleb1' or equivalently b = @circleb1

Data Types: double | char | function_handle

Mesh points, specified as a 2-by-Np matrix of points, where Np is the number of points in the mesh. For a description of the (p,e,t) matrices, see Mesh Data.

Typically, you use the p, e, and t data exported from the PDE app, or generated by initmesh or refinemesh.

Example: [p,e,t] = initmesh(gd)

Data Types: double

Mesh edges, specified as a 7-by-Ne matrix of edges, where Ne is the number of edges in the mesh. For a description of the (p,e,t) matrices, see Mesh Data.

Typically, you use the p, e, and t data exported from the PDE app, or generated by initmesh or refinemesh.

Example: [p,e,t] = initmesh(gd)

Data Types: double

Mesh triangles, specified as a 4-by-Nt matrix of triangles, where Nt is the number of triangles in the mesh. For a description of the (p,e,t) matrices, see Mesh Data.

Typically, you use the p, e, and t data exported from the PDE app, or generated by initmesh or refinemesh.

Example: [p,e,t] = initmesh(gd)

Data Types: double

Stiffness matrix, specified as a sparse matrix or full matrix. Generally, you obtain K from a previous call to assema or assempde. For the meaning of stiffness matrix, see Elliptic Equations.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Mass matrix, specified as a sparse matrix or full matrix. Generally, you obtain M from a previous call to assema or assempde. For the meaning of mass matrix, see Elliptic Equations.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Finite element f representation, specified as a vector. Generally, you obtain F from a previous call to assema or assempde. For the meaning of this representation, see Elliptic Equations.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Neumann boundary condition matrix, specified as a sparse matrix or full matrix. Generally, you obtain Q from a previous call to assemb or assempde. For the meaning of this matrix, see Elliptic Equations.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Neumann boundary condition vector, specified as a sparse vector or full vector. Generally, you obtain G from a previous call to assemb or assempde. For the meaning of this vector, see Elliptic Equations.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Dirichlet boundary condition matrix, specified as a sparse matrix or full matrix. Generally, you obtain H from a previous call to assemb or assempde. For the meaning of this matrix, see Algorithms.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Dirichlet boundary condition vector, specified as a sparse vector or full vector. Generally, you obtain R from a previous call to assemb or assempde. For the meaning of this vector, see Algorithms.

Example: [K,M,F,Q,G,H,R] = assempde(model,c,a,f)

Data Types: double
Complex Number Support: Yes

Subdomain labels, specified as a vector of positive integers. For 2-D geometry only. View the subdomain labels in your geometry using the command

pdegplot(g,'SubdomainLabels','on')

Example: sdl = [1,3:5];

Data Types: double

Output Arguments

collapse all

PDE solution, returned as a vector.

  • If the PDE is scalar, meaning only one equation, then u is a column vector representing the solution u at each node in the mesh. u(i) is the solution at the ith column of model.Mesh.Nodes or the ith column of p.

  • If the PDE is a system of N > 1 equations, then u is a column vector with N*Np elements, where Np is the number of nodes in the mesh. The first Np elements of u represent the solution of equation 1, then next Np elements represent the solution of equation 2, etc.

To obtain the solution at an arbitrary point in the geometry, use pdeInterpolant.

To plot the solution, use pdeplot for 2-D geometry, or see Plot 3-D Solutions and Their Gradients.

Stiffness matrix, returned as a sparse matrix. See Elliptic Equations.

u1 = Kc\Fc returns the solution on the non-Dirichlet points. To obtain the solution u at the nodes of the mesh,

u = B*(Kc\Fc) + ud

Generally, Kc, Fc, B, and ud make a slower but more accurate solution than Ks and Fs.

Load vector, returned as a vector. See Elliptic Equations.

u = B*(Kc\Fc) + ud

Generally, Kc, Fc, B, and ud make a slower but more accurate solution than Ks and Fs.

Dirichlet nullspace, returned as a sparse matrix. See Algorithms.

u = B*(Kc\Fc) + ud

Generally, Kc, Fc, B, and ud make a slower but more accurate solution than Ks and Fs.

Dirichlet vector, returned as a vector. See Algorithms.

u = B*(Kc\Fc) + ud

Generally, Kc, Fc, B, and ud make a slower but more accurate solution than Ks and Fs.

Finite element matrix for stiff-spring approximation, returned as a sparse matrix. See Algorithms.

To obtain the solution u at the nodes of the mesh,

u = Ks\Fs.

Generally, Ks and Fs make a quicker but less accurate solution than Kc, Fc, B, and ud.

Load vector corresponding to the stiff-spring approximation for Dirichlet boundary condition, returned as a vector. See Algorithms.

To obtain the solution u at the nodes of the mesh,

u = Ks\Fs.

Generally, Ks and Fs make a quicker but less accurate solution than Kc, Fc, B, and ud.

Stiffness matrix, returned as a sparse matrix. See Elliptic Equations.

K represents the stiffness matrix alone, unlike Kc or Ks, which are stiffness matrices combined with other terms to enable immediate solution of a PDE.

Typically, you use K in a subsequent call to a solver such as assempde or hyperbolic.

Mass matrix. returned as a sparse matrix. See Elliptic Equations.

Typically, you use M in a subsequent call to a solver such as assempde or hyperbolic.

Load vector, returned as a vector. See Elliptic Equations.

F represents the load vector alone, unlike Fc or Fs, which are load vectors combined with other terms to enable immediate solution of a PDE.

Typically, you use F in a subsequent call to a solver such as assempde or hyperbolic.

Neumann boundary condition matrix, returned as a sparse matrix. See Elliptic Equations.

Typically, you use Q in a subsequent call to a solver such as assempde or hyperbolic.

Neumann boundary condition vector, returned as a sparse vector. See Elliptic Equations.

Typically, you use G in a subsequent call to a solver such as assempde or hyperbolic.

Dirichlet matrix, returned as a sparse matrix. See Algorithms.

Typically, you use H in a subsequent call to a solver such as assempde or hyperbolic.

Dirichlet vector, returned as a sparse vector. See Algorithms.

Typically, you use R in a subsequent call to a solver such as assempde or hyperbolic.

More About

collapse all

Reduced Linear System

This form of the finite element matrices eliminates Dirichlet conditions from the problem using a linear algebra approach. The finite element matrices reduce to the solution u = B*(Kc\Fc) + ud, where B spans the null space of the columns of H (the Dirichlet condition matrix representing hu = r). R is the Dirichlet condition vector for Hu = R. ud is the vector of boundary condition solutions for the Dirichlet conditions. u1 = Kc\Fc returns the solution on the non-Dirichlet points.

See Systems of PDEs for details on the approach used to eliminate Dirichlet conditions.

Stiff-Spring Approximation

This form of the finite element matrices converts Dirichlet boundary conditions to Neumann boundary conditions using a stiff-spring approximation. Using this approximation, assempde returns a matrix Ks and a vector Fs that represent the combined finite element matrices. The approximate solution u is u = Ks\Fs.

See Elliptic Equations. For details of the stiff-spring approximation, see Systems of PDEs.

Algorithms

assempde performs the following steps to obtain a solution u to an elliptic PDE:

  1. Generate the finite element matrices [K,M,F,Q,G,H,R]. This step is equivalent to calling assema to generate the matrices K, M, and F, and also calling assemb to generate the matrices Q, G, H, and R.

  2. Generate the combined finite element matrices [Kc,Fc,B,ud]. 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 boundary condition solutions for the Dirichlet conditions.

  3. Calculate the solution u via

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

assempde uses one of two algorithms for assembling a problem into combined finite element matrix form. A reduced linear system form leads to immediate solution via linear algebra. You choose the algorithm by the number of outputs. For the reduced linear system form, request four outputs:

[Kc,Fc,B,ud] = assempde(_)

For the stiff-spring approximation, request two outputs:

[Ks,Fs] = assempde(_)

For details, see Reduced Linear System and Stiff-Spring Approximation.

As explained in Elliptic Equations, the full finite element matrices and vectors are the following.

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

  • M is the mass matrix, 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. See Systems of PDEs.

Introduced before R2006a

Was this topic helpful?