# Buckling of a ﬂush-mounted plate in simple shear ﬂow

Μηχανική

18 Ιουλ 2012 (πριν από 6 χρόνια και 7 μέρες)

1.015 εμφανίσεις

Arch Appl Mech (2006) 76:549–566
DOI 10.1007/s00419-006-0061-5
ORIGINAL
H.Luo ∙ C.Pozrikidis
Buckling of a ﬂush-mounted plate in simple shear ﬂow
Received:15 February 2006/Accepted:28 June 2006/Published online:12 August 2006
Abstract The buckling of an elastic plate with arbitrary shape ﬂush-mounted on a rigid wall and deforming
under the action of a uniformtangential load due to an overpassing simple shear ﬂow is considered.Working
under the auspices of the theory of elastic instability of plates governed by the linear von Kármán equation,an
eigenvalue problemis formulated for the buckled state resulting in a fourth-order partial differential equation
with position-dependent coefﬁcients parameterized by the Poisson ratio.The governing equation also describes
the deformation of a plate clamped around the edges on a vertical wall and buckling under the action of its own
weight.Solutions are computed analytically for a circular plate by applying a Fourier series expansion to derive
an inﬁnite system of coupled ordinary differential equations and then implementing orthogonal collocation,
and numerically for elliptical and rectangular plates by using a ﬁnite-element method.The eigenvalues of
the resulting generalized algebraic eigenvalue problemare bifurcation points in the solution space,physically
representing critical thresholds of the uniformtangential load above which the plate buckles and wrinkles due
to the partially compressive developing stresses.The associated eigenfunctions representing possible modes
of deformation are illustrated,and the effect of the Poisson ratio and plate shape is discussed.
Keywords Plate buckling · Membrane wrinkling · Elastic instability · Linear von Kármán equation ·
Finite-element method · Generalized eigenvalue problem
1 Introduction
Buckling and wrinkling of biological and fabricated membranes under the action of a compressive load arises
in several physical contexts including shear ﬂowpast liquid capsules enclosed by polymerized membranes [1,
2],instability of compressed hard ﬁlms on soft substrates [3],and undulation instability of lamellar phases
consisting of multiple molecular sheets in channel ﬂow [4].Although scaling laws describing the topology
of the wrinkled shapes have been derived and classical results from the theory of plates and shells have been
invoked to estimate critical buckling thresholds,a detailed theoretical analysis of the conditions under which
buckling and wrinkling occur under speciﬁc circumstances is not available.
In this paper,we consider the buckling of a membrane patch modeled as an elastic plate deforming under the
action of a simple shear ﬂowapplied on the upper side.The termmembrane refers to a biological object rather
than to a thin shell with inﬁnitesimal bending stiffness,as is traditionally deﬁned in mechanics.Accordingly,
we shall refer to the membrane patch interchangeably as a plate.The shear ﬂow imparts to the membrane a
uniform shear stress that can be distributed over the cross section to yield a uniform body force tangential to
the undeformed shape.Physiological situations where this occurs include instances where a biological cell or
vesicle is captured on a surface [5] or an endothelial cell is subjected to capillary blood ﬂow[6].Fung and Liu
H.Luo · C.Pozrikidis (
B
)
Department of Mechanical and Aerospace Engineering University of California,San Diego,La Jolla,CA 92093-0411,USA
E-mail:cpozrikidis@ucsd.edu
550 H.Luo,C.Pozrikidis
[7] discuss the mechanics of the endotheliumand suggest that the main effect of an overpassing shear ﬂow is
to generate tensions over the exposed part of the cell membrane,while the cell interior is virtually unstressed.
In an idealized depiction,the exposed membrane is a thin elastic patch anchored around its edges on the
endothelium wall and connected to the basal lamina by side walls.In an alternative physical interpretation,
the uniform body force may be attributed to the weight of a uniform vertical plate clamped on a rigid wall.
The in-plane deformation of the membrane generates in-plane tensions whose precise distribution depends
on the material properties,shape of the anchoring rim,and assumed boundary conditions around the rim.
Elementary mechanics indicates that the upstreampart of the plate is stretched,while the downstreamportion
is compressed.Compression raises the possibility of buckling and wrinkling when the uniform load exceeds
critical thresholds representing bifurcation points in the solution space.The computation of these bifurcation
points and associated eigenfunctions is the main objective of our work.
In Sect.2,the relevant eigenvalue problemis formulated based on the linear von Kármán equation for elas-
tic plates.Though solutions of similar eigenvalue problems are available for plates and shells with rectangular
and circular shapes,with a fewexceptions,previous studies have addressed situations where the second deriv-
atives on the right-hand side of the von Kármán equation are multiplied by constant coefﬁcients representing
spatially uniform in-plane stresses,and this considerably simpliﬁes the analysis.In Sects.3 and 4,analytical
and numerical solutions are presented for plates with circular,elliptical,and rectangular shapes,and the effect
of the plate shape on the threshold levels for instability is discussed.
2 Theoretical model
We consider an elastic membrane ﬂush-mounted on a rigid wall with the edge clamped around the rim,as
illustrated in Fig.1.The upper surface of the membrane is exposed to an overpassing shear ﬂow along the
x-axis with velocity u
x
= Gz,where G is the shear rate,and the z-axis is normal to the wall,as shown in Fig.1.
The lower surface of the membrane is in contact with a stationary ﬂuid medium that is unable to withstand
shear stress.The shear ﬂowimparts to the upper surface of the membrane a uniformhydrodynamic shear stress,
τ = μG,where μ is the ﬂuid viscosity.In the context of thin-shell theory for a zero-thickness membrane,
the shear stress can be smeared from the upper surface into the whole cross section of the membrane.When
this is done,the shear stress effectively amounts to an in-plane body force uniformly distributed over the cross
section with components
b
x
=
τ
h
=
μG
h
,b
y
= 0,(1)
where h is the membrane thickness.Justiﬁcation for smearing the shear stress is provided in the Appendix
where the numerical solution of an analogous two-dimensional problemobtained by the ﬁnite-element method
is discussed.The emergingproblemalsodescribes the bucklingof a homogeneous clampedplate ﬂushmounted
on a vertical wall and deforming under the action of its own weight.
We shall assume that the in-plane stresses developing due to the in-plane deformation in the absence of
buckling,σ
i j
,are related to the in-plane strains,
i j
,by the linear constitutive equation

