. 66
( 67 .)


I= , (18.68)
where b is the width of the beam and h the height of the beam and E the Young™s
modulus of the material. The ratio h/L of height over length equal to 0.1 is chosen.
Using the meshes depicted in Fig. 18.3, the results of Table 18.1 are obtained,
which presents the ratio of uL and the computed displacement uc at x = L.
It is clear that the displacement of the beam, using a mesh of linear trian-
gles is much too small. The poor performance of the linear triangle can easily
be understood; because of the linear interpolation of the displacement ¬eld u, the
associated strains computed from
Infinitesimal strain elasticity problems

Figure 18.3
Triangular and quadrilateral element distribution for the bar problem.

µ = Bu (18.69)

are constant per element. The bi-linear quadrilateral element, on the other hand,
is clearly enhanced. A typical shape function, of for example the ¬rst node in an
element, is given by (with respect to the local coordinate system)
N1 = (1 ’ ξ ) (1 ’ ·) . (18.70)
N1 = (1 ’ ξ ’ · + ξ ·) , (18.71)
which means that an additional non-linear term is present in the shape functions.
Therefore a linear variation of the stress ¬eld within an element is represented.
Two remarks have to be made at this point:
• The numerical analysis in this example was based on plane stress theory, while in the
chapter the equations were elaborated for a plane strain problem. How this elaboration
is done for a plane stress problem is discussed in Exercise 18.1.
• Strictly speaking, a concentrated force in one node of the mesh is not correct. Mesh
re¬nement in this case would lead to an in¬nite displacement of the node where the
force is acting. This can be avoided by applying a distributed load over a small part of
the beam.


18.1 For Hooke™s law the Cauchy stress tensor σ is related to the in¬nitesimal
strain tensor µ via

σ = Ktr( µ) I + 2Gµ d .

(a) What is tr(µ) for the plane strain case?
(b) What is tr(µ) for the plane stress case?
325 Exercises

What are the non-zero components of µ and σ for the plane stress
What are the non-zero components of µ and σ for the plane strain
(e) Consider the plane stress case. Let
⎡ ¤ ⎡ ¤
σxx µxx
⎢ ⎥ ⎢ ⎥
σ = ⎣ σyy ¦ ∼ = ⎣ µyy ¦.

σxy µxy
What is the related H matrix for this case, using
σ = H µ?
∼ ∼

(f) For the plane strain case the matrix H is given by
⎡ ¤ ⎡ ¤
4 ’2
110 0
⎢ ⎥ G⎢ ⎥
H = K ⎣ 1 1 0 ¦ + ⎣ ’2 4 0 ¦.
000 0 0 3
Rewrite this matrix in terms of Young™s modulus E and Poisson™s ratio
ν. Is the resulting matrix linearly dependent on E?
18.2 Consider the mesh given in the ¬gure below. The mesh consists of two
linear triangular elements and a linear elasticity formulation applies to this
element con¬guration. The solution u is given by

uT = [ u1 w1 u2 w2 u3 w3 u4 w4 ] ,

where u and v denote the displacements in the x- and y-direction, respec-
(a) What is a possible top array of this element con¬guration assuming
equal material and type identi¬ers for both elements?
(b) What is the pos array for this element con¬guration?

4 3


1 2

(c) What is the dest array?
(d) Based on the pos array the non-zero entries of the stiffness matrix
can be identi¬ed. Visualize the non-zero entries of the stiffness
Infinitesimal strain elasticity problems

Suppose that the boundary nodes in the array usercurves are
stored as
1 2];
The solution u is stored in the array sol. How are the displace-

ments in the x-direction extracted from the sol array along the ¬rst
(f) Let the solution array sol and the global stiffness matrix q be given.
Suppose that both displacements at the nodes 1 and 2 are suppressed.
Compute the reaction forces in these nodes. Compute the total reac-
tion force in the y-direction along the boundary containing the nodes
1 and 2.
18.3 Consider the bi-linear element of the ¬gure below. It covers exactly the
domain ’1 ¤ x ¤ 1 and ’1 ¤ y ¤ 1. Assume that the plane strain
condition applies.
(a) Compute the strain displacement matrix B for this element in point
( x, y) = ( 1 , 1 ).
(b) Let the nodal displacements for this element be given by

uT = .
0 0 0 0 1 0 1 0

What are the strains in ( x, y) = ( 1 , 1 )?

(“1, 1) (1, 1)

4 3


1 2
(“1, “1) (1, “1)

If G = 1 and K = 2 what are the stress components σxx , σyy and σxy
at ( x, y) = ( 0, 0)? What is the value of σzz ?
18.4 Given the solution vector sol and the pos array,
(a) Extract the solution vector for a given element ielem.
(b) Compute the strain at an integration point.
(c) Compute the stress at an integration point.
327 Exercises

