ES 247 Fracture Mechanics http://imechanica.

org/node/7448 Zhigang Suo

THE J INTEGRAL

Fracture Mechanics, with or without Field Theory

For a crack in an elastic body subject to a load, the elastic energy stored in
the body is a function of two independent variables: the displacement of the load,
( )
and the area of the crack, U Δ, A . This function of two independent variables
can be determined by alternative methods:
1. Experimental method. Measure the load-displacement curve for a body
containing a pre-cut crack. During the measurement, the crack does not
extend. Integrate the load-displacement curve to obtain the elastic
( )
energy U Δ, A for the body with the fixed area of crack. Repeat the
procedure for the body with a pre-cut crack of a different area. This
method invokes no field theory of elasticity.
2. Computational method. Solve the boundary-value problem for a body
containing a pre-cut crack of a given area. Calculate the energy density at
every material particle in the body. Integrate the energy density overall
all material particles to obtain the elastic energy for the body with the
specific area of crack. Repeat the procedure for the body with a pre-cut
crack of a different length. This method requires that the body be
modeled by a field theory of elasticity.

The energy release rate is defined as

G=−
∂U Δ, A
.
( )
∂A
This definition of the energy release rate assumes that the body is elastic, but
invokes no field theory of elasticity. Indeed, the energy release rate can be
determined experimentally by measuring the load-displacement curves of
identically loaded bodies with cracks of different areas. No field need be
measured.

The J integral. Many materials, however, can be modeled with a field

theory of elasticity. When a material is modeled by such a field theory, we can
represent the energy release rate in terms of the field in the body:
( )
G = ∫ WN 1 − Ti Fi1 dL .
The integral is known as the J integral (J for James R. Rice). This lecture
describes the J integral, along with examples of calculation. Uses of the J integral
are often better appreciated in the context of individual applications, which we
will describe in later lectures. In current engineering practice, the field in an
elastic body is commonly determined by using finite element method. Once the
field is determined, the finite element code calculates the J integral, which gives
the energy release rate.

Nonlinear Field Theory of Elasticity

The J integral can be developed for both linear and nonlinear theory of
elasticity. The nonlinear theory will be used in class, and the linear theory of will
be used in a homework problem.

4/8/16 1
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

For a detailed development of the field P

theory of elasticity, see a
http://imechanica.org/node/7794. The main
ingredients of the theory are as follows.

Uniaxial deformation. We proceed A

with our subject incrementally, beginning with l
the simplest structure: a bar. When the bar is L
not subject to any force, the cross-sectional
area is A and the length is L . We will call this
state the reference state. The bar is then
subject to an axial force P, and deforms to a
new state, cross-sectional area a and length l. P
We will call this state the current state. The reference
experimentalist records the force as a function state current
of the length. For an elastic body, the force- state
length curve is reversible upon loading and
unloading.
Define the stretch, λ , as the length of the bar in the current state divided
by the length of the bar in the reference state:
length in current state l
λ= = .
length in reference state L
When dealing with large deformation, we must be specific about the area
used in defining the stress. Define the nominal stress, s, as the force applied to
the bar in the current state divided by the cross-sectional area of the bar in the
reference state:
force in the current state P
s= = .
area in the reference state A
The nominal stress is also known as the engineering stress, or the first Piola-
Kirchhoff stress.
When the bar elongates from length l to length l + δl , the force P does
work Pδl . Recall the definitions of stress and strain, P = sA and l = λL .
Consequently the work done by the force is Pδl = ALsδλ . Since AL is the volume
of the bar in the reference state, we note that
increment of work in the current state
sδλ = .
volume in the reference state
We say that the nominal stress and the stretch are work-conjugate.
We assume that the bar is made of an elastic material. In the current
state, denote the elastic energy of the bar by U. The work done by the force
equals the elastic energy, δ U = Pδ l . That is, the elastic energy is the area under
the force-length curve.
Define the nominal density of elastic energy, W, as the elastic energy in
the bar in the current state divided by the volume of the bar in the reference
state:
elastic energy inthe currentstate U
W= = .
volumein thereferencestate AL
( )
The nominal density of elastic energy is a function of the stretch, W λ .