σ
xx
σ
yy
σ
xy

=
E
1 −ν
2

1 ν 0
ν 1 0
0 0 1 −ν

·


xx

yy

xy

,(2)
x
y
z
θ
u = Gz
x
Membrane
Fig.1 Shear ﬂow past a membrane patch modeled as an elastic plate ﬂush-mounted on a plane wall
Buckling of a ﬂush-mounted plate in simple shear ﬂow 551
where

kl
=
1
2

∂v
k
∂x
l
+
∂v
l
∂x
k

,(3)
(v
x
,v
y
) is the tangential displacement of membrane point particles in the x–y plane,E is the membrane
modulus of elasticity,and ν is the Poisson ratio.The upper limit,ν = 0.5,corresponds to an incompressible
material.Although the Poisson ratio is positive for the vast majority of materials encountered in practice,
zero or slightly negative values have been reported for wrinkled membranes [8].Physically,smoothing of the
wrinkles under uniaxial extension leads to expansion in the lateral direction associated with a negative Poisson
Force equilibriumrequires the differential balances
∂σ
xx
∂x
+
∂σ
yx
∂y
+b
x
= 0,
∂σ
xy
∂x
+
∂σ
yy
∂y
+b
y
= 0,(4)
subject to the boundary conditions v
x
= 0 and v
y
= 0 around the clamped rimof the plate.
As a speciﬁc application,we consider an elliptical plate whose axes are aligned in the x and y directions.
The rim is described by the equation (x/a
x
)
2
+ (y/a
y
)
2
= 1,where a
x
and a
y
are the two semi-axes.The
solution of the plane-stress problemcan be found by inspection and is given by
v
x
= V

1 −
x
2
a
2
x

y
2
a
2
y

,v
y
= 0,(5)
where
V =
τ
Eh
(1 −ν
2
) a
2
x
a
2
y
(1 −ν) a
2
x
+2 a
2
y
(6)
is a constant with dimensions of length.The associated in-plane stresses are given by
σ
xx
=
E
1 −ν
2
∂v
x
∂x
= −2
EV
(1 −ν
2
) a
2
x
x = −
τ
h
2 a
2
y
(1 −ν) a
2
x
+2 a
2
y
x,
σ
xy
=
E
2(1 +ν)
∂v
x
∂y
= −
EV
(1 +ν) a
2
y
y = −
τ
h
(1 −ν) a
2
x
(1 −ν) a
2
x
+2 a
2
y
y,(7)
σ
yy
= ν σ
xx
,
independent of the modulus of elasticity,E.For a circular plate of radius a,a
x
= a
y
= a,we obtain the
simpliﬁed expressions
v
x
=
τ
Eh
1 −ν
2
3 −ν
(a
2
−x
2
− y
2
),v
y
= 0,(8)
and associated stresses
σ
xx
= −
2
3 −ν
τ
h
x,σ
xy
= −
1 −ν
3 −ν
τ
h
y,σ
yy
= ν σ
xx
.(9)
These expressions conﬁrmthat the stream-wise component of the in-plane normal stress,σ
xx
,is positive (ten-
sile) on the upstream half,and negative (compressive) on the downstream half of the plate.The transverse
component of the normal stress,σ
yy
,is also positive or negative depending on the sign of the Poisson ratio.
Compression raises the possibility of buckling and wrinkling when the shear stress,τ,exceeds a critical
threshold.To compute the transverse deﬂection along the z axis upon inception of buckling,z = f (x,y),we
552 H.Luo,C.Pozrikidis
work under the auspices of linear elastic stability of thin plates and shells and derive the linear von Kármán
equation,

4
f ≡ ∇
2

2
f =

4
f
∂x
4
+2

4
f
∂x
2
∂y
2
+

4
f
∂y
4
=
h
E
B

σ
xx

2
f
∂x
2
+2 σ
xy

2
f
∂x∂y

yy

2
f
∂y
2
−b
x
∂ f
∂x
−b
y
∂ f
∂y

,(10)
where E
B
is the bending modulus (e.g.,[11] p.18;[12],p.305).
This is a fourth-order differential equation with position-dependent coefﬁcients multiplying the second
derivatives on the right-hand side.Since the membrane is clamped around the rim,the deﬂection satisﬁes the
homogeneous Dirichlet and Neumann boundary conditions f = 0 and ∂ f/∂n = 0,where ∂/∂n denotes the
normal derivative around the rimin the x–y plane.
Substituting the expressions for the in-plane shear stresses in (10),and nondimensionalizing lengths by a
characteristic membrane surface length,a,we derive the dimensionless parameter ˆτ = (τa
3
)/E
B
,expressing
the strength of the shear ﬂow relative to the developing bending moments.Equation (11) admits the trivial
solution,f = 0,for any value of ˆτ,and nontrivial eigensolutions at a sequence of discrete eigenvalues.The
computation of these eigenvalues and corresponding eigenfunctions for a speciﬁed membrane shape is the
main objective of our analysis.
3 Fourier series solution for a circular plate
We begin by considering the buckling of a circular plate of radius a.Substituting expressions (9) for the
in-plane stresses and expression (1) for the body force in (10),we obtain

4
f = −
α
a
3

x

2
f
∂x
2
+(1 −ν) y

2
f
∂x∂y
+ν x

2
f
∂y
2
+
3 −ν
2
∂ f
∂x

,(11)
where
α ≡
2τa
3
E
B
(3 −ν)
(12)
is a dimensionless parameter.The solution is to be found subject to the homogeneous Dirichlet and Neumann
boundary conditions,f = 0 and ∂ f/∂r = 0 at r = a.
An eigensolution of (11) can be expressed as a Fourier series with respect to the plane polar angle θ,deﬁned
such that x =r cos θ and y =r sin θ,as shown in Fig.1,in the form
f (r,θ) =
1
2
p
0
(r) +

n=1

p
n
(r) cos nθ +q
n
(r) sin nθ

=

n=−∞
F
n
(r) exp(−inθ),(13)
where i is the imaginary unit,p
n
(r),q
n
(r) are real functions,
F
n
(r) ≡
1
2

p
n
(r) +i q
n
(r)

(14)
for n ≥ 0,F
n
(r) = F

−n
(r),and an asterisk denotes the complex conjugate.A straightforward computation
yields

4
f =

n=−∞

n
(r) exp(−inθ),(15)
where

n
(r) ≡ F

n
+
2
r
F

