Non-linear analysis of double cantilever beam

Abstract

Non-linear fracture analysis of the Double Cantilever Beam (DCB) specimen is performed by using the J-integral approach. For this purpose, methods of Mechanics of materials are applied. Closed form analytical solutions of the J-integral are obtained by using two non-linear beam models (ideally elastic–plastic model and model with linear hardening). It is assumed that the material has the same properties in tension and compression. The fracture behaviour is analyzed at each typical stress level in the DCB in the process of increasing of the external load. The influence of the material non-linearity on the fracture behaviour is discussed. The effect of slope of the stress–strain curve in the portion with linear hardening on the fracture behaviour is evaluated. The analysis revealed that the J-integral value increases when the slope increases. The analytical approach developed in the present paper is very useful for parametric analyses, since the simple formulae derived capture the fundamentals of non-linear fracture in the DCB configuration.

Keywords

Non-linear fracture beam theory analytical solution

1. Introduction

Due to their superior specific strength and stiffness, fibre reinforced composite laminates are ideal for weight critical structures. With increasing use of laminates in load-bearing structural components, the need to understand and predict the failure behaviour of these materials has also increased.

Delamination, i.e. the separation of the layers in laminated materials, is the most severe type of damage which reduces considerably the strength and the stiffness of the structure and may, in some cases, lead to complete failure [1,2]. The most common causes of delamination are defects induced in a result of careless manufacturing. Similar defects may also arise during lifetime of the laminated structure due to low-velocity impact events, cyclic loading, thermo-mechanical loads, and in-service incidents such as tool drops or careless handling during maintenance. Delamination cracks may be very difficult to detect since they are often embedded within the structure, which makes them especially insidious. Therefore, the delamination fracture mechanisms should be carefully considered in evaluation of laminated structures for durability and damage tolerance.

The literature on the delamination problem is wide [3–19]. Usually specimens of beam type are used for characterization of delamination fracture behavior of laminated composite structures (one of the most commonly used configurations is the Double Cantilever Beam (DCB) specimen). The delamination is studied in terms of the strain energy release rate by applying the concepts of linear-elastic fracture mechanics. For this purpose, various analytical solutions are obtained using the linear-elastic beam theory. The basic assumption of these solutions is that the relationship between the stresses and strains can be expressed by the Hook’s law. However, in the cases of composites with high fracture toughness, the applicability of linear-elastic fracture mechanics is limited, because the plastic deformation usually begins prior to onset of macroscopic crack growth from the initial crack tip position. In such cases, the fracture analysis should take into account the material non-linearity.

The objective of the present paper is to perform a non-linear fracture analysis of DCB specimen. The integration path independent J-integral is used as a fracture characterizing parameter. Two non-linear stress–strain curves are considered in the analysis. Closed form analytical solutions of J-integral are obtained using a model which is based on Mechanics of materials. Effects of stress–strain curves parameters on non-linear fracture are evaluated.

Fig. 1.

The DCB configuration.

2. Non-linear fracture analysis

The DCB configuration analyzed in the present paper is illustrated schematically in Fig. 1. There is a crack in the mid-plane of the beam. The external load consists of two moments applied at the free end of the crack arms. The fracture behavior is analyzed by using the J-integral approach. The J-integral is written as [20] $J = \int_{Γ} [u_{0} cos α - (p_{x} \frac{\partial u}{\partial x} + p_{y} \frac{\partial v}{\partial x})] d s,$ (1) where Γ is a contour of integration going around the crack from one crack face to the other in the counter clockwise direction, $u_{0}$ is the strain energy density, α is the angle between the outwards normal vector to the contour of integration and the crack direction, $p_{x}$ and $p_{y}$ are the components of the stress vector, u and v are the components of the displacement vector with respect to the crack tip coordinate system, $x y$ , and $d s$ is a differential element along the contour Γ.

The J-integral is invariant with respect to the integration contour [20]. Therefore, we decide the integration contour, Γ, to coincide with the beam contour as shown in Fig. 1. In this way, the J-integral solution is substantially facilitated. For instance, the J-integral value is non-zero only in segments A and B of the contour (Fig. 1). Due to the symmetry conditions, it is enough to solve the J-integral in segment A and to double the result obtained.

Fig. 2.