4/8/16 2
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

For the bar made of an elastic material, the work done by the force equals
the elastic energy, δ U = Pδ l . Dividing both sides by the volume of the bar in the
reference state, AL, we obtain that
δ W = sδλ .
According to calculus, this expression is equivalent to

s=
dW λ ( ).
dλ
Once we experimentally determined the force-length curve, we divide the force by
A and divide the length by L. This procedure scales the force-length curve to the
stress-stretch curve. The area under the stress-stretch curve gives the nominal
( )
density of elastic energy, W λ .

A body is a sum of many small parts. We next generalize the above

procedure to a body of an arbitrary shape and subject arbitrary load. The field
theory regards the body as a sum of many small pieces. Each small piece is called
a material particle, and undergoes a homogeneous deformation. The deformation
of the body is in general inhomogeneous—that is, the deformation varies from
one part of the body to another part.
In what follows, we will first describe the homogeneous deformation of an
individual small piece, and then describe the inhomogeneous deformation of a
body.

Homogeneous deformation in three dimensions. We first

consider a small part of the body. We list the three ingredients of the theory:
geometry of deformation, balance of forces, and conditions of thermodynamic
equilibrium.
Geometry of deformation. In the undeformed state, the part is a unit
cube, with edges coinciding with the coordinates. Subject to forces on the faces of
the part, the part undergoes a homogeneous deformation. In the deformed state,
the part deforms into a parallelepiped.
The deformation maps the three orthogonal edges of the cube to three
edges of the parallelepiped. The edges of the parallelepiped are three vectors,
noted as Fi1 , Fi2 , and Fi 3 . The first subscript indicates the component of each
vector, and the second subscript differentiates the three edges. The nine
quantities FiK together are known as the deformation gradient.

si3
si2
Fi3 si1

Fi2
Fi1
1 si1

1 1
si3 si2

undeformed state deformed state

4/8/16 3
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

Balance of forces. Acting on the three pairs of parallel faces of the

parallelepiped are forces si1 , si2 , and si 3 . The first subscript indicates the
component of each vector, and the second subscript differentiates the three pairs
of parallel faces. The nine quantities siK together are known as the nominal
stress.
Let N be the unit vector normal to an element of an area in the
undeformed body. Define the nominal traction T by the force acting on the
element of the surface in the deformed body divided by the area of the element in
the undeformed body. The balance of forces relates the nominal traction to the
nominal stress:
Ti = siK N K .

Condition of thermodynamic equilibrium. We consider a small part of

the body in thermodynamic equilibrium with external forces and with an
environment of constant temperature. We will drop temperature from the list of
variables. Let W be the Helmholtz free energy of the parallelepiped. We say a
material is elastic if the Helmholtz free energy is a function of the deformation
( )
gradient, W F . This function is an input to the field theory of elasticity. W is
also known as the nominal density of elastic energy.
When the part is in a state of thermodynamic equilibrium, the nominal
stress relates to the deformation gradient as

siK =
∂W F ( )
.
∂FiK
This equation of states generalizes Hooke’s stress-strain relation.

Inhomogeneous deformation. We call each small part of the body a

material particle. As the body deforms, each material particle moves in space.
Now focus on a particular material particle in the body. When the body is in the
undeformed state, the material particle is at a place in the space, of coordinate X.
When the body is in a deform state, the same material particle moves to a
( )
difference place in space, of coordinate x. The function x X fully describes the
deformation of the body. The deformation gradient relates to the deformation as

FiK =
∂xi X
.
( )
∂X K
When the deformation of the body is homogeneous, the deformation gradient is
the same for all material particles in the body. In general, the deformation of the
body is inhomogeneous, so that the deformation gradient varies from one
material particle to another. That is, the deformation gradient as a function of
( )
material particle, FiK X .
In general, the stress in the deformed body is inhomogeneous, and we
( )
write the stress as a function of material particles, siK X . The balance of forces
requires that the field of stress satisfy

4/8/16 4
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

( ) =0.
∂siK X
∂X K
Here we neglect body force and inertia.

Elastic energy of the body. When a body is in equilibrium with a set