n

1 +2n
2
r
2
F

n
+
1 +2n
2
r
3
F

n
+n
2
n
2
−4
r
4
F
n
,(16)
Buckling of a ﬂush-mounted plate in simple shear ﬂow 553
and a prime denotes a derivative with respect to r.To express the right-hand side of (11) in a Fourier series as
well,we use the Cartesian-to-plane-polar transformation rules
∂ f
∂x
= cos θ
∂ f
∂r

sin θ
r
∂ f
∂θ
,

2
f
∂x
2
= cos
2
θ

2
f
∂r
2

sin 2θ
r

2
f
∂r∂θ
+
sin
2
θ
r
∂ f
∂r
+
sin 2θ
r
2
∂ f
∂θ
+
sin
2
θ
r
2

2
f
∂θ
2
,
(17)

2
f
∂y
2
= sin
2
θ

2
f
∂r
2
+
sin 2θ
r

2
f
∂r∂θ
+
cos
2
θ
r
∂ f
∂r

sin 2θ
r
2
∂ f
∂θ
+
cos
2
θ
r
2

2
f
∂θ
2
,

2
f
∂x∂y
=
1
2

sin 2θ

2
f
∂r
2
+2
cos 2θ
r

2
f
∂r∂θ

sin 2θ
r
∂ f
∂r
−2
cos 2θ
r
2
∂ f
∂θ

sin 2θ
r
2

2
f
∂θ
2

.
Using Euler’s formula to express the cosines and sines in terms of the complex exponential exp(−iθ),and
substituting the Fourier expansion,we ﬁnd
∂ f
∂x
=
1
2

n=−∞


F

n
+n
F
r

e

+

F

n
−n
F
r

e
−iθ

exp(−inθ),

2
f
∂x
2
=
1
2

n=−∞

H
n
+ R
n
e
2iθ
+T
n
e
−2iθ

exp(−inθ),
(18)

2
f
∂y
2
=
1
2

n=−∞

H
n
− R
n
e
2iθ
−T
n
e
−2iθ

exp(−inθ),

2
f
∂x∂y
= −
i
2

n=−∞

R
n
e
2iθ
−T
n
e
−2iθ

exp(−inθ),
where
H
n
≡ F

n
+
F

n
r
−n
2
F
n
r
2
,R
n

1
2

F

n
+(2n −1)
F

n
r
+n(n −2)
F
n
r
2

,
(19)
T
n

1
2

F

n
−(2n +1)
F

n
r
+n(n +2)
F
n
r
2

.
If the function F
n
is real,in which case F
n
= F
−n
,the functions R
n
and T
n
are related by R
n
= T
−n
.
Using the preceding expressions,we ﬁnd that the term enclosed by the square brackets on the right-hand
side of (11) takes the simple form
1
2

n=−∞


n
e

+
n
e
−iθ

exp(−inθ) =
1
2

n=−∞


n+1
+
n−1

exp(−inθ),(20)
where

n
=r F

n
+

3 +ν
2
+(1 −ν)n

F

n
+n

1 +ν
2
−νn

F
n
r
,
(21)

n
=r F

n
+

3 +ν
2
−(1 −ν)n

F

n
−n

1 +ν
2
+νn

F
n
r
.
Substituting (20) and (15) in (11),and equating corresponding Fourier coefﬁcients,we derive an inﬁnite
tridiagonal systemof ordinary differential equations,

n
= −
α
2a
3
(

n+1
+
n−1
),(22)
for n = 0,±1,±2,....Approximate eigenvalues can be computed by truncating the systemat a ﬁnite level,
n = ±N,and solving the eigenvalue problem deﬁned by the retained ordinary differential equations.In the
554 H.Luo,C.Pozrikidis
case of eigensolutions with left-to-right symmetry with respect to the z–x plane,the Fourier series onlyinvolves
cosine terms,the component functions F
n
are real,F
n
= F
−n
,and

−n
=
n
.The general system(22) then
reduces to

0
= −
α
a
3

1
,
n
= −
α
2a
3
(

n+1
+
n−1
),(23)
for n = 1,2,...,N.On the other hand,in the case of eigensolutions that are antisymmetric with respect to the
z–x plane,the Fourier series only involves sine terms,the component functions F
n
are imaginary,F
n
= −F
−n
,

−n
= −
n
,and
0
= 0.
To satisfy the boundary conditions F = F

= 0 around the rim,we set
F
i
= (a −r)
2
H
i
(r),(24)
and approximate the modulating functions H
i
(r) with Mth-degree polynomials in r,
H
i
(r) = c
i 0
+c
i 1
r +· · · +c
i M
r
M
,(25)
where M is a speciﬁed polynomial order.Substituting (24) and (25) into the governing differential equations,
computing the derivatives analytically by differentiating the polynomials,and enforcing the resulting equations
at the Chebyshev collocation points
r
j
= cos


j −
1
2

π
M +1