Stress–strain curve of an elastic–perfectly plastic material.

2.1. Ideally elastic–plastic model

We assume that the mechanical response of the beam can be described by the ideally elastic–plastic diagram for both tension and compression (Fig. 2). When the external load magnitude is relatively low, the beam deforms in linear-elastic stage, i.e. the Hook’s law is valid. In this case, for the components of the J-integral in segment A (Fig. 1), we have $p_{x} = \frac{M}{I} y_{1}, p_{y} = 0,$ (2) where $I = \frac{b h^{3}}{12}$ is the principal moment of inertia of the upper crack arm cross-section (here b is the width of beam cross-section). The coordinate, $y_{1}$ , varies in the interval $[- \frac{h}{2}; + \frac{h}{2}]$ , and $d s = - d y_{1}$ , $cos α = - 1$ .

The partial derivative, $\frac{\partial u}{\partial x}$ , is obtained as $\frac{\partial u}{\partial x} = - ε_{x} = - \frac{p_{x}}{E} = - \frac{M}{E I} y_{1},$ (3) where E is the longitudinal modulus of elasticity. The following formula for strain energy density is used: $u_{0} = \frac{p_{x}^{2}}{2 E} = \frac{1}{2 E} {[\frac{M}{I} y_{1}]}^{2} = \frac{M^{2}}{2 E I^{2}} y_{1}^{2} .$ (4)

By substitution of (2), (3) and (4) in (1), we obtain $J = 2 \int_{\frac{h_{1}}{2}}^{- \frac{h_{1}}{2}} {[\frac{M^{2}}{2 E I^{2}} y_{1}^{2} \cdot (- 1) - \frac{M}{I} y_{1} (- \frac{M}{E I} y_{1})]} (- d y_{1}) = \frac{12 M^{2}}{E b^{2} h^{3}} .$ (5)

It is known that at linear-elastic behaviour of the material, the J-integral value coincides with the strain energy release rate. This fact is used here to verify solution (5). Indeed, (5) coincides with the expression for the strain energy release rate in DCB [21].

Fig. 3.

Plastic zones and diagrams of normal stresses and linear strains in the crack arms.

The crack growth will begin when $J = J_{C},$ (6) where $J_{C}$ is the fracture toughness. From (5) and (6), we obtain the critical magnitude of the external moment (at the onset of crack growth) $M_{C} = \frac{J_{C} E b^{2} h^{3}}{12} .$ (7)

However, if the fracture toughness is high, the plastic deformation will begin prior to the onset of crack growth. When the moment gradually increases, plastic deformation first will begin in the remotest edges of the crack arms, where the normal stress will attain the yield stress limit, $f_{y}$ : $\frac{M_{PL}}{W} = f_{y},$ (8) where $W = \frac{b h^{2}}{6}$ (9) is the second modulus of the crack arm cross-section. From (8) and (9), we obtain the critical moment $M_{PL} = f_{y} \frac{b h^{2}}{6} .$ (10)

If $M_{C} > M_{PL}$ , the plastic deformation will begin before the onset of crack growth from the initial crack tip position. Therefore, the material non-linearity has to be taken into account in the fracture analysis. We assume that plastic zones are formed in the crack arms as illustrated in Fig. 3. The magnitude of the moment is $M_{ep}$ ( $M_{ep}$ is higher than $M_{PL}$ ). The value of the J-integral in segment A is obtained by integration along the elastic and plastic zones, i.e. $J = J_{el} + J_{pl},$ (11) where the J-integral in the elastic zone of the crack arm cross-section is written as $J_{el} = \int_{\frac{h_{el}}{2}}^{- \frac{h_{el}}{2}} {u_{0_{el}} cos α - [p_{x_{el}} \frac{\partial u}{\partial x}]} d s,$ (12) where $h_{el}$ is the height of the elastic zone (Fig. 3).

The normal stresses in the elastic zone are distributed linearly $p_{x_{el}} = \frac{M_{ep}}{I} y_{1} .$ (13)

The partial derivative, $\frac{\partial u}{\partial x}$ , is written as: $\frac{\partial u}{\partial x} = - ε_{x} = - \frac{1}{ρ_{ep}} y_{1},$ (14) where the curvature of the crack arm in elastic–plastic stage of deformation can be obtained using the following formula from Mechanics of materials [22] $\frac{1}{ρ_{ep}} = \frac{M_{PL}}{E I \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}} .$ (15)