18.5 The m-¬le demo_bend in the directory twode analyses the pure bending
of a single element (for quadrilaterals) or two elements (for triangles). The
analysis is based on plane stress theory. The geometry is a simple square
domain of dimensions 1 — 1. Along the left edge the displacements in the
x-direction are set to zero, while at the lower left corner the displacement in
the y-direction is set to zero to prevent rigid body motions. The two nodes
at the right edge are loaded with a force F of opposite sign, to represent
a pure bending moment. Investigate the stress ¬eld for various elements:
linear triangle, bi-linear quadrilateral, and their quadratic equivalents.
Explain the observed differences. To make these choices use itype and
itype = 1 : quadrilateral element
itype = 20 : triangular element
norder = 1 :(bi-)linear element
norder = 2 :(bi-)quadratic element
18.6 In a shearing test a rectangular piece of material is clamped between a
top and bottom plate, as schematically represented in the ¬gure. This
experiment is generally set-up to represent the so-called ˜simple-shear™
con¬guration. In the simple-shear con¬guration the strain tensor is given as

µ = µxy ex ey + µxy ey ex

using the symmetry of the strain tensor. As a consequence, the stress-strain
relation according to Hooke™s law reduces to

σxy = 2Gµxy .

Hence, measuring the clamp forces and the shear displacement provides


a direct means to identify the shear modulus G. However, the ˜simple-
shear™ state is dif¬cult to realize experimentally, since the con¬guration
of the ¬gure does not exactly represent the simple-shear case. This may be
analysed using the m-¬le demo_shear.
(a) Analyse the shear and the simple-shear case using this m-¬le. What
is the difference in boundary conditions for these two cases?
Infinitesimal strain elasticity problems

(b) Why is the simple_shear=0 case not equal to the exact simple-
shear case?
What ratio /h is required to measure the modulus G within 10 %
accuracy? How many elements did you use to obtain this result?
Explain the way the m-¬le computes the modulus G.

[1] Adams, R. A. (2003) Calculus: A Complete Course (Addison, Wesley, Longman).
[2] Brekelmans, W. A. M., Poort, H. W. and Slooff, T. J. J. H. (1972) A new method to
analyse the mechanical behaviour of skeletal parts. Acta Orthop. Scandnav. 43,
[3] Brooks, A. N. and Hughes, T. J. R. (1990) Streamline upwind/Petrov-Galerkin
formulations for convection dominated ¬‚ows with particular emphasis on the
incompressible Navier-Stokes equations. Computer Methods in Applied Mechanics
and Engineering Archive “ Special edition, 199“259.
[4] Carslaw, H. S. and Jaeger, J. C. (1980) Conduction of Heat in Solids (Clarendon
[5] Cacciola, G. R. C. (1998) Design, simulation and manufacturing of ¬bre reinforced
polymer heart valves. Ph.D. thesis, Eindhoven University of Technology.
[6] Frijns, A. J. H., Huyghe, J. M. R. J. and Janssen, J. D. (1997) A validation of the
quadriphasic mixture theory for intervertebral disc tissue. Int. J. Engng. Sci.
35(15), 1419“29.
[7] Fung, Y. C. (1990) Biomechanics: Motion, Flow, Stress, and Growth
[8] Fung, Y. C. (1993) Biomechanics: Mechanical Properties of Living Tissues, 2nd
edition (Springer-Verlag).
[9] Hill, A. V. (1938) The heat of shortening and the dynamic constants in muscle.
Proc. Roy. Soc. London 126, 136“65.
[10] Hughes, T. J. R. (1987) The Finite Element Method (Prentice Hall).
[11] Huyghe, J. M. R. J., Arts, T. and Campen, D. H. van (1992) Porous medium
¬nite element model of the beating left ventricle. Am.J.Physiol. 262
[12] Huxley, A. F. (1957). Muscle structure and theory of contraction. Prog. Biochem.
Biophysic. Chem., 255“318.
[13] Mow, V. C., Kuei, S. C. and Lai, W. M. (1980) Biphasic creep and stress
relaxation of articular cartilage in compression. J. Biomech. Engng. 102 73“84.
[14] Oomens, C. W. J. (1985) A mixture approach to the mechanics of skin and
subcutis. Ph.D. thesis, Twente University of Technology.

[15] Oomens, C. W. J., Maenhout, M., Oijen, C. H. van, Drost, M. R. and
Baaijens, F. P. T. (2003) Finite element modelling of contracting skeletal muscle.
Phil. Trans. R. Soc. London B 358, 1453“60.
[16] Sengers, B. G., Oomens, C. W. J., Donkelaar, C. C. van, and Baaijens, F. P. T.
(2005) A computational study of culture conditions and nutrient supply in cartilage
tissue engineering. Biotechnology Progress 21, 1252“1261.
[17] Sengers, B. G. (2005) Modeling the development of tissue engineered cartilage.


. 66
( 67 .)