,(26)
for j = 1,2,...,M +1,we obtain a generalized eigenvalue problemof the form
A· c = α B· c,(27)
where the vector c contains the polynomial coefﬁcients c
i j
.The bandedmatrices AandBarise fromthe orthogo-
nal collocation.Speciﬁcally,thematrixAconsists of 2N+1diagonal blocks withdimensions (M+1)×(M+1),
and the matrix B consists of 2N super-diagonal and sub-diagonal blocks with same dimensions.
The numerical task is to compute as many eigenvalues as possible,beginning with the smallest eigenvalue
and moving upward.Physically,the smallest eigenvalue represents the minimum shear stress for buckling.
Numerical experimentation showed that the best results are obtained by ﬁrst inverting the matrix A,and then
solving the standard eigenvalue problem,D· c = λc,for λ ≡ 1/α,concentrating on the largest eigenvalue,λ,
and moving downward,where D ≡ A
−1
· B.Since the matrix A arises from the discretization of the elliptic
biharmonic operator,it is nonsingular and well conditioned.By contrast,the matrix B is poorly conditioned.
The computation of the eigenvalues was carried out using the Matlab eig function.
Eigenvalues come in degenerate positive–negative pairs representing identical buckling states that arise
when the shear ﬂow is oriented toward the positive or negative direction of the x-axis.All the corresponding
eigenvectors,c,turn out to be real,which means the eigenfunctions are either symmetric or antisymmetric
with respect to the z–x plane,designated by S or A,respectively.Symmetry or antisymmetry can be detected
by comparing the signs of the coefﬁcients of the polynomials H
−n
and H
n
,and was veriﬁed by solving (23).
Multiple eigenvalues with distinct eigensolutions arise only at speciﬁc values of the Poisson ratio,as will be
discussed in the next section.
By way of illustrating the performance of the numerical method and demonstrating the convergence of
the numerical results,in Table 1a we list numerical values for the smallest eigenvalue for ν = 0.5,against the
Fourier truncation level,N,and polynomial degree,M.The corresponding mode of deformation is designated
as symmetric mode,S1,as shown in Fig.2.The results in Table1a suggest that accuracy up to the ﬁfth decimal
point for this mode can be achieved with a moderate Fourier truncation level,N = 6,and polynomial degree
truncation level,M = 25.Table1b illustrates the convergence of the sixth smallest eigenvalue corresponding
to the symmetric mode,S4,and Table1c illustrates the convergence of the twelfth smallest eigenvalue corre-
sponding to the symmetric mode,S7.The associated eigenfunctions are shown in Fig.2.Larger eigenvalues
require higher truncation levels,resulting in matrices of large size whose eigenvalues can be identiﬁed with
sufﬁcient accuracy only up to a point.In the next section,results will be presented for those eigenvalues that
could be extracted with conﬁdence to shown accuracy.
Buckling of a ﬂush-mounted plate in simple shear ﬂow 555
Table 1 a Convergence of the smallest eigenvalue,α = 71.5795 (S1 mode) for ν = 0.5,listed with respect to the Fourier
truncation level,N,and expansion order,M.b Convergence of the eigenvalue α = 275.18 (S4 mode) for ν = 0.5.c Convergence
of the eigenvalue α = 452 (S7 mode) for ν = 0.5
M (a) N = 1 N = 2 N = 3 N = 4 N = 5 N = 6
14 78.77035 95.12050 58.18914 80.58670 79.79018 79.77308
17 78.77057 73.47277 70.94721 70.82917 70.81285 70.81218
18"72.05807 71.61643 71.60297 71.60484 71.60510
21"72.07656 71.60049 71.57975 71.57926 71.57926
25"72.07594 71.60076 ” 71.57951 ”
M (b) N = 5 N = 6 N = 7 N = 8 N = 9 N = 10
22 286.02569 282.62062 280.65468 280.08687 280.00046 280.02944
25 282.33248 277.94532 275.83596 275.29922 275.22545 275.21879
30"277.91176 275.80062 275.26349 275.18968 275.18355
M (c) N = 5 N = 6 N = 7 N = 8 N = 9 N = 10
22 405.02099 517.40551 483.33009 476.57590 461.19203 449.06499
25 446.92262 436.43061 469.63320 457.21055 453.31981 451.89132
30 450.28448 442.29069 464.69797 454.91600 452.68277 451.94615
3.1 Results and discussion
Table2 displays the computed eigenvalues for three Poisson ratios,listed in order of increasing magnitude.
Concentrating on the ﬁrst column for ν = 0.5 corresponding to an incompressible membrane,we note the
interlaced appearance of symmetric (S) and antisymmetric (A) modes.A left-to-right break of symmetry in
the solution space occurs for the antisymmetric modes.The component functions,F
n
,and corresponding
eigenfunctions are plotted in Fig.2.All functions F
n
except for F
0
,tend to zero at the origin,r = 0,which
is necessary for the deﬂected shape to be single-valued.A truly incompressible membrane will not be able
to exhibit transverse deﬂection without reducing the area of the base enclosed by the clamped circular edge.
However,the change in area due to the deﬂection is of second order with respect to the maximum deﬂected
height,and the results for ν = 0.5 are relevant to a nearly incompressible material.
To interpret the computed eigenvalues from a physiological standpoint,we reverse the deﬁnition of α
and ﬁnd that the critical shear rates where buckling occurs are given by G = α E
B
(3 −ν)/(2μa
3
).Taking
E
B
 1 ×10
−12
dyncm,which is typical of a biological membrane,α = 5μm,and μ= 1.2cp = 1.2mPas,we
ﬁnd G = α (1 −ν/3) s
−1
.For ν = 0.5,buckling in the S1 mode will occur when G = 71.6 (1 −1/6) s
−1

60.0s
−1
.In human circulation,the shear stress varies in the range 1–2Pa through all branches,corresponding to
G ∼ 100s
−1
.This estimate exceeds the S1 buckling threshold but not subsequent thresholds.In ﬂowthrough
a cylindrical capillary of radius b,the wall shear rate is related to the mean velocity,U
m
,by G = 4U
m
/b.
Taking b = 10 μm,we ﬁnd that,as U
m
is raised from 0.001 to 1cm/s,the wall shear rate increase from 4 to
4,000s
−1
,which includes several buckling modes.
To establish a further point of reference for the numerical results,we consider a circular membrane that
is compressed around the edges with a uniform force per unit length,P
r
.Classical analysis shows that the
symmetric S1 buckling mode occurs at the critical eigenvalue β
S1
≡ P
r
a
2
/E
B
 14.68 ([12],p.368) which