The strain energy density is found as $u_{0_{el}} = \frac{1}{2} p_{x_{el}} ε_{x} = \frac{M_{PL} M_{ep}}{2 E I^{2} \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}} y_{1}^{2} .$ (16) By substitution of (13)–(16) in (12), we obtain $J_{el} = \frac{6 M_{PL} M_{ep}}{E b^{2} h^{6} \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}} h_{el}^{3} .$ (17) The height of the elastic zone is determined using a formula from Mechanics of materials [22] $h_{el} = h \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}} .$ (18) By substitution of (18) in (17), we obtain $J_{el} = \frac{6 M_{ep} (3 M_{PL} - 2 M_{ep})}{E b^{2} h^{3}} .$ (19)

The J-integral in the plastic zones of the crack arm cross-section is written as $J_{pl} = 2 \int_{\frac{h_{el}}{2}}^{\frac{h}{2}} {u_{0_{pl}} cos α - [p_{x_{pl}} \frac{\partial u}{\partial x}]} d s,$ (20) where Fig. 2 takes into account the fact that there are two plastic zones in the crack arm cross-section (Fig. 3).

Fig. 4.

Strain energy density in the plastic zones of the crack arms.

The normal stresses in the plastic zones are equal to the yield stress limit, i.e. $p_{x_{pl}} = f_{y} .$ (21)

The strain energy density in the plastic zones of the crack arm cross-section is equal to the area enclosed by the stress–strain curve [23–26] (Fig. 4): $u_{0_{pl}} = \frac{1}{2} f_{y} ε_{y} + f_{y} (ε - ε_{y}) = \frac{f_{y}}{2} (\frac{2 M_{PL}}{E I \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}} y_{1} - \frac{f_{y}}{E}) .$ (22) By substitution of (14), (15), (21) and (22) into (20), we obtain $J_{pl} = \frac{f_{y}^{2}}{2 E} (h - h_{el}) .$ (23) Further, (18) is substituted in (23), i.e. $J_{pl} = \frac{f_{y}^{2} h}{2 E} (1 - \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}) .$ (24) By taking into account (10), we write (24) in the form $J_{pl} = \frac{18 M_{PL}^{2}}{E b^{2} h^{3}} (1 - \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}) .$ (25) Substitution of (19) and (25) in (11) yields $J = \frac{6 M_{ep} (3 M_{PL} - 2 M_{ep})}{E b^{2} h^{3}} + \frac{18 M_{PL}^{2}}{E b^{2} h^{3}} (1 - \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}) .$ (26) The final solution is obtained by doubling (26) in view of the symmetry conditions, i.e. $J = \frac{12 M_{ep} (3 M_{PL} - 2 M_{ep})}{E b^{2} h^{3}} + \frac{36 M_{PL}^{2}}{E b^{2} h^{3}} (1 - \sqrt{3 - 2 \frac{M_{ep}}{M_{PL}}}) .$ (27) It should be mentioned that formula (27) is valid if $M_{PL} < M_{ep} < M_{S},$ (28) where $M_{S} = \frac{b h^{2}}{4} f_{y}$ (29) is the limit bending moment at the plastic collapse of the crack arms.

Fig. 5.

Stress–strain curve with linear hardening.

2.2. Elastic–plastic model with linear hardening

We analyze the fracture behaviour of the DCB configuration using a stress–strain curve with linear hardening (Fig. 5). At low magnitudes of the external moments, the J-integral solution is expressed by formula (5).

Fig. 6.

Plastic zones and diagrams of normal stresses and linear strains in crack arms.

