Dr A. Cookson Department of Mechanical Engineering
2020/21
Page 1 of 4
Assignment 1: Static MATLAB-Based FEM Modelling
ME40064: System Modelling & Simulation
ME50344: Engineering Systems Simulation

You must develop a FEM-based simulation tool to solve the static diffusion-reaction equation
in a 1D mesh, which your consultancy firm can re-use for future modelling projects. This tool
must be developed following good software development practice and must be validated
against standard test cases, accompanied by appropriate documentary proof.
Part 1: Solving Laplace’s Equation using FEM

Local Element Matrix Functions [5%]

a. Create a MATLAB function to calculate the local 2-by-2 element matrix for the diffusion
operator, for any element, eN, in the finite element mesh.

The inputs to the function are:
• the diffusion coefficient, D
• the local element number, eN
• the mesh data structure

You should use the mesh data structure to access each element’s Jacobian, J, and its local
nodal positions x0 and x1 . Create the mesh data structure using the function
OneDimLinearMeshGen.m.

Your function should pass the unit test, CourseworkOneUnitTest.m, which is
available on the unit Moodle page. In this script you will see a suggested name for your
function, as well as the required arguments for your function.

Include the following in your report:
• a legible screenshot that shows your function passes the test
• the function’s source code

b. Create a MATLAB function to calculate the local 2-by-2 element matrix for the linear
reaction operator, for any element, eN, in the finite element mesh.

The inputs to the function are:
• the reaction coefficient,
• the local element number, eN
• the mesh data structure

Write unit tests for this function. Explain tests are designed to check and why these tests
are sufficient to show that your function works correctly.

Include the following in your report:
• a legible screenshot that shows your function passes the test
Dr A. Cookson Department of Mechanical Engineering
2020/21
Page 2 of 4
• the function’s source code
• the unit test’s source code

Laplace’s Equation [35%]

Using these functions, now develop a finite element code to solve Laplace’s equation, using
additional unit tests to verify any new functions, as appropriate. Demonstrate that the code is
correct using the following test cases.

Solve Laplace’s equation:



2
2
= 0

Eq. 1
for the four-element mesh of the domain between x=0 and x=1, shown in Figure 1.

i. Solve for the case of the two Dirichlet boundary conditions:

( = 0) = 2
( = 1) = 0

Compare your numerical solution to the analytical solution (Eq. 2) for Laplace’s
equation (Eq. 1) with these boundary conditions:



= 2(1 − )

Eq. 2
How accurate is your solution? Why is this?

ii. Run your code again for Laplace’s equation (Eq. 1), this time with the following
Neumann boundary condition and Dirichlet boundary condition:



( = 0) = 2

( = 1) = 0

Explain the effect that changing this boundary condition has had on the solution and
therefore show that your code is correct.

Figure 1: Uniformly sized four-element mesh for the domain x=0 to x=1. Global
mesh nodes are numbered N1 to N5.
Dr A. Cookson Department of Mechanical Engineering
2020/21
Page 3 of 4
Part 2: Verification of Diffusion-Reaction Equation Solver

Extend your code to be able to solve the full steady state diffusion-reaction equation with
source terms. You should aim to:
• provide unit tests for additional functions
• test different mesh resolutions & domain sizes, where appropriate
• test different boundary conditions and combinations thereof
• explore a range of parameter values

Diffusion-Reaction Equation [10%]

Check that your code handles reaction terms correctly, by using it to solve the diffusion-
reaction equation:



2
2
+ = 0

Eq. 3
for the following parameter values:
= 1, = −9

and Dirichlet boundary conditions:
( = 0) = 0
( = 1) = 1

Note, that for these parameters & boundary conditions, Eq. 3 has the following analytical
solution:


() =
3
6 − 1
(3 − −3)

Eq. 4
Use a range of different mesh resolutions and evaluate how well your solution converges to
the analytical function (Eq. 4) as you increase the number of elements in the mesh.

Poisson’s Equation [10%]

Now verify that your code can correctly solve Poisson’s equation (Eq. 5), that is, equations of
the form:



2
2
=

Eq. 5
Where D is the diffusion coefficient and f is a spatially constant source term. You should test
your code for a range of values of D, f and boundary conditions.

Reaction-Diffusion Equation with Spatially Varying Source [10%]

Solve the reaction-diffusion equation (Eq. 6) with a spatially varying source term:
Dr A. Cookson Department of Mechanical Engineering
2020/21
Page 4 of 4



2
2
+ =

Eq. 6
Where D is the diffusion coefficient, = , and and k are constants.

i. Include in your report the derivation of the analytical function for the local element
vector, which represents this spatially varying source term.
ii. Implement the modified source term in your code and investigate & explain how the
results vary for different values of D, and k, and boundary conditions.
Additional Marking Criteria

Presentation [10%]
• Clear graphical presentation of equations, text, diagrams, and plots
• Standard of written English and explanations
• Proper use of references, figure numbering & captions, and table headings
• Clear layout of report

Code Quality [10%]
• Correctness and elegance
• Good software practice e.g. use of functions, data structures etc
• Readability – use of meaningful variable names and code commenting

Verification [10%]
• Thoroughness of software testing and verification
• Analysis of accuracy of numerical solutions
SUBMISSION GUIDELINES

You may structure your report as a set of answers to these questions – there is no
requirement to write this in a lab report format. Please note, you do not need to re-explain the
entirety of the finite element method in order to derive. However, your report must be self-
contained and therefore must not assume that the reader knows the content in this document.

• You must include all your MATLAB source code as text in the Appendices – failure to
do so will cause you to lose marks i.e.
- Do not paste your code into the document as an OLE item or as an image.
- Do not upload archived/zipped/compressed folders of these source files.
• Do not use MATLAB’s symbolic algebra toolbox – you will lose marks.
• Word limit of 1500 words (not including source code).
• You are reminded of the university’s policies regarding academic plagiarism: your code
and your report must be your own work.

Submit your work via the online submission link on the unit Moodle page.

Deadline: 4pm on Friday 6th November 2020.

欢迎咨询51作业君