is nearly ﬁve times smaller than the present eigenvalue for the corresponding mode,α
S1
= 71.58.This large
difference is expected in light of the partial compression of the membrane patch in shear ﬂow presently con-
sidered.The uniformly compressed plate also exhibits antisymmetric buckling modes associated with a break
of radial symmetry.The ﬁrst antisymmetric mode,A1,occurs at the critical eigenvalue β
A1
 28.87 ([13],see
e.g.,[11],p.230) which is also nearly ﬁve times less than the present eigenvalue for the corresponding mode.
The results in Table2 show that the Poisson ratio can have an important effect on the spectrum of eigen-
values.The most striking effect is that,as ν is reduced from0.5 to 0.25,the relative position of symmetric and
antisymmetric eigenmodes is altered.For example,when ν = 0.25,the eigenvalue of the S2 mode is lower
than that of the A1 mode,and is thus expected to appear ﬁrst in an experiment where the shear rate is gradually
increased fromzero.Because of mode crossing,there is a critical Poisson ratio where the eigenvalues of the A1
and S2 modes are identical.Our computations showthat this occurs when ν = 0.42703 and α = 138.308.The
eigenfunctions of this double eigenvalue are arbitrary superpositions of symmetric and antisymmetric modes,
and may thus have an arbitrary orientation is space.As the Poisson ratio is further decreased,the ordering of
556 H.Luo,C.Pozrikidis
Fig.2 Fourier modes and buckled shapes for ν = 0.5,corresponding to the eigenvalue:a α = 71.5795 (S1 mode),bα = 132.128
(A1 mode),c α = 137.29 (S2 mode),dα = 208.409 (S3 mode).e α=212.107 (A2 mode),f α = 275.18 (S4 mode),g α = 300.13
(A3 mode),h α = 318.1 (S5 mode),i α=364.5 (A4 mode),j α = 403 (S6 mode),k α = 429 (A5 mode),l α = 452 (S7 mode)
the eigenvalues changes even further.Families of eigenfunctions for Poisson ratios ν = 0.25 and 0 are shown
in Figs.3 and 4.Comparing these shapes with those shown in Fig.2 for ν = 0.5,we ﬁnd that the low modes,
such as S1,S2,A1,and A2,are virtually identical.The Poisson ratio affects the shape of higher modes by
changing elastic properties of the membrane.
Because in the mathematical model the upper and lower surface of the membrane are indistinguishable,
the eigenfunctions can be ﬂipped with respect to the x–y plane in order to facilitate the comparison.In the
context of the zero-thickness membrane model,the transverse deformation is indeterminate.In reality,the
inner surface of a cell membrane is supported by the cytoskeleton and the outer surface of a membrane in
Buckling of a ﬂush-mounted plate in simple shear ﬂow 557
Fig.2 (contd.)
the endotheliumis coated with the glycocalyx and exposed to the shear ﬂow.Numerical solutions of a model
problemdiscussed in the Appendix suggest that the indeterminacy is removed for a ﬁnite-thickness membrane.
4 Finite-element solution for arbitrary membrane shapes
Adual ﬁnite-element method was implemented for describing the buckling of a non-circular plate.The numeri-
cal procedure involves two steps:determination of the in-plane membrane tensions in the absence of transverse
deﬂection,and solution of the eigenvalue problem based on the linear von Kármán equation.For elliptical
membranes aligned with the ﬂow,the analytical solution presented in Sect.2 was used to generated the in-plane
558 H.Luo,C.Pozrikidis
Fig.2 (contd.)
tensions.For more general shapes,the tensions were found by performing a plane-stress analysis using the
Galerkin ﬁnite-element method with six-node quadratic descendant elements (e.g.,[14]).
Grid generation was carried out by successively subdividing a parental element structure into smaller ele-
ments.The results of the plane-stress code were conﬁrmed to reproduce available analytical solutions.An
alternative of this dual implementation would be a uniﬁed approach wherein the eigenvalue problemis solved
at one stage in the context of shell buckling.The present formulation reduces the number of unknowns at the
expense of solving a fourth-order differential equation.
Since the von Kármán equation involves the fourth-order biharmonic operator on the left-hand side,
Hermitian conforming elements must be used to ensure continuity of the solution and its gradient across
Buckling of a ﬂush-mounted plate in simple shear ﬂow 559
Table 2 Eigenvalues for a circular membrane at three Poissonratios;Sdenotes a symmetric mode,andAdenotes anantisymmetric
mode
Mode ν = 0.5 ν = 0.25 ν = 0
1 71.5795000 (S1) 77.4095411 (S1) 83.4655475 (S1)
2 132.128774 (A1) 140.1698392 (S2) 141.6740559 (S2)
3 137.294138 (S2) 154.5573454 (A1) 178.89886 (A1)
4 208.409 (S3) 225.54982 (A2) 233.834 (A2)
5 212.107 (A2) 242.6051 (S3) 263.6254 (S3)
6 275.18 (S4) 295.922 (S4) 320.74 (S4)
7 300.13 (A3) 340.92 (A3) 365.4 (A3)
8 318.1 (S5) 346.67 (S5) 372.55 (S5)
9 364.5 (A4) 413.8 (S4) 442 (A4)
10 403 (S6) 442 (S6)
11 429 (A5) 472 (S7)
12 452 (S7) 487 (A5)
Fig.3 Eigenfunctions for a circular membrane with Poisson ratio ν = 0.25
the element edges,adding to the complexity of the formulation.We have adopted the Hrieh-Clough-Tocher
(HCT) triangular element deﬁned by three vertex nodes,three edge nodes,and one hidden interior node (e.g.,
[14]).In the formulation,each element is subdivided into three subelements,and the solution over each subel-
ement is approximated with a complete cubic polynomial in x and y involving 10 coefﬁcients,for a total of 30
coefﬁcients per element.After imposing constraints the C
1
continuity condition,the HCT element is endowed
with 12 degrees of freedomconsisting of the values of the solution and its Cartesian or directional derivatives.
The ﬁnite-element grid is identical to that employed for solving the plane-stress problem.The input to the
buckling code includes the components of the in-plane stress tensor at the element vertex nodes.
560 H.Luo,C.Pozrikidis
Fig.4 Eigenfunctions for a circular membrane with Poisson ratio ν = 0
To develop the Galerkin ﬁnite-element method for the buckling problem,we multiply (10) by each one of
the global interpolation functions associated with the available degrees of freedom,and integrate the product
over the patch surface,D,to obtain

D
φ
i
(x,y) ∇
4
f dx dy =

D
φ
i
(x,y) w(x,y) dx dy,(28)
where w stands for the right-hand side of (10).Integrating by parts on the left-hand side to reduce the order of
the biharmonic operator,∇
4
f,we write

D
φ
i

4
f dx dy =

D
φ
i
(x,y) ∇ ·

∇(∇
2
f )

dx dy
=

D
∇ ·

φ
i
∇(∇
2
f )

dx dy −

D
∇φ
i
· ∇(∇
2
f ) dx dy
=

C
φ
i
n · ∇(∇
2
f ) dl −

D
∇φ
i
· ∇(∇
2
f ) dx dy,(29)
where C is the boundary of D,n is the unit vector in the x–y plane that is normal to C and points outward
from D,and l is the arc length along C.Integrating by parts once more the last integral in (29),we ﬁnd
Buckling of a ﬂush-mounted plate in simple shear ﬂow 561

D
φ
i

4
f dx dy =

C
φ
i
n · ∇(∇
2
f ) dl


C

2
f (n · ∇φ
i
) dl +

D

2
φ
i

2
f dx dy.(30)
Finally,we substitute this expression in (28),and rearrange to obtain