If $M = M_{ep} > M_{PL},$ (30) the material non-linearity has to be accounted for in the fracture analysis. When $M_{ep} > M_{PL}$ , two symmetric plastic zones are formed in the crack arms (Fig. 6). Here, again the J-integral solution can be obtained using formula (11). In the elastic zone of the crack arm cross-section, the J-integral solution is expressed by (12). The normal stresses in the elastic zone are obtained by (13). The partial derivative again can be found by formulae (14). However, the curvature of crack arm, $1 / ρ_{hr}$ , is different. We will use the following equation in order to find the curvature [22]: $M_{ep} = \frac{M_{PL}}{2} {2 \frac{E_{1}}{E} \frac{ρ_{y}}{ρ_{hr}} + (1 - \frac{E_{1}}{E}) [3 - {(\frac{ρ_{hr}}{ρ_{y}})}^{2}]},$ (31) where $E_{1}$ is the modulus of elasticity in regime of hardening (Fig. 5). The curvature radius, $ρ_{y}$ , when the normal stress in the remotest edges of the crack arms attains the material yield stress limit is obtained by the following formula from Mechanics of materials [22]: $ρ_{y} = \frac{E h}{2 f_{y}} .$ (32) Equation (31) is derived by considering the equilibrium of the elementary forces in a rectangular cross-section in elastic–plastic stage of deformation [22]. We solve (31) with respect to the curvature, $1 / ρ_{hr}$ , using the method for solution of cubic equation [27]. The result is $\frac{1}{ρ_{hr}} = \frac{1}{ρ_{y}} (Q_{1} + Q_{2} - \frac{P}{3}),$ (33) where $\begin{array}{rclr} Q_{1} = \sqrt[3]{- {(\frac{P}{3})}^{3} - \frac{n - 1}{4 n} + \sqrt{Q}}, & (34) \\ Q_{2} = \sqrt[3]{- {(\frac{P}{3})}^{3} - \frac{n - 1}{4 n} - \sqrt{Q}}, & (35) \\ Q = \frac{2 P^{3}}{27} (\frac{n - 1}{4 n}) + {(\frac{n - 1}{4 n})}^{2}, & (36) \\ P = \frac{3 M_{PL} (1 - n) - 2 M_{ep}}{2 M_{PL} n}, & (37) \\ n = \frac{E_{1}}{E} . & (38) \end{array}$

It should be mentioned that at $n = 1$ , from (33)–(38), we obtain $\frac{1}{ρ_{hr}} = \frac{M_{ep}}{E I},$ (39) which is the curvature of the crack arm in linear-elastic stage of deformation.

We replace $1 / ρ_{ep}$ by $1 / ρ_{hr}$ in (14) to obtain the partial derivative, $\frac{\partial u}{\partial x}$ , i.e. $\frac{\partial u}{\partial x} = - ε_{x} = - \frac{1}{ρ_{hr}} y_{1} .$ (40)

The strain energy density is found as $u_{0_{el}} = \frac{1}{2} p_{x_{el}} ε_{x} = \frac{M_{ep}}{2 I} \frac{1}{ρ_{hr}} y_{1}^{2} .$ (41) By substitution of (13), (40) and (41) in (12) and solving the integral, we obtain $J_{el} = \frac{M_{ep} h_{el}^{3}}{24 I} \frac{1}{ρ_{hr}} .$ (42) The height of the elastic zone is determined using a formula from Mechanics of materials [22] $h_{el} = h \frac{ρ_{hr}}{ρ_{y}},$ (43) where $ρ_{y}$ is obtained by formula (32).

Fig. 7.

Strain energy density in the plastic zones of the crack arms in stage of hardening.

In the plastic zones of the crack arm cross-section, the J-integral is expressed by (20). The normal stresses in the plastic zones are obtained as (Fig. 5) $p_{x_{pl}} = f_{y} + E_{1} (ε - ε_{y}) .$ (44)