of applied forces, the field in the body is in general inhomogeneous. The
Helmholtz free energy of the body is the sum of energy of all material particles:
∫
U = WdV .
The energy density is of the body in the deformed state, but the integration is
carried over the volume of the body in the reference state.

A composite system that does not receive work from the rest of
the world. Consider a body subject to a constant force P, with Δ being the
displacement of the force. We can always picture the constant force as a hanging
weight. The potential energy of the weight is −PΔ . The body and the weight
together form a composite system. The composite exchanges with the
environment by heat, but not by work. The Helmholtz free energy of the
composite is
Π = U − PΔ .
This quantity is known as the potential energy of the body in mechanics.
The general idea is to form a composite system that does not receive any
work from the rest of the world. We have considered two examples:
1. When a body is loaded by a fixed displacement, the body itself is
composite system, and the free energy of the system is U.
2. When a body is loaded by a hanging weight, the body and the weight
together form a composite system, and the free energy of the system is
U − PΔ .
We can also consider yet another example. When a body in series with a spring is
loaded by a fixed displacement, the body and the spring together form a
composite system, and the free energy of the system is the sum of that of the body
and that of the spring.

Energy Release Rate

Express energy release rate in terms of field variables. Model a

body with the nonlinear theory of elasticity. Under the plane-strain or plane-
stress conditions, we represent the body by a region in the plane, and the
boundary of the body by a curve in the plane. Consider a sample of unit thickness.
Each point on the boundary is subject to constant nominal traction Ti . The
potential energy per unit thickness of the body is
Π = W (F )dA − Ti (X )xi (X )dL
∫ ∫
The first integral extends over the body, and the second integral extends over the
boundary of the body.
Now compare two specimens, one having a crack of length C in the
undeformed state, and other having a crack of length C + δC in the undeformed
state. The cracks are in the direction of coordinate X 1 . The two specimens are of
the same general shape. In the deformed states, the cracks in the two specimens

4/8/16 5
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

are traction-free, and external boundaries of

the specimens are subject to the same
tractions. The two specimens have different
( ) (
potential energies, Π C and Π C + δC . )
In the undeformed state, we translate C
one specimen relative to the other by a
distance δC , so that the tips of the cracks in
the two specimens coincide. Consider an
element of the boundary dL. Let N be the unit
vector normal to the element. As one
specimen translates relative to the other by a
distance δC , the element dL sweeps across a C + δC
region in the shape of a parallelogram. The
area of the parallelogram is δCN 1dL . In the
deformed states of the two specimens, the
parallelogram contributes a difference in the N
elastic energies in the two specimens,
−WN 1δCdL . dL
δC
Consider a set of material particles. In
the undeformed state, the set of particles forms
a straight segment δC . In the deformed state,
the set of material particles forms a segment
represented by vector Fi1δC . The applied force acting on the element is Ti dL .
Consequently, the applied force acting on the element dL gives a difference in the
potential energies of the two specimens, Fi1δCTi dL .
Combining the above contributions, we obtain the difference in the
potential energies of the two specimens:
( ) ( )
Π C + δC − Π C = − ∫ WN 1 δ C dL + ∫ Ti Fi1 δ C dL .
The first term is due to the elastic energy of the body, and the second term is due
to the potential energy of the applied traction.
By definition, the energy release rate is

G=−
(
Π C + δC − Π C ) ( )
.
δC
The partial derivative means that we vary the length of the crack, but keep the
loading conditions fixed. Comparing the above two expressions, we obtain that
(
G = ∫ WN 1 − Ti Fi1 dL . )
This equation expresses the energy release rate in terms of field variables.

Properties of the J Integral

The J integral. Let us focus on the integral:

(
J = ∫ WN 1 − Ti Fi1 dL . )
The path of integration is a curve, drawn in the body in the undeformed state,
connecting two material particles. The unit vector N is normal to the path of
integration.

4/8/16 6
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

Recall the relations in the theory of

( ) ( )
elasticity: FiK = ∂xi X /∂X K , siK = ∂W F /∂FiK
and Ti = siK N K . The J integral can also be
written as
B
⎡ ∂W (F ) ∂xi (X )⎤
∫
J = ⎢WN 1 − N K
⎣ ∂FiK ∂X 1 ⎦
⎥dL .
The integral is a functional of x (X ), the field
A
that describes the deformation of the body.
The dummy variable of integration is X, the
coordinates of material particles in the
undeformed body.