D

2
φ
i

2
f dx dy = −

C
φ
i
n · ∇(∇
2
f ) dl
+

C

2
f n · ∇φ
i
dl +

D
φ
i
(x,y) w(x,y) dx dy.(31)
The contour integrals on the right-hand side of (31) are non-zero only if φ
i
is a boundary mode associated
with a boundary node.If φ
i
is a mode associated with an interior node,these integrals are zero,yielding the
simpliﬁed equation

D

2
φ
i

2
f dx dy =

D
φ
i
(x,y) w(x,y) dx dy.(32)
Since the rim of the membrane is assumed to be clamped,the transverse displacement and its gradient are
known along the boundary,C,the corresponding projections are excluded fromthe ﬁnite-element formulation,
and the simpliﬁed formulation applies.
Assume that the ﬁnite-element expansion involves N
m
G
modes associated with the available degrees of
freedom,
f (x,y) =
N
m
G

j =1
h
G
j
φ
j
(x,y),(33)
where the coefﬁcients h
G
j
represent either the values of the solution or its spatial derivatives at the unique
global nodes.Inserting this expansion in (32),we derive a generalized eigenvalue problemof the form
K
i j
h
G
j
= ˆτ R
i j
h
G
j
,(34)
involving N
m
G
unknowns,where ˆτ = (τa
3
)/E
B
is the dimensionless membrane tension,and summation over
j is implied on the left-hand side.The global bending stiffness matrix,K
i j
,can be compiled in the usual way
fromthe element bending stiffness matrices,
A
(l)
i j


E
l

2
ψ
i