The strain energy density in the plastic zones of the crack arm cross-section is equal to the area enclosed by the stress–strain curve [23–26] (Fig. 7), i.e. $u_{0_{pl}} = \frac{1}{2} f_{y} ε_{y} + \frac{1}{2} [2 f_{y} + E_{1} (ε - ε_{y})] (ε - ε_{y}) .$ (45) By substitution of $ε = \frac{1}{ρ_{hr}} y_{1}$ and $ε_{y} = \frac{f_{y}}{E}$ in (45), we obtain $u_{0_{pl}} = \frac{f_{y}^{2}}{2 E} + [f_{y} + \frac{E_{1}}{2} (\frac{y_{1}}{ρ_{hr}} - \frac{f_{y}}{E})] (\frac{y_{1}}{ρ_{hr}} - \frac{f_{y}}{E}) .$ (46) We substitute (40), (44) and (46) into (20) and solve the integral. The result is $J_{pl} = \frac{f_{y}^{2}}{2 E} (h - h_{el}) (1 - \frac{E_{1}}{E}) + \frac{E_{1}}{24 ρ_{hr}^{2}} (h^{3} - h_{el}^{3}) .$ (47) By substitution of (42) and (47) in (11), we obtain $J = \frac{M_{ep} h_{el}^{3}}{24 I} \frac{1}{ρ_{hr}} + \frac{f_{y}^{2}}{2 E} (h - h_{el}) (1 - \frac{E_{1}}{E}) + \frac{E_{1}}{24 ρ_{hr}^{2}} (h^{3} - h_{el}^{3}) .$ (48) The final solution is obtained by doubling (48) in view of the symmetry conditions, i.e. $J = \frac{M_{ep} h_{el}^{3}}{12 I} \frac{1}{ρ_{hr}} + \frac{f_{y}^{2}}{E} (h - h_{el}) (1 - \frac{E_{1}}{E}) + \frac{E_{1}}{12 ρ_{hr}^{2}} (h^{3} - h_{el}^{3}),$ (49) where $1 / ρ_{hr}$ and $h_{el}$ are obtained be formulae (33) and (43), respectively.

Fig. 8.

The J-integral value plotted against the proportion, $E_{1} / E$ , at two magnitudes of the external moment, $M_{ep}$ .

It should be mentioned that at $E_{1} = E$ and $E_{1} = 0$ Eq. (49) transforms into the linear-elastic solution (5) and ideally elastic–plastic solution (27), respectively. These transformations are indications for the correctness of Eq. (49), since at $E_{1} = E$ and $E_{1} = 0$ the model with linear hardening transforms into linear-elastic and ideally elastic–plastic model, respectively.

The influence of the stress–strain curve slope in the portion with hardening on the non-linear fracture behaviour is analyzed. For this purpose, the J-integral is calculated by formula (49) at different values of the proportion, $E_{1} / E$ . The following dimensions of the crack arm cross-section are used: $b = 0.02 m$ and $h = 0.01 m$ . Calculations are performed for $M_{ep} = 40 Nm$ and $M_{ep} = 45 Nm$ , i.e. the external moment magnitude is higher than $M_{PL}$ (refer to formula (10)). Thus, the beam deforms in elastic–plastic stage. It is assumed also that $E = 5.378 \times 10^{11} Pa$ and $f_{y} = 1.034 \times 10^{8} Pa$ . The results of these calculations are presented graphically in Fig. 8. The diagrams in Fig. 8 indicate that increase of the proportion, $E_{1} / E$ , leads to increase of the J-integral value.

3. Conclusions

Non-linear fracture in the DCB configuration is studied theoretically using two models (ideally elastic–plastic model and model with linear hardening) which are based on the Mechanics of materials. The fracture analysis is carried-out using the J-integral approach. The integration contour is chosen to coincide with the beam contour in order to facilitate obtaining of the J-integral closed form analytical solutions. The influence of the slope of the stress–strain curve in the portion of linear hardening on the non-linear fracture is also analyzed. It is found that the increase of the slope leads to increase of the J-integral value. The J-integral solutions obtained in the present article can be applied to calculate the critical J-integral value with taking into account the material non-linearity using experimentally determined critical load at the onset of crack growth from the initial crack tip position in the DCB specimen. The analytical approach applied in the present study is very useful for parametric analyses, since it captures with simple formulae the fundamental non-linear fracture behaviour of the DCB configuration.

References

[1] Bolotin

V.V.

, Delaminations in composite structures: Its origin, buckling, growth and stability, Composites: Part B 27 (1996), 129–145.

[2] Pagano

N.J.

and Schoeppner

G.A.

, Delamination of polymer matrix composites: Problems and assessment, in: Comprehensive Composite Materials, Vol. 2 Kelly

and Zweben

, eds, Elsevier Science, Oxford, 2000, pp. 433–528.

[3] De Morais

and Pereira

, Mixed mode II + III interlaminar fracture of carbon/epoxy laminates, Composites Science and Technology 68 (2008), 2022–2027.

[4] Benzeggagh

M.L.

and Kenane

, Measurement of mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites with mixed-mode bending apparatus, Composites Science and Technology 56 (1996), 439–449.