Yet another alternative expression of the J integral. A material

particle occupies a place of coordinate X in the undeformed state, and occupies a
place of coordinate x in the deformed state. The displacement of the material
particle is
( )
u= x X −X.
Thus
( ) = ∂x ( X) − δ
∂ui X i
.
i1
∂X 1 ∂X 1
Inserting this expression into the definition of the J integral, we obtain that
# ∂ui X &
J = ∫ %WN 1 − Ti
( )
( dL .
% ∂X (
$ 1 '
In reaching this expression, we have used the following relation:
∫ T1 dL = 0 .
This relation results from the balance of forces.

Divergence theorem. In the following development, we will need a

result in calculus: the divergence theorem. Consider a region in a plane. Let
( )
f X be a smooth function defined in the region, and N K be the unit vector
normal to the curve around the region. The divergence theorem is
( ) dA =
∂f X
∫ ∂ XK ∫ fN K
dL .

The integral on the left is over the area of a region, and integral on the right is
over the curve surrounding the region.

The J integral is path-independent.

The two material particles A and B can be
connected by numerous paths. For any two B
paths, so long as no singularity exists between
them, the J integrals along the two paths are
identical.
A

4/8/16 7
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

To prove this statement, reverse the direction of one path, so that the two
paths form a closed contour. This closed contour is used as a path of integration.
We need to prove that the J integral vanishes if the closed contour encloses no
singularity. The proof invokes the divergence theorem, as well as the field theory
of elasticity. We show that the second part of the J integral equals the first part:
∫ Ti Fi1 dL = ∫ N K siK Fi1 dL recall Ti = siK N K

=∫
(
∂ siK Fi1 ) dA divergence theorem
∂ XK
⎛ ∂s ∂Fi1 ⎞
= ∫ ⎜ iK Fi1 + siK dA
⎝ ∂ XK ∂ X K ⎟⎠
∂Fi1 ∂siK
= ∫ siK dA recall =0
∂ XK ∂ XK

= ∫ siK
∂2 x i
dA recall FiK =
( )
∂xi X
.
∂ X1 ∂ X K ∂ XK

=∫
∂W F ( ) ∂2 x i
dA recall siK =
∂W F( )
∂FiK ∂ X1 ∂ X K ∂FiK

=∫
( )
∂W F ∂FiK
dA
∂FiK ∂ X 1

=∫
∂W F ( ) dA
∂ X1
( )
= ∫ W F N 1 dL divergence theorem

Undeformed body Deformed body

The J integral associated with a traction-free crack. Consider a

crack in a body subject to a load. The faces of the crack are traction-free. In the
undeformed body, a path of integration starts from a material particle on the
lower face of the crack, surrounds the tip of the crack, and ends at a material

4/8/16 8
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

particle on the upper face of the crack. Because the path encloses the singular tip
of the crack, the J integral does not vanish. We can make two statements:
1. Given the two material particles on the faces of the crack, the J integral is
independent of the path.
2. The J integral is the same independent of the choice of the material
particles on the two faces of the crack.

Use the J integral to calculate energy release rate in the finite-

element method. For a crack in an elastic body subject to a load, the energy
release rate is defined as the decrease in the elastic energy of the body associated
with extension of the crack by unit area, while the load is held fixed. We have
used this definition directly to calculate the energy release rate for a few cases. In
each case, we somehow obtain the elastic energy in the body as a function of the
area of the crack, and then take the partial derivative of this function with respect
to the area of the crack.
This procedure, however, is difficult to apply when we solve the
boundary-value problem by using the finite-element method. Let us imagine how
we might go about it. For a given area of the crack, we use the finite-element
method to determine the field in the body, and then integrate to obtain the elastic
energy stored in the body. We then change the area of the crack slightly, solve a
new boundary-value problem, and then obtain the elastic energy in the body. We
then take the ratio between the difference in the elastic energy and the difference
in the area of crack. Here are two obvious difficulties:
• In each boundary-value problem, the field is singular around the front of
the crack.
• The elastic energies in the two bodies are large quantities compared to
their difference. To calculate the energy release rate this way would
require us to determine the field very accurately.

Both of the above difficulties are avoided if we use the J integral to

calculate the energy release rate. The integrand involves the field in the body of a
fixed crack. The path of integration needs to enclose the tip of crack, but can be
chosen to be far away from the tip of the crack.

Examples

A long crack in a strip of a material pulled by two rigid grips.

Rivlin and Thomas (1953) described this experimental setup. They obtained the
energy release rate directly from its definition. The result is
G = HW (λ ),
where H is the width of the undeformed strip, and λH is the width of the
deformed stripe. The quantity HW (λ ) can be determined from experimentally
measured load-displacement curve of a strip with no crack, pulled by two rigid
grips.

4/8/16 9
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

H λH

The same result can be obtained by using the J integral along the dotted
lines indicated in the undeformed strip (Rice, 1968). Use the J integral in the
form
# ∂ui X &
%
J = ∫ WN 1 − Ti
( )
( dL .
% ∂X (
$ 1 '
The only nonzero contribution comes from the vertical line in the strip ahead of
the crack.

Layered materials. In the above, the material is taken to be

homogeneous. When the two material particles are of the same deformation
gradient, they have the same free energy function—that is, the nominal density of
free energy is a function the deformation gradient:
W = W (F ) .
The function does not depend on X explicity.
For a layered material, which is homogenous in the X 1 -direction, but
inhomogeneous in the X 2 -direction, the energy density takes the form
W = W (F, X 2 ).
Going through the same steps, you can confirm that the J integral is still path-
independent.
Finally consider a generally inhomogeneous material, for which the
energy density takes the form
W = W (F, X 1 , X 2 ).
Going through the same steps, you will find that the J integral is path-dependent.

film

substrate

Detaching highly stretchable materials. Jinda Tang and Jianyu Li

are studying the delamination of highly stretchable materials, such as elastomers
and gels, using a specific experimental setup. Attach a film to a substrate, and
leave a detached region as a pre-existing crack. Both materials are elastic and
highly stretchable. When a force pulls the substrate, the tip of the crack blunts.
When the force reaches some critical value, the crack runs. The length of the
crack is large compared to the thicknesses of the film and the substrate, so that
the crack runs in a steady state, and the energy release rate is independent of the
length of the crack. We wish to relate the energy release rate of the crack to the
applied force.

4/8/16 10
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

First measure the stress-stretch relations of the film and the substrate
separately. When a uniaxial force pulls the film, the nominal stress S f in the film

()
relates to the stretch λ through a curve S f λ , and the area under this curve is
the nominal density of the free energy in the film as a function of the stretch,
()
W f λ . Similarly, when a uniaxial force pulls the substrate, the nominal stress
Ss in the substrate relates to the stretch λ through a curve Ss λ , and the area()
under this curve is the nominal density of the free energy in the film as a function
( )
of the stretch, Ws λ .
When the laminate is in the undeformed state, the thickness of the film is
H f , the thickness of the substrate is H s , and the width of the film and substrate
is B. In the deformed state, a force P pulls the substrate. Far behind the tip of the
crack, the film is stress-free, and the substrate is in a state of uniaxial stress, with
the stretch λ ! determined by
( )
H s Ss λ ! = P .
Far ahead the tip of the crack, the film and substrate are attached and both have
the same stretch λ !! , determined by
( ) ( )
H f S f λ !! + H s Ss λ !! = P .
Now compare two specimens. Both specimens have the total length L in
the reference state. One specimen has a crack of length C in the undeformed
state, and has potential energy
Π C = C $% BH sWs λ " − P λ "&' + L − C $% BH f W f λ "" + BH sWs λ "" − P λ ""&'
( ) ( ) ( ) ( ) ( )
in the deformed state.
The other specimen has a crack of length C + δC in the undeformed state,
and has potential energy
Π C + δC = C + δC $% BH sWs λ " − P λ "&' + L − C − δC $% BH f W f λ "" + BH sWs λ "" − P λ ""&'
( ) ( ) ( ) ( ) ( ) ( )
in the deformed state.
By definition, the energy release rate is the reduction in the potential
energy associated with unit increase of the crack area:

G=−
(
Π C + δC − Π C
,
) ( )
BδC
giving that
P
( ) ( ) ( )
G = H f W f λ !! + H sWs λ !! − H sWs λ ! + λ ! − λ !! .
B
( )
In writing the potential energy of each specimen, we have assumed that
the detached substrate is in a homogeneous state of stretch λ ! , and the attached
laminate is in a homogeneous state of λ !! . In the region around the tip of the
crack, however, the field is inhomogeneous. Because the crack is in a steady
state, the inhomogeneous field in this region is identical in the two specimens.
Consequently, this inhomogeneous field does not affect the calculation of the
energy release rate.
We can also use the J integral to obtain the above expression for the
energy release rate.

4/8/16 11
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

These relations together relate the energy release rate to the applied force
()
P. We measure the stress-stretch curves of the two materials, S f λ and Ss λ . ()
The areas under these curves give the nominal densities of energy of the two
() ( )
materials, W f λ and Ws λ . The applied force P determines the stretch λ ! in

( )
the detached substrate by H s Ss λ ! = P , and determines the stretch λ !! in the

( ) ( )
attached laminate by H f S f λ !! + H s Ss λ !! = P .

A thin, compliant film on a substrate. In a limiting case that the

substrate is much stiffer than the substrate, the attached film does not constrain
the substrate, so that the stretch in the detached substrate is nearly the same as
the stretch in the attached laminate, λ ! ≈ λ !! . In this limiting case, the energy
release rate reduces to
G = H fWf λ . ()
When the stretch is of order unity, the free energy density is on the order of the
elastic modulus, E. For a hydrogel, a representative value is E = 104 Pa. For the
film of thickness H f = 0.1 mm, the energy release rate is on the order of 1 J/m2.
For given materials, the laminate will not detach if the film is thin and stretch is
small.

A stiff film on a compliant substrate. In this limiting case, the

attached laminate does not deform, λ !! = 1, so that the energy release rate is
entirely due to the detached substrate:
G = H s #$ Ss λ ! λ ! −Ws λ ! %& .
( ) ( )
For experimental details, see Jingda Tang, Jianyu Li, Joost J. Vlassak,
Zhigang Suo. Adhesion between highly stretchable materials. Soft Matter 12,
1093-1099 (2016). http://www.seas.harvard.edu/suo/papers/350.pdf

Historical Notes
Eshelby studied the change in the potential energy associated with the
movement of a singularity. He called this change the force on the singularity, and
expressed this force as an integral. He did not relate his idea to the Griffith
theory of cracks. Rice and Cherepanov independently discovered the J integral
within the context of fracture mechanics.

• J.D. Eshelby, The force on an elastic singularity. Phil. Trans. Roy. Soc.
London A 244, 87-111 (1951).
• J.D. Eshelby, The continuum theory of lattice defects. Solid State Physics
3, 79-144 (1956).
• J.R. Rice, A path independent integral and the approximate analysis of
strain concentration by notches and cracks. J. Appl. Mech. 35, 379-386
(1968).
http://esag.harvard.edu/rice/015_Rice_PathIndepInt_JAM68.pdf
• J.R. Rice, Mathematical analysis in the mechanics of fracture, Chapter 3
of Fracture: An Advanced Treatise (Vol. 2, Mathematical Fundamentals)
(ed. H. Liebowitz), Academic Press, N.Y., 1968, pp. 191-311.
http://esag.harvard.edu/rice/018_Rice_MathAnalMechFract_68.pdf

4/8/16 12
ES 247 Fracture Mechanics http://imechanica.org/node/7448 Zhigang Suo

• G.P. Cherepanov, Crack propagation in continuous media. Journal of

Applied Mathematics and Mechanics 31, 503 (1967).
• G.P. Cherepanov, Cracks in solids. International Journal of Solids and
Structures 4, 811-831 (1968).
• S.J. Chang, Path-independent integral for rupture of perfectly elastic
materials. Journal of Applied Mathematics and Physics 23, 149-152
(1972).

4/8/16 13