2
ψ
j
dx dy,(35)
where E
l
denotes the lth element.The matrix R
i j
can be compiled from corresponding element matrices.To
eliminate spurious eigenvalues,the ﬁnal algebraic systemis condensed by eliminating the constrained bound-
ary degrees of freedom.The best results are obtained by transforming the generalized eigenvalue probleminto
a standard eigenvalue problem for the matrix D ≡ K
−1
· R,whose eigenvalues are λ = 1/ˆτ.Since we are
interested in the lowest eigenvalues,ˆτ,representing the threshold in the shear ﬂowwhere buckling occurs,we
concentrate on the highest inverse eigenvalues,λ,which we compute using the Matlab eig function.
The plate buckling code was validated by comparing the results with available solutions.In a ﬁrst test,the
buckling of an entirely clamped square plate with side length a was considered,subject to compression along
all four edges with a uniform force per unit length,p.The in-plane stress ﬁeld admits the simple isotropic
representation σ
xx
= −p/h,σ
xy
= 0,σ
yy
= −p/h,corresponding to the null body force,b
x
= 0,b
y
= 0.
Dimensionless eigenvalues,ˆp ≡ pa
2
/E
B
,computed with two-element discretizations are shown in Table3a.
562 H.Luo,C.Pozrikidis
Table 3 a Dimensionless eigenvalues,ˆp,for a clamped square plate compressed along all four edges;N
E
is the number of
elements,N
F
is the number of condensed degrees of freedom,and an asterisk denotes a double eigenvalue.b Dimensionless
eigenvalues,ˆτ,for a clamped circular membrane in shear ﬂow,and Poisson ratio,ν = 0.5;S indicates a symmetric mode,and A
indicates an antisymmetric mode
Finite-element method Finite-element method Separation of variables [15]
N
E
=128,N
F
=323 (a) N
E
=512,N
F
=1,411
52.824 52.391 52.358
93.437* 92.263* 92.157*
132.157 128.684 128.302
156.788 154.391 154.226
171.083 167.469 167.066
198.552* 190.702* 189.790*
– 247.112* 246.416*
– 270.520 269.503
Mode (b) Finite-element method Finite-element method Fourier series
N
E
=64,N
F
=163 N
E
=256,N
F
=707
S1 96.478 90.794 89.474
A1 182.621 168.501 165.161
S2 185.222 174.222 171.618
S3 290.422 267.233 260.511
A2 271.613 265.133
S4 352.741 343.979
A3 387.715 375.162
The third column shows numerical results obtained by the method of separation of variables [15].The ﬁnite-
element solution with 512 elements and 1,411 degrees of freedomcaptures with remarkable accuracy the ﬁrst
eight eigenvalues up to the second or third signiﬁcant ﬁgure.
In the second test,the buckling of a circular plate of radius a discussed in Sect.3 was considered.Dimen-
sionless eigenvalues,ˆτ,computed using the ﬁnite-element method for two-element discretizations are shown
in Table3b for Poisson ratio,ν = 0.5.The numerical results with the ﬁne discretization comprised of 256
elements and 707 degrees of freedom are in good agreement with the highly accurate results based on the
Fourier series.The difference escalates from1.5%for the ﬁrst eigenvalue to 3%for the seventh eigenvalue.
Results for non-circular shapes are expected to carry a comparable amount of error.
4.1 Results and discussion
Table4 documents the effect of the aspect ratio of an elliptical membrane on the spectrumof eigenvalues for
ν = 0.5 and 0.0.The membrane with aspect ratio 2:1 is oriented in the direction of the ﬂow,and the membrane
with aspect ratio 1:2 is oriented transversely to the ﬂow.In all cases,the reduced critical shear stress is deﬁned
as ˆτ ≡ τa
3
/E
B
;the equivalent membrane radius,a,is related to the membrane surface area by A = πa
2
.The
most signiﬁcant effect of the membrane aspect ratio is an increase in the threshold for the onset of buckling
with respect to the circular shape.A second effect is a change in the order of appearance of symmetric and
antisymmetric modes.In particular,as the membrane becomes more elongated in the direction of the ﬂow,the
symmetric modes aggregate at the lowest thresholds.On the other hand,as the membrane elongates laterally,
antisymmetric modes are favored near the lowest thresholds.Figure5 illustrates the ﬁrst six buckling modes as
they arise fromthe ﬁnite element solution for ν = 0.5.The shapes for ν = 0 are qualitatively similar.Table5
Table 4 Dimensionless eigenvalues,ˆτ,for a clamped elliptical membrane with Poisson ratio,ν = 0.5 (ﬁrst three columns) and
0 (last three columns),computed with N
E
=256 elements and N
F
=707 degrees of freedom
a
x
/a
y
= 2.0 a
x
/a
y
= 1.0 a
x
/a
y
= 0.5 a
x
/a
y
= 2.0 a
x
/a
y
= 1 a
x
/a
y
= 0.5
122.987 (S1) 89.474 (S1) 158.843 (S1) 235.514 (S1) 125.198 (S1) 180.467 (S1)
164.837 (S2) 165.161 (A1) 219.173 (A1) 265.156 (S2) 212.511 (S2) 277.782 (A1)
262.202 (S3) 171.618 (S2) 295.629 (S2) 421.099 (S3) 268.348 (A1) 383.081 (S2)
277.502 (A1) 260.511 (S3) 360.295 (S3) 514.724 (S4) 350.751 (A2) 413.458 (S3)
339.901 (A2) 265.133 (A2) 382.463 (A2) 545.558 (A1) 395.438 (S3) 502.843 (A2)
368.626 (A3) 343.975 (A2) 466.346 (A3) 558.298 (A2) 481.110 (S4) 588.446 (A3)
448.904 (A4) 300.13 (A3) 489.130 (S4) 714.967 (S5) 548.100 (A3) 650.265 (S4)
Buckling of a ﬂush-mounted plate in simple shear ﬂow 563
S1 S2
–1.5
–1
–0.5
0
0.5
1
1.5
1
0.5
0
0.5
1
–0.5
0
0.5
1
y
x
z
–1.5
–1
–0.5
0
0.5
1
1.5
1
0.5
0
0.5
1
–0.3
–0.2
–0.1
0
0.1
0.2
0.3
0.4
y
x
z
S3 A1
–1.5
–1
–0.5
0
0.5
1
1.5
–1
–0.5
0
0.5
1
0
y
x
z
–1.5
–1
–0.5
0
0.5
1
1.5
1
0.5
0
0.5
1
–0.6
–0.4
–0.2
0
0.2
0.4
0.6
y
x
z
A2 A3
–1.5
–1
–0.5
0
0.5
1
1.5
1
0.5
0
0.5
1
–3
–2
–1
0
1
2
3
y
x
z
–1.5
–1
–0.5
0
0.5
1
1.5
1
0.5
0
0.5
1
–0.6
–0.4
–0.2
0
0.2
0.4
0.6
y
x
z
Fig.5 Buckled shapes of an elliptical membrane with aspect ratio 2:1,for Poisson ratio ν = 0.5
summarizes the effect of the membrane aspect ratio on the lowest threshold for instability,corresponding to the
S1 mode,for ν = 0.5.Recalling that the numerical error is on the order of 1%,we can state with conﬁdence
that,given the membrane surface area,the circular shape has the lowest buckling threshold.
All computations presented thus far correspond to the aligned elliptical shape,wherein the in-plane stress
ﬁeld of the undeﬂected shape is a known linear function of x and y.Next,we consider a case where the in-plane
stress ﬁeld is determined fromthe ﬁnite-element solution.Figure6 shows the distribution of the in-plane stress,
σ
xx
,on a square membrane with edge length 2a,together with the ﬁrst ﬁve buckled shapes occurring at the
critical shear stresses,ˆτ = 73.294,129.577,141.174,203.630,and 212.307.These results demonstrate that the
presence of corners does not have a profound effect on the overall features of the buckling modes.
564 H.Luo,C.Pozrikidis
Table 5 Effect of the aspect ratio of an elliptical membrane aligned with the ﬂowon the lowest threshold for buckling instability
for ν = 0.5
Aspect ratio ˆτ Aspect ratio ˆτ
1.0 89.474 1.0 89.474
1.2 89.223 0.9 94.524
1.4 92.850 0.8 101.222
1.6 100.103 0.7 112.246
1.8 110.262 0.6 129.959
2.0 122.987 0.5 158.843
–1
–0.5
0
0.5
1
1
0.5
0
0.5
1
–1
–0.8
–0.6
–0.4
–0.2
0
0.2
0.4
0.6
0.8
1
y
x
–1
–0.8
–0.6
–0.4
–0.2
0
0.2
0.4
0.6
0.8
1
–1
–0.5
0
0.5
1
0
0.2
0.4
0.6
y
x
z
–1
–0.5
0
0.5
1
–1
–0.5
0
0.5
1
–1
–0.8
–0.6
–0.4
–0.2
0
0.2
0.4
0.6
0.8
1
y
x
z
–1
–0.5
0
0.5
1
–1
–0.5
0
0.5
1
–0.5
–0.4
–0.3
–0.2
–0.1
0
0.1
0.2
0.3
0.4
0.5
y
x
z
–1
–0.5
0
0.5
1
–1
–0.5
0
0.5
1
–1.5
–1
–0.5
0
0.5
1
1.5
y
x
z
–1
–0.5
0
0.5
1
–1
–0.5
0
0.5
1
–0.3
–0.2
–0.1
0
0.1
0.2
0.3
y
x
z
Fig.6 Distribution of the in-plane stress,σ
xx
,and the ﬁrst ﬁve buckled modes of a square membrane for Poisson ratio ν = 0.5
Buckling of a ﬂush-mounted plate in simple shear ﬂow 565
κ
τ
Fig.7 Schematic illustration of the solution diagramin the shear stress – center-point curvature solution plane for a ﬂat membrane
(solid lines),and for a curved membrane (broken line)
5 Discussion
The computed eigenvalues represent bifurcation points in a certain solution plane,such as the plane of the
hydrodynamic shear stress,τ,and curvature at the center of a circular or elliptical membrane,κ.Branches of
possible shapes originating from these eigenvalues are schematically drawn with solid lines in Fig.7.In the
linearized analysis presently undertaken,only shapes near the bifurcation points are described to ﬁrst order
with respect to the membrane deﬂection amplitude.To describe substantially deformed shapes we must carry
out a nonlinear post-buckling analysis based on the full von Kármán equation (e.g.,[16]).The graphs in Fig.7
indicate that,for a given shear stress,there may be one,two,or a higher number of buckled states.Which one
will occur depends on the stability of the eigenfunctions,which is presently unknown.
The simplifying assumption of a perfectly ﬂat undeformed membrane patch is responsible for the sudden
appearance of buckled shapes at the lowest critical shear stress.In blood ﬂow through a capillary over an
endothelial cell,because the height of the cell can be as high as 10%the lateral diameter,the cell membrane
possesses a certain amount of curvature in the undeformed conﬁguration [10].
The solution diagramfor a curved membrane originates froma point above the origin of the κ axis,as shown
with the broken line in Fig.7.As the shear stress is raised,the resting shape is continuously deformed and
then undergoes a rapid transition near the bifurcation points.In this light,the results of the idealized analysis
for a ﬂat patch are useful in that they provide estimates for the conditions under which sudden changes in the
membrane shape are most likely to occur.
A second consequence of the ﬂat resting-shape approximation is that the shear stress is constant over
the entire upper surface of the membrane,corresponding to the simple shear ﬂow.A curved shape causes a
disturbance ﬂow with an associated perturbation shear stress that is highest near the elevated center of the
membrane and lowest around the depressed region around the rim.Moreover,when a membrane patch is
naturally curved or buckles in shear ﬂow,the normal hydrodynamic stress varies over its area as in hydro- and
shear stress and pressure are ﬁrst-order effects with respect to the cell height,neglecting them is consistent
with the small-deﬂection buckling analysis based on the linearized von Kármán equation.
Acknowledgments Support for this research has been provided by the National Science Foundation.
Appendix
Deformation of a sheared elastic slab
To support the argument that a shear stress over the upper surface of a thin membrane is tantamount to a
parallel body force distributed over the cross section of the membrane,we consider a two-dimensional model
wherein a rectangular slab of an elastic material is subjected to a shear stress on the upper surface,while the
traction is required to be zero on the lower surface.The deformation is restricted by rigid side walls.The slab
material is assumed to obey the constitutive equations of linear elasticity involving a modulus of elasticity and
the Poisson ratio.
566 H.Luo,C.Pozrikidis
(a) (b)
–1
–0.8
–0.6
–0.4
–0.2
0
0.2
0.4
0.6
0.8
1
–0.15
–0.1
–0.05
0
0.05
0.1
0.15
x
y
1
0.5
0
–0.5
–1
–0.15
–0.1
–0.05
0
0.05
0.1
0.15
–0.06
–0.04
–0.02
0
0.02
0.04
0.06
x
y
Fig.8 Deformation of an elastic slab subjected to a shear stress on the upper surface.a Finite-element nodes before (crosses)
and after (circles connected by lines) deformation.b Corresponding distribution of the horizontal normal stress,σ
xx
Theequations of planestress for themodel problemweresolvedusingastandardﬁnite-element methodwith
six-node quadratic elements.Figure8a shows the ﬁnite-element nodes before (crosses) and after deformation
(circles connected by lines),for a slender slab with aspect ratio equal to 10 and Poisson ratio ν = 0.25.To
avoid the occurrence of corner singularities,the distribution of shear stress on the upper surface is parabolic,
reaching a maximum at the center and dropping to zero at the two ends,σ
xy
∼ 1 − x
2
.The results reveal
that the shear stress causes a sinusoidal deformation with the crest occurring at the left half and the trough
occurring at the right half of the slab.
Figure8b shows the corresponding distribution of the horizontal normal stress,σ
xx
.Even though the aspect
ratio of the slab is only moderately small,σ
xx
exhibits a mild variation across the slab cross section far from
the side walls.If the shear stress on the upper surface were smeared uniformly fromthe upper surface into the
whole cross section of the slab,σ
xx
would be constant over each cross section.The variation of σ
xx
with x
arising fromthe ﬁnite element solution is consistent with the theoretical prediction based on a zero-thickness
membrane model,σ
xx
∼ x −x
3
/3.
References
1.Walter,A.,Rehage,H.,Leonhard,H.:Shear-induceddeformations of microcapsules:shape oscillations andmembrane folding.
Coll Surf A Physicochem.Eng.Asp.183,123–132 (2001)
2.Finken,R.,Seifert,U.:Wrinkling of microcapsules in shear ﬂow.J.Phys.Condens.Matter 18,L185–L191 (2006)
3.Huang,Z.,Hong,W.,Suo,Z.:Evolution of wrinkles in hard ﬁlms on soft substrates.Phys.Rev.E 70,030601(R) (2004)
4.Zilman,A.G.,Granek,R.:Undulation instability of lamellar phases under shear:a mechanismfor onion formation?Eur.Phys.
J.B 11,593–608 (1999)
5.Kamm,R.D.:Cellular ﬂuid mechanics.Annu.Rev.Fluid Mech.34,211–232 (2002)
6.Liu,S.Q.:Biomechanical basis of vascular tissue engineering.Crit.Rev.Biomed.Eng.27,75–148 (1999)
7.Fung,Y.C.,Liu,S.Q.:Elementary mechanics of the endotheliumof blood vessels.J.Biomech.Eng.ASME 115,1–12 (1993)
8.Pieper,G.,Rehage,H.,Barthès-Biesel,D.:Deformation of a capsule in a spinning drop apparatus.J.Coll.Interf.Sci.202,
293–300 (1988)
9.Boal,H.,Seifert,U.,Shillock,J.C.:Negative Poisson ratio in two-dimensional networks under tension.Phys.Rev.E 48,
4274–4283 (1993)
10.Schmid-Schoenbein,G.W.,Kosawada,T.,Skalak,R.,Chien,S.:Membrane model of endothelial-cells and leukocytes – a
proposal for the origin of a cortical stress.J.Biomech.Eng.ASME 117,171–178 (1995)
11.Bloom,F.,Cofﬁn,D.:Handbook of Thin Plate Buckling and Postbuckling.Boca Raton:Chapman &Hall/CRC (2001)
12.Timoshenko,S.P.,Gere,J.M.:Theory of Elastic Stability,2nd edn.New York:McGraw Hill (1961)
13.Mossakowski,J.:Buckling of circular plates with cylindrical orthotropy.Arch.Mech.Stos.(Poland) 12,583–596 (1960)
14.Pozrikidis,C.:Introduction to Finite and Spectral Element Methods using Matlab.Boca Raton:Chapman &Hall/CRC(2005)
15.Muradova,A.D.:Numerical techniques for linear and nonlinear eigenvalue problems in the theory of elasticity.In:Proceed-
ings of the 12th Computational Techniques and Applications Conference CTAC-2004 (Melbourne,2004).Special Edition of
the Electronic Supplement of the ANZIAMJ,46(E),pp.C426–C438 (2005)
16.Dossou,K.,Pierre,R.:A Newton–GMRES approach for the analysis of the postbuckling behavior of the solutions of the
von Kármán equations.SIAMJ.Sci.Comput.24,1994–2012 (2003)