[5] Carlsson

L.A.

and Pipes

R.B.

, Experimental Characterization of Advanced Composite Materials, Technomic Publishing Co. Inc., 1997.

[6] Kevin

, Interlaminar fracture toughness: The long and winding road to standardization, Composites: Part B 29 (1998), 57–62.

[7] Marat-Mendes

and de Freitas

, Characterization of the Edge Crack Torsion (ECT) test for the measurement of the mode III interlaminar fracture toughness, Engineering Fracture Mechanics 76(18) (2009), 2799–2809.

[8] Szekrenyes

and Uj

, Comparison of some improved solutions for mixed-mode composite delamination coupons, Composite Structures 72 (2006), 321–329.

[9] Liu

, Mechanics and Mechanisms of Fracture: An Introduction, ASM International, 2005.

10.

[10] Gong

X.-J.

and Benzeggagh

, Mixed mode interlaminar fracture toughness of unidirectional glass/epoxy composite, in: Composite Materials: Fatigue and Fracture, Vol. 5 Martin

R.H.

, ed., ASTM STP, Vol. 1230, American Society for Testing and Materials, Philadelphia, 1995, pp. 100–123.

11.

[11] Szekrenyes

and Uj

, Application of elastic foundation beams in composite fracture mechanics – Analysis of the ELS specimen, in: Proceedings of the Fourth Conference on Mechanical Engineering, Vol. I, Springer Hungarica, Budapest, 2004, pp. 191–195.

12.

[12] Hodgkinson

J.M.

, Mechanical Testing of Advanced Fiber Composites, Imperial College, London, UK, 2000.

13.

[13] Suemasu

et al., Double notched split cantilever test method to measure the mixed mode II and III interlaminar toughness, in: 17th International Conference on Composite Materials, Edinburgh, UK, 2009.

14.

[14] Hashemi

Kinloch

and Williams

J.G.

, The effects of geometry, rate and temperature on mode I, mode II and mixed-mode I/II interlaminar fracture toughness of carbon-fibre/poly(ether-ether ketone) composites, Journal of Composite Materials 24 (1990), 918–956.

15.

[15] Ogihara

and Matsuda

, Evaluation of mode III interlaminar fracture toughness of laminated composites, in: 16th International Conference on Composite Materials, Kyoto, Japan, 2007.

16.

[16] Szekrenyes

, Interlaminar fracture analysis in the G_I-G_II-G_III space using prestressed transparent composite beams, Journal of Reinforced Plastics and Composites 30(19) (2011), 1655–1669.

17.

[17] Blackman

B.R.K.

and Kinloch

A.J.

, Fracture tests for structural adhesive joints, in: Fracture Mechanics Testing Methods for Polymers, Adhesives and Composites Pavan

Moore

D.R.

and Williams

J.G.

, eds, Elsevier Science, Oxford, UK, 2001.

18.

[18] Martin

R.H.

, Interlaminar fracture characterization, Key Eng. Mater. 120–121 (1996), 329–346.

19.

[19] Davies

Blackman

B.R.K.

and Brunner

A.J.

, Standard test methods for delamination resistance of composite materials: Current status, Appl. Compos. Mater 5 (1998), 345–364.

20.

[20] Broek

, Elementary Engineering Fracture Mechanics, Springer, 1986.

21.

[21] Hutchinson

J.W.

and Suo

, Mixed mode cracking in layered materials, Advances in Applied Mechanics 29 (1992), 63–191.

22.

[22] Malinin

N.N.

, Applied Theory of Plasticity and Creep, Mashinostroenie, Moscow, 1968.

23.

[23] Hoff

N.J.

, The Analysis of Structures, Wiley, New York, 1956.

24.

[24] Nadai

, Theory of Flow and Fracture of Solids, Vol. 2, McGraw-Hill, New York, 1963.

25.

[25] Seely

F.B.

and Smith

J.O.

, Advanced Mechanics of Materials, Wiley, New York, 1967.

26.

[26] Washizu

, Variational Methods in Elasticity and Plasticity, Pergamon Press, Oxford, 1974.

27.

[27] Korn

and Korn

, Mathematical Handbook for Scientists and Engineers, Nauka, Moscow, 1970.