Dr A. Cookson Department of Mechanical Engineering
20192020
ME50344: Engineering Systems Simulation
University of Bath
Assignment 1: Static MATLABBased FEM Modelling
Summary
You must develop a FEMbased simulation tool to solve the static diffusion reaction equation in a 1D mesh, which your consultancy firm can reuse 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.
The specific tasks of this deliverable are as follows:
Part 1: Software VerificationAnalytical Testing
a. Create a MATLAB function to calculate the local 2by2 element matrix for the diffusion operator, for any element, eN, in the finite element mesh. 2
The inputs to the function are:
the diffusion coefficient, D
the local element number, eN
the mesh data structure through which you can access each elements
Jacobian, J and nodal positions x0 and x1
Create the mesh using the function OneDimLinearMeshGen.m , and store
it using the mesh data structure, as provided in Tutorial 2.
Your function should pass the unit test, CourseworkOneUnitTest.m, which is available from the courses 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 functions source code
b. CreateaMATLABfunctiontocalculatethelocal2by2elementmatrixforthe linear reaction operator, for any element, eN, in the finite element mesh. 3
The inputs to the function are:
the reaction coefficient, lambda
Page 1 of 6
Dr A. Cookson Department of Mechanical Engineering 20192020 University of Bath
the local element number, eN
the mesh data structure through which you can access each elements
Jacobian, J and nodal positions x0 and x1
Create your own unit test for this function and explain in your report what your unit test has been 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
the functions source code
the unit tests source code
c. Using these functions, develop a finite element solver for the static reaction diffusion equation, using additional unit tests to verify any new functions as appropriate.
Demonstrate that the code is correct using the following two analytical test cases:
Analytical Test OneLaplaces equation: Solve Laplaces equation:
2c x2
for the fourelement mesh of the domain between x0 and x1, shown in Figure 1.
i. Solve for the case of the two Dirichlet boundary conditions 25:
c2 at x0 c0 at x1
Compare your numerical solution to the analytical solution for this equationboundary conditions, which is:
c21x Page 2 of 6
0
Dr A. Cookson Department of Mechanical Engineering
20192020
How accurate is your solution? Why is this?
University of Bath
ii. Run your code again for Laplaces equation, this time with the following Dirichlet boundary condition and Neumann boundary condition 5:
c2 at x0 x
c0 at x1
Explain the effect that changing this boundary condition has had on the
solution 5.
Figure 1: Uniform 4element mesh defined between x0 and x1.
Analytical Test TwoDiffusionReaction equation
iii. Now check your code solves reaction terms correctly, by using it to
Page 3 of 6
solve the diffusionreaction equation 10: D 2cc0
x2
for the following parameters:
D1, 9
and Dirichlet boundary conditions:
c0 at x0 c1 at x1
Dr A. Cookson Department of Mechanical Engineering 20192020 University of Bath
For these parametersboundary conditions that equation has the following analytical solution:
cxe3 e3xe3x e61
Use several different mesh resolutions and check that your solution approaches the analytical value as you increase the number of elements.
Part 2: ModellingSimulation Results
Thermophoresis is a phenomenon in which particles of different types diffuse differently through a medium depending on the local temperature gradient that they are exposed to. The particles diffuse at a rate that is negatively proportional to the local temperature gradient i.e. the steeper the negative temperature gradient, the faster the particles will diffuse. In biomaterial processing this effect can be harnessed to create a spatially varying concentration of these particles, which causes a spatially varying stiffness to be created in the biomaterial. This is useful for growing cell cultures under biologically realistic conditions.
In Part 2 you will use your code to model steadystate heat transfer through a material that is filled with small diameter heating channels. By varying the temperature of the heating liquid, the flow rate, or the spatial distribution of temperature in these channels, different temperature profiles can be obtained for the thermophoresisbased manufacturing process. Your task is to characterise the behaviour of this system for its different parameters.
Figure 2 shows a crosssection of the material:
x0.01m x
x0m
TINNER
TL Q TOUTER
Figure 2: Crosssection of the material
Page 4 of 6
Dr A. Cookson Department of Mechanical Engineering 20192020 University of Bath
Given that the boundary conditions are uniform in all directions, we can reduce the problem to a 1D approximation for heat transfer through the material:
k2T QTL T0 x2
a. Use your code to solve this model and hence compute the temperature distribution T in Kelvin for different combinations of values of Q the liquid flow rate and TL liquid temperature in Kelvin in the ranges 20:
Q0.5 to 1.5
TL 294.15Kto322.15Ki.e.21to49degC
The other material parameters and boundary conditions are fixed as: k1.01e05
DirichletBCatx0mis: TOUTER 323.15K50degC Dirichlet BC at x0.01 m is: TINNER293.15 K 20 deg C.
i. Plot the results from your parameter space study and explain the effect each parameter has on both the temperature distribution and on the gradients of the temperature.
ii. If some particles are initially placed at x0 and are allowed to diffuse for a short amount of time, describe how you think their distribution, and hence the variation in stiffness, might change with Q and TL, based on the effect these have on the temperature gradients in the material?
iii. Assess the effect of different mesh resolutions and explain how this affects the accuracy of the results and your conclusions.
b. Itissuggestedthatvaryingthetemperatureoftheliquidfromonechannelto the next could provide greater control of the temperature distribution. This variation in temperature and its effect on heating performance can be represented by the following linear function of x for the source term 10:
TL14x
This means the equation to be solved is now:
k2T QTL14xT0 x2
Page 5 of 6
Dr A. Cookson Department of Mechanical Engineering 20192020 University of Bath
Remember that x can be written as the sum of the two basis functions, within each local element:
xx0 0x1 1
i. Derive the analytical expression for the local element vector that represents this modified source term and include it in your report
ii. Implement the modified source term in your code and investigateexplain how it changes the behaviour of this system
Presentation 10
Clear graphical presentation of equations, text, and results
Quality of English, proper use of references, figures, and tables
Code Quality 10
Correctness, elegance, and readability, of MATLAB code
Evidence of software verification and other good development practice
SUBMISSION GUIDELINES
Structure your report as a set of answers to these questionsthere is no requirement to write this in a lab report format. You also do not need to re explain the entirety of the finite element method. However, your report must be selfcontained 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 Appendicesfailure to do so will cause you to lose marks. Do not paste your code into the document as an OLE item or as an image.
Do not upload archivedzippedcompressed folders of these source files.
Do not use MATLABs symbolic algebra toolbox.
Your code should use meaningful variable names and include comments,
in line with good practice.
Word limit of 2000 words not including source code.
Submit your work using the online submission function on the units Moodle page.
Deadline: 4pm on Wed 6th November 2019.
Page 6 of 6
Reviews
There are no reviews yet.