The purpose of this tutorial is to demonstrate the solution of solid mechanics problems using unstructured meshes. We focus primarily on two-dimensional meshes, generated using Jonathan Shewchuk's open-source mesh generator
Triangle, based on the output from the open-source drawing program xfig. An example three-dimensional problem is described here; see also 3D Problems later in this document.
The solid mechanics problem studied here also serves as a "warm-up problem" for the corresponding fluid-structure interaction problem in which the solid object is immersed in (and loaded by) a viscous fluid.
Here is a sketch of the problem: A slightly strange-looking elastic solid is loaded by gravity and by a pressure load, acting on its upper face.
We start by drawing the outline of the solid as a polyline in xfig:
(Note that we draw the solid upside-down because of the way xfig orients its coordinate axes.)
We save the figure as a *.fig file and convert it to a *.poly file using
oomph-lib's conversion code
fig2poly (a copy of which is located in
This creates a file called
solid.fig.poly that can be processed using Triangle. For instance, to create a quality mesh with a maximum element size of 0.025 we use
Here is a plot of the mesh, generated by
showme distributed with
*.poly, *.ele and
*.node files generated by
Triangle can now be used as input to
The animation shown below illustrates the solid's deformation. The first frame shows the undeformed, stress-free reference configuration; the second frame shows the deformation induced by gravity. Subsequent frames illustrate the deformation in response to the (additional) increasing "suction" applied at the upper face.
The blue markers show the position of pinned nodes.
We create the mesh by multiple inheritance from
TriangleMesh and the
SolidMesh base class:
The constructor calls the constructor of the underlying
TriangleMesh, and, as usual, sets the Lagrangian coordinates to the current nodal positions, making the current configuration stress-free.
TriangleMesh constructor associates each polyline in the
xfig drawing with a distinct
oomph-lib mesh boundary. Hence the boundary nodes are initially located on the same, single boundary. To facilitate the application of boundary conditions, we divide the single boundary into three:
We loop over all nodes in the mesh and identify nodes on the lower (pinned) boundary by their y-coordinate. We remove the node from the boundary 0 and re-allocate it to the new boundary 1:
Similarly, we identify all nodes on the upper boundary and re-assign them to boundary 2, before re-generating the various boundary lookup schemes that identify which elements are located next to the various mesh boundaries:
As usual we define the various problem parameters in a global namespace. We define Poisson's ratio and prepare a pointer to a constitutive equation.
Next we define the gravitational body force
and the pressure load to be applied at the upper boundary
The driver code is straightforward. We specify an output directory and instantiate a constitutive equation. (Recall that the single-argument constructor to the
GeneralisedHookean constitutive law implies that all stresses are non-dimensionalised on Young's modulus ).
We create the
Problem object using a displacement formulation of the equations and output the initial configuration
Finally, we perform a straightforward parameter study, applying a constant gravitational load and slowly increasing the suction (negative pressure) on the upper boundary.
The parameter study is then repeated for a pressure-displacement formulation with and without an incompressibility constraint, see the source code for details.
Problem class has the usual member functions and provides storage for the two sub-meshes: the bulk mesh of 2D solid elements and the mesh of 1D traction elements that will be attached to the upper boundary.
We start by building the bulk mesh, using the files created by
Next we create traction elements, attaching them to the "bulk" solid elements that are adjacent to boundary 2. We also specify the load function.
We add both meshes as sub-meshes to the
Problem and build the global mesh
Next we apply the boundary conditions at the lower boundary where we suppress the displacements in both directions. We document the position of the pinned nodes to allow us to check that the boundary IDs were identified correctly – see Comments and Exercises for a further discussion of this issue.
Finally, we complete the build of the solid elements by specifying their constitutive equation and the body force before assigning the equation numbers.
The post-processing routine outputs the deformed domain shape and the applied traction. In the spirit of continuing paranoia we also document the domain boundaries; see Comments and Exercises.
This tutorial demonstrates that the use of unstructured meshes for solid mechanics problems is extremely straightforward. The only aspect that requires some care (and not just for solid mechanics applications) is the correct identification/assignment of domain boundaries when the mesh is generated with xfig. The fact that we documented the mesh boundaries and the positions of the pinned nodes in the driver code suggests (correctly!) that we managed to get both assignments (slightly) wrong when we first wrote the driver code. The manual identification of nodes on domain boundaries is tedious and therefore error prone and, as usual, it pays to be as a paranoid as possible! Ignore this advice at your own risk...
We also welcome any improvements to our conversion code fig2poly.cc that would allow the specification of domain boundaries from within xfig. It's probably not possible but if you have any ideas how to go about this (either by hijacking information that can be generated by xfig or via some other GUI interface), let us know.
Unstructured meshes may also be used for the 3D solid mechanics problems. The overall procedure is very similar to that documented above, we will not provide a detailed discussion of the corresponding 3D driver code
which computes the deformation of the hollow cube used in the Mesh Generation with Tetgen Tutorial.
Here is an animation of the cube's deformation when it is subjected to gravity, acting in the negative -direction, and a "suction" force applied on the upper face, while its right face is held in a fixed position.
An alternative three-dimensional problem: the deformation of a bifurcating tube is described here.
A pdf version of this document is available.