Tension failure analysis for bolted joints using a closed-form stress solution

Journal Pre-proofs Tension failure analysis for bolted joints using a closed-form stress solution M. Nguyen-Hoang, W. Becker PII: DOI: Reference:

S0263-8223(19)33896-6 https://doi.org/10.1016/j.compstruct.2020.111931 COST 111931

To appear in:

Composite Structures

Received Date: Revised Date: Accepted Date:

17 October 2019 10 January 2020 11 January 2020

Please cite this article as: Nguyen-Hoang, M., Becker, W., Tension failure analysis for bolted joints using a closedform stress solution, Composite Structures (2020), doi: https://doi.org/10.1016/j.compstruct.2020.111931

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2020 Published by Elsevier Ltd.

Tension failure analysis for bolted joints using a closed-form stress solution M. Nguyen-Hoang∗, W. Becker Technische Universität Darmstadt, Fachgebiet Strukturmechanik, Franziska-Braun-Straße 7, 64287 Darmstadt, Germany

Abstract Bolted joints are widely used to connect safety-critical parts in the aeronautical industry and require precise structural assessment tools. Focus of this paper is the development of an overall analytical calculation method to determine the stresses of a bolted joint with quasi-isotropic composite material, which is then used in a tension failure analysis. The stress solution is obtained by means of the Airy stress function. To predict crack initiation in the tension failure analysis use is made of the Theory of Critical Distances. Crack initiation occurs in the net section plane and its assessment requires precisely calculated net section stresses. For this purpose, the stress boundary conditions in load direction at the straight free edges must be fulfilled, which is reached by supplementing the field modelling the load introduction with auxiliary functions. This technique is a physically motivated means to deal with symmetric finite geometry problems involving straight free edges. However stress boundary conditions perpendicular to the load direction are not covered and some inaccuracies in the stress solution may arise. A failure analysis is conducted to investigate their impact on the predicted failure load using literature values based on Finite Fracture Mechanics as reference. Keywords: Bolted joint, Stress analysis, Tension failure, Airy stress function, Theory of Critical Distances

1. Introduction Bolted joints are widely used in the aeronautical industry for reasons of inexpensive manufacturing and the ability to disassemble. Nevertheless drilling a hole introduces a stress-concentrator requiring precise means for stress and failure assessment to ensure both structural integrity and a lightweight optimal design. In the frame of structural analysis five common failure modes need to be considered (Camanho and Matthews, 1997): net section, shear-out, bearing, cleavage and pull-through failure shown in Fig. 1.

Tension/Net section

Shear-out

Pullthrough Bearing

Cleavage

Fig. 1. Failure modes of bolted joints.

Which mode is triggered depends on various factors such as joint geometry, stacking sequence and fibre orientation of the laminate as well as clamping pressure. Corresponding tests ∗ Corresponding

author Email address: [email protected] (M. Nguyen-Hoang) Preprint submitted to Elsevier

were performed by Collings (1977) and Kretsis and Matthews (1985). A review containing mechanical modelling including stress and failure analysis as well as experimental validation is given by Camanho and Lambert (2006). The connection should be designed such that bearing failure involving a gradual failure process, experimentally investigated by Collings (1982) and Xiao and Ishikawa (2005), is likely to be triggered. Contrary, net section failure leads to instantaneous destruction of the joint and must not occur in safety-critical parts. Therefore, precise stress and failure analysis tools are vital for prevention. The present work is dedicated to provide an analytical means to predict the failure load of this mode. Crack initiation and the failure load are influenced by many factors such as width to diameter ratio w/d (Fig. 2), hole size effect and lay-up (Collings, 1977; Hart-Smith, 1980; Kretsis and Matthews, 1985). The stress state is generally three-dimensional. Capturing all failure mechanism factors is complex especially if intending an analytical modelling beneficial in terms of computational effort. To reduce the complexity a two-dimensional in-plane idealisation shall be implemented incorporating the main relevant effects of the failure mechanism. This is applicable if the bolted joint connection contains an overall symmetrical setup with respect to the midplane so that no secondary bending occurs. Particularly the number of plates must be uneven and the plate material behaviour without bending extension coupling. The stress field of the two-dimensional idealisation needs to be determined. Herein, especially the net section stresses in load direction should be precisely modelled since these mainly cause crack initiation. The stress analysis of the pin-loaded hole was analytically performed by many scientists with different level of complexity such as contact modelling January 13, 2020

with or without friction, material anisotropy or finite dimensions. An extensive review is given by Matthews (1987) as well as Camanho and Matthews (1997). Among the solutions there are various idealising the bolt contact by a sinusoidal radial stress distribution along half of the hole boundary. If the material is quasi-isotropic or even isotropic the plane elasticity problem can be solved using the Airy stress function, which was performed by Bickley (1928) for infinite dimensions providing a good approximation of the circumferential stresses at the hole boundary. In Knight (1935) the stress field of both the infinite and the finite width joint with a frictionless load introduction model was obtained. The requirement of vanishing shear stresses at the hole boundary is slightly violated. Furthermore the load transfer of a finite dimensions joint cannot be adequately modelled only considering a finite width while neglecting finite height. For in-plane problems involving an anisotropic material many solutions were developed using the Lekhnitskii complex potential formalism (Lekhnitskii, 1968). The infinite geometry joint without bending extension coupling was treated by De Jong (1976) and further extended for any arbitrary load direction by De Jong and Vuil (1981). A solution for the finite problem was developed by De Jong (1977) using finite width correction factors. Drawbacks of this method are the violated stress boundary conditions at the straight free edges. A similar approach was implemented in Ogonowski (1981) and Echavarría et al. (2007). The former reveals good agreement in the circumferential stresses at the hole boundary in comparison to the finite element (FE) solution for rather large dimensions e/d = 9 and w/d > 4 (ref. Fig. 2) whereas the latter only analyses stress concentration factors at the hole boundary. These may serve for local failure criteria but are insufficient if intending to use those enabling to model the hole size effect requiring the net section stress decay. Kratochvil and Becker (2010) treated a finite dimensions bolted joint by approximating the rectangular plate geometry as an ellipse showing good correlation to FE in the circumferential stresses for [±45°] s laminates and w/d = 3 as well as for [0°/90°]s laminates and w/d = 5. Bending extension coupling was covered by Grüber et al. (2007) revealing good agreement to FE solution for w/d = {15, 45}, a ratio in which finite width effects are slight or even have decayed and tension failure is not expected. The calculus is extended to multilayered composites bolted joints with interference-fit under thermal load (Grüber et al., 2018). The solutions summarised above reveal agreement for rather large dimensions if validated using FE at all. To the authors’ knowledge for w/d ≤ 4 there exists no analytical approach to accurately calculate the stress field especially the net section stresses essential for tension failure assessment, which is not only validated using FE but also further established in a subsequent failure analysis. These gaps shall be treated by the present paper, in which the stress field of a bolted joint with finite dimensions is determined using the Airy stress function. The calculus is dedicated to quasi-isotropic laminates and can be extended to the general anisotropic case once the methodology is proven to be feasible. Based on the present stress solution a failure analysis of the net section mode likely occurring for small widths w/d ≤ 4 (Hart-

Smith, 1980) is conducted. In a first step the stress function describing the load introduction is obtained while neglecting stress free boundary conditions at the straight edges. For a full load transfer through the net section plane at least the stress free boundary conditions in load direction need to be fulfilled. This is achieved by superimposing virtual auxiliary functions created by mirroring the load introduction field enabling full elimination of non-zero tractions in load direction at the straight edges while maintaining the stress boundary conditions at the hole. Those tractions perpendicular to the load direction at the straight edges do not play a role in the equilibrium in load direction and are thus not covered by the present calculus. To the authors’ knowledge this approach is used for the first time to deal with symmetric finite domain problems containing straight free edges. The resulting net section stresses will transfer the whole external load but nevertheless may be differently shaped in comparison to the FE solution obeying all stress boundary conditions. The criticality of those potential stress deviations shall be further quantified by conducting a failure analysis using the Theory of Critical Distances (TCD, Whitney and Nuismer (1974); Taylor (2007)). The TCD postulates crack initiation if the net section stress in load direction equals the plain strength of the material at or averaged over a certain characteristic hole distance, which is assumed to be invariant with respect to the hole diameter. The resulting failure stresses are compared to the experiment and to the predictions using Finite Fracture Mechanics (FFM). In the frame of FFM (Weißgraeber et al., 2015) a crack of the a priori unknown finite length ∆a instantly initiates if both a stress and energy criterion are fulfilled. These criteria yield ∆a and thus the location where the stress criterion is evaluated being dependant on the defect size neglected by the TCD. Hence investigating to which extent the use of a characteristic distance in the frame of the TCD is justifiable shall be another focus of the present paper. 2. Determination of the stress field The stress field for a quasi-isotropic pinned-hole problem with finite dimensions w/d and e/d as shown in Fig. 2 is determined. Herein the quantity Py denotes an external force per plate thickness applied in the vertical direction whereas d represents both bolt and hole diameter. This setting is idealised as a boundary value problem, whose stresses are determined by means of the Airy stress function. Fundamentals of the Airy stress function Let us assume a two-dimensional plane stress or plane strain problem with linear elastic and isotropic material behaviour. Then the governing equations equilibrium, compatibility and Hooke’s law can be reduced to one single differential equation involving the Airy stress function F. For a plane stress problem with vanishing body forces the function F needs to satisfy

2

Py

y

(Waszczak and Cruse, 1971). The complete stress boundary conditions read  2F y    sin(ϕ) for 0 ≤ ϕ ≤ π, −  π R σr (R, ϕ) =      (4) 0 for π ≤ ϕ ≤ 2π,

e

r ϕ x

∅d

τrϕ (R, ϕ) = 0

for 0 ≤ ϕ ≤ 2π,

where R = d/2 denotes the hole radius. The hole boundary conditions shall be fulfilled by FLI representing the load introduction but also by the final finite bolted joint solution F. w

Stress function FLI modelling the load introduction The field FLI obeying the tractions along the hole boundary in Eq. (4) is developed. The stress-free straight edges are neglected at this stage. As in Knight (1935) two partial functions are constructed by (5) FLI = F1 + F2 .

Fig. 2. Geometry of bolted joint.

the following biharmonic equation (Timoshenko and Goodier, 1951; Sadd, 2005): ∆∆F = 0

with

 2  ∂ F ∂2 F    + 2    ∂y  ∂x2 ∆=  2  ∂ 1 ∂2 ∂ 1      ∂r2 + r ∂r + r2 ∂ϕ2

The first partial function represents a full sine in the radial stresses over the complete hole boundary whereas the second a cosine Fourier series expansion of |σr1 (R, ϕ)| belonging to F1 . The radial stresses σr1 (R, ϕ) and σr2 (R, ϕ) have the same behaviour at the upper circular boundary whereas at the lower they have a reversed sign and therefore superimposed lead to cancellation. This is illustrated in Fig. 3. Contrary to Knight (1935) both partial functions fulfil vanishing shear stresses along r = R.

in cartesian coordinates, (1) in polar coordinates.

Functions F fulfilling the biharmonic equation are called biharmonic functions and were addressed by Michell (1899). The plane stress components are derived using σx =

∂2 F , ∂y2

σr =

1 ∂2 F 1 ∂F ∂2 F + 2 2 , σϕ = 2 , r ∂r r ∂ϕ ∂r

σy =

∂2 F , ∂x2

∂2 F , (2) ∂x∂y   ∂ 1 ∂F =− . (3) ∂r r ∂ϕ

+

τ xy = − τrϕ

=

F1 full sine

Stress functions F are now to be chosen in such a way that the corresponding stresses fulfil the given stress boundary conditions.

F2 cosine Fourier series of magnitude of full sine

FLI contact idealisation

Fig. 3. Radial stresses at the hole boundary. The tractions of F1 and FLI yield a non-zero force resultant contrary to F2 .

According to Bickley (1928) the stress function F1 and the corresponding stresses calculated using Eq. (3) read

Overview of the calculus At first a bolt contact idealisation is performed providing the hole boundary stress conditions to satisfy. A stress function FLI fulfilling them is then developed and represents the load introduction, whose force flux is distributed in all directions reaching stress free conditions at infinity and therefore describes the pinned hole with infinite geometry. In the finite dimensions case, however, the force flux is directed towards the clamp (ref. Fig. 2) through the net section plane. This load transfer is modelled by supplementing FLI with auxiliary functions eliminating the stresses in load direction at the straight edges while keeping the hole boundary conditions unchanged.

1 r sin ϕ + b13 sin ϕ, R r   1 1 1 σr1 = −2a15 + b12 − 2b13 3 sin ϕ, r r r   1 1 σϕ1 = b12 + 2b13 3 sin ϕ, r r   1 1 τrϕ1 = −b12 + 2b13 3 cos ϕ. r r F1 =a15 rϕ cos ϕ + b12 r ln

Bolt contact idealisation The bolt plate contact problem is idealised by a sinusoidal function in the radial direction and with the assumption of no friction and thus vanishing shear stresses at the hole boundary

(6) (7) (8) (9)

The involved coefficients are determined by b12 1 = (1 − ν), a15 2 3

(10)

ensuring single-valued circumferential displacements uϕ (Bickley, 1928; Timoshenko and Goodier, 1951) with ν being the Poisson’s ratio of the plate material, σr1 (R, π/2) = −

1 Fy π R

τrϕ2 = − 2

σr2 (R, ϕ)/σ0 =

with

r

−2

2π 1 1 σr1 (R, ϕ) cos nϕ dϕ π σ0 0   4 1 + cos nπ    for even n, − π2 1 − n2 =     0 for odd n.

n(2n + 1)A2,n

 R 2n+2

 R 2n  r

(19)

f2,n =

f2,0 , 2

−2[n(2n + 1)A2,n + (n + 1)(2n − 1)B2,n ] = f2,2n .

(20)

(21) (22)

Further ensuring τrϕ2 (R, ϕ) = 0 in Eq. (18) results in (2n + 1)A2,n + (2n − 1)B2,n = 0.

n=1

+ (n + 1)(2n − 1)B2,n

r

sin 2nϕ σ0 . (18)

f2,0 f2,n cos nϕ, + 2 n=1

b2 =

r

N    R 2  R 2n+2

σϕ2 = − b2 n(2n + 1)A2,n +2 r r n=1  R 2n  + (n − 1)(2n − 1)B2,n cos 2nϕ σ0 , r

2n − 1 B2,n , 2n + 1 f2,2n 1 . = 2 n(2n − 1) − (n + 1)(2n − 1)

A2,n = − B2,n

(24)

Note that F2 does not depend on ν and the corresponding tractions at the hole boundary incorporate a vanishing force resultant. Also refer to Fig. 3 for illustration. The stress field FLI is now completely determined. In order to eliminate its remaining tractions in load direction at the straight edges and thus leading the load transfer through the net section plane towards the clamp auxiliary functions are superimposed. Finite geometry modelling using auxiliary stress functions Let us investigate the free body diagram in Fig. 4 giving an idea about the stresses of the final solution. For the transfer of the whole external load through the net section plane its normal stresses σy (x, 0) must equilibrate Fy alone, which requires vanishing tractions in load direction at the straight free edges. This is reached by superimposing virtual neighbouring auxiliary plates created by mirroring the load introduction field 11 FLI = FLI and arranging it as illustrated in Fig. 5. For elimina11 tion of the normal stresses σ−y (x, e) at the horizontal edge FLI needs to be vertically copied, its loading reversed to tension and rotated by π along the hole boundary creating the auxiliary

+

cos 2nϕ σ0 ,

(23)

Using Eq. (22) and (23) the Airy stress coefficients can be calculated by

N   R 2n  R 2n−2  r  F2 = R2 b2 ln + + B2,n A2,n cos 2nϕ σ0 , R n=1 r r (15) with the stress components

b2

 R 2n 

and equate the coefficients of Eq. (16) and (19) while taking into account that uneven f2,n are zero yielding

being multi-valued due to the linear term in ϕ and thus modelling a non-zero traction. However if the force resultant of the tractions along each of all boundaries is zero as in doubly symmetric problems then the stress function is independent of the elastic constants (Michell, 1899) such as the circular hole under uniform tension in an infinite domain (Kirsch, 1898; Timoshenko and Goodier, 1951; Sadd, 2005). The magnitude of the radial stresses σr1 (R, ϕ) shall now be expanded by a Fourier cosine series. The corresponding stress function has the general form

σr2 =

+

(14)

=0

N 

r

N∗

(12)

 σ 0 = sin ϕ cos ϕ − ϕ π

 R 2

 R 2n+2

Let us represent |σr1 (R, ϕ)| by the cosine Fourier series

leading to vanishing shear stresses at the hole boundary. NorF malising with respect to the reference stress σ0 = dy , F1 then is given by r 1 1 1 R rϕ cos ϕ+ (1−ν) r ln sin ϕ+ (1−ν) R2 sin ϕ σ0 . F1 = π 2 R 4 r (13) Note that the stress function and thus the stress state depend on ν. This is true for any boundary value problem containing tractions with non-zero force resultant (Bickley, 1928) such as the current problem or the infinite plate loaded by a single force (Timoshenko and Goodier, 1951). The corresponding stress function must contain first order terms of sin ϕ and/or cos ϕ so that at least one of the loading functions is multi-valued. This is exemplarily shown for F1 in the current bolted joint problem with

ϕ py1 (ϕ) = σr1 (R, ϕ) ˆ sin ϕˆ + τrϕ1 (R, ϕ) ˆ cos ϕˆ dϕˆ



n=1

(11)

τrϕ1 (R, ϕ) = 0,

n(2n + 1)A2,n

+ n(2n − 1)B2,n

transferring half of the external load and eventually

0

N  

(16)

(17) 4

σy (x, e) = 0

τ xy (−w/2, y) ≈ 0

with

τ xy (w/2, y) ≈ 0

[x j ] = [ x1 x2 x3 x4 ...] = [x x − w x + w x − 2w ...] ,

(26)

[yi ] = [y1

(27)

y2 ] = [y y − 2e] .

The superimposed field of the load introduction is obtained by sup

FLI =

σr (R, ϕ)

y

n nx y =2 k =2 n

Fki j ,

(28)

k=1 i=1 j=1

 sup

=Fk

sup

x

with Fk being a superimposed partial field. The quantity n x represents the total number of horizontally and ny the total number of vertically aligned plates. For the specific bolted joint problem n x is regarded as sufficiently large if there is less than 1 % inaccuracy in the load transfer due to remaining shear stresses at the vertical edges in the superimposed load introsup duction field FLI . This is expressed by

σy (x, 0) Fig. 4. Free form body. Only stresses relevant for equilibrium in load direction are shown.

T ... τ−xy

τ+xy

T τ−xy

τ+xy

23 FLI

τ−xy

21 FLI

22 FLI

σ+y

σ+y

σ−y

σ−y

σ−y

τ−xy

C τ+xy

13 FLI

τ−xy

 ny =2 nx,min  2  e  i j j  τ xy1 (±w/2, y) + τixy2 (±w/2, y) dy  ≤ 0.01, Fy 0 i j (29) and the smallest n x obeying Eq. (29) eventually leads to n x,min . Note that this mirror technique may be applied to determine the stress field of other plane problems involving free straight edges in a finite domain if the setup is symmetric. This would only require to modify the Fourier series of the load introduction field being capable to model any traction expressible by a continuously differentiable function. Also note that the approach to scale the net section stress decay of the infinite problem such that its integrated stresses carry the external load over the width of the finite problem as performed for open-holes by Tan (1988) fails. In doing so the terms of the Fourier series representing the net section stresses σy,LI (x, 0) of the load introduction field are equally scaled. This approach coincidentally yields good results for open-holes involving only a few terms to model the infinite stress field. Enhancing the approach by scaling them differently is rather hard to perform since there is no clue how these shall be scaled without knowing the solution a priori. The mirror technique however provides a systematic and physically motivated means for load transfer modelling no matter how complex the load introduction is.

τ+xy ...

σ+y C ... τ+xy

T

C τ+xy

11 FLI

τ−xy ... 12 FLI

11 . C = Fig. 5. Schematic how to arrange virtual auxiliary plates of FLI = FLI ˆ Compression, T =ˆ Tension.

21 plate FLI . For vanishing shear stresses τ∓xy (±w/2, y) at the two 11 12 straight vertical edges FLI is horizontally mirrored creating FLI 13 and FLI . Since the left auxiliary plate does not only cancel the shear stresses at x = −w/2 but also affect those at the opposite edge x = w/2 and vice versa for the right auxiliary plate and its opposite edge at x = −w/2 there will be a remaining shear traction at both vertical edges. Its magnitude depends on the number of horizontally aligned auxiliary plates in use. If a periodic row of auxiliary plates is superimposed the shear tractions will vanish. Note that this shear stress cancellation approach is limited to symmetric problems. As the horizontal auxiliary plates will induce further normal stresses σ−y at y = e they need 11 to be vertically copied as FLI ensuring full cancellation. The auxiliary plates for the partial fields F1 and F2 of the load introduction field FLI are calculated using

  Fki j (x, y) = (−1)i+1 · Fk (−1)i+1 x j , (−1)i+1 yi , Tension/ Compression

Correction functions to mitigate deviations in the hole boundary conditions Beyond the desired effect of stress elimination at the free edges the auxiliary plates also affect the hole boundary conditions fulfilled by the unmirrored load introduction field FLI . The arising deviations in the tractions along the hole boundary are mitigated by superimposing a correction field F3 expanding the deviations in a Fourier series. Since this additionally introduced field disturbs the zero stress boundary conditions in the load direction at the free straight edges it needs to be mirrored in the same way as FLI , which again yields slightly violated hole boundary conditions. Iterating this procedure by superimposing further correction functions Fk>3 and mirroring them

(25)

load shifting along hole

5

Input parameters w/d, e/d, ν

leads to vanishing deviations in the hole boundary conditions. The full solution F that correctly models the load introduction as well as fulfils the boundary conditions in load direction at the straight free edges with the smallest computational effort is eventually assembled by

Assemble load introduction field FLI

n y =2 n x,min n k,min

Determine n x,min using

F=

i

j

Fki j .

(30)

k

LI χLI τ xy ≤ χτ xy,min = 0.01

The task of each partial field Fk is listed in Tab. 1. Table 1 Stress fields assembling the full solution.

stress field F1 F2 FLI = F1 + F2 Fk≥3

task at hole boundary represents full sine in radial stresses Fourier series expansion of |σr1 (R, ϕ)| sinusoidal bolt load introduction sup correction of deviating stresses in Fk−1

Mirror and superimpose auxiliary field

The minimum number of required correction functions nk,min is regarded as sufficiently large if the relations 

 R  π  σr (R, ϕ) sin ϕ dϕ  = 1, χσr = Fy  0

2π  R σr (R, ϕ) 2 sin ϕ dϕ ≤ 0.01, ψσr = (31) Fy π

2π  R ψτrϕ = τrϕ (R, ϕ) 2 cos ϕ dϕ ≤ 0.01, Fy 0

Hole BCs violated? • χ σr =

2

+

N∗

+

τ

τrϕk (R, ϕ)/σ0 =

fk,0rϕ 2

gσk,nr sin nϕ, (33)

n=1

N∗

N∗

N∗

rϕ gk,n sin nϕ, (34)

n=1

σr fk,n cos nϕ +

Full load transfer through net section? • 0.98 ≤ χσy ≤ 1 ?

τ

No

Increase n x,min and rebuild correction functions Fk≥3

Yes Full stress field Fig. 6. Flowchart of the bolted joint calculation method.

Table 2 Deviating stresses at the hole boundary to be expanded.

correction σdev rk (R, ϕ) field  sup sup  F3 σr1 + σr2 − σr1 − σr2 r=R  sup  F4 σr3 − σr3 r=R

n=1

fk,nrϕ cos nϕ +

Introduce additional field Fk≥3 to correct hole BCs

No

A too small χσy means that the shear stresses elimination at the vertical edges is insufficient which can be cured by increasing the number of horizontal auxiliary plates n x,min . Note that after its increase the iterative correction function routine must be recalculated since the deviations at the hole boundary depend on the number of auxiliary plates. The flowchart in Fig. 6 summarises the overall calculus. Let us proceed with the determination of the correction functions Fk≥3 . Their radial and shear stresses at the hole edge expand the deviations caused by the previously superimposed auxiliary fields enabling their mitigation. The expansions for the specific bolted joint setting symmetric to the y-axis read σr fk,0

Yes

• ψτrϕ ≤ 0.01 ?

with the load transfer ratios based on those stresses of the full superimposed solution F are true. A load transfer value χσr = 1 indicates that the external load is fully transferred by the radial sinusoidal stresses at the hole boundary. The quantities ψσr and ψτrϕ are means to assess undesirable non-zero oscillations at r = R. Inaccurate load transfer values concerning the hole tractions indicate that the number of correction functions nk,min shall be increased. A final check of the force flux through the net section area is done by quantifying the load transfer ratio

w/2 2 χσy = σy (x, 0) dx with 0.98 ≤ χσy ≤ 1. (32) Fy R

σrk (R, ϕ)/σ0 =

1 ?

• ψσr ≤ 0.01 ?

τ

n=1

6

τdev rϕk (R, ϕ) sup sup  − τrϕ1 − τrϕ2 r=R  sup  τrϕ3 − τrϕ3 r=R



taking into account that for uneven n the coefficients

2π 1 1 σr = σdev fk,n rk (R, ϕ) cos nϕ dϕ, π σ0 0

2π 1 1 σdev gσk,nr = rk (R, ϕ) sin nϕ dϕ, π σ0 0

2π 1 1 τ τdev fk,nrϕ = rϕk (R, ϕ) cos nϕ dϕ, π σ0 0

2π 1 1 τrϕ = τdev gk,n rϕk (R, ϕ) sin nϕ dϕ. π σ0 0

 R 3  R τ cos ϕ − τrϕk = − dkσr cos ϕ + 2dk rϕ r r N   R 2n+2

n(2n + 1)Ak,n + −2 r n=1  R 2n  + n(2n − 1)Bk,n sin 2nϕ + r N   R 2n+3

+ 2(2n + 1) + (n + 1) Ck,n r n=1  R 2n+1  + n Dk,n cos(2n + 1) ϕ σ0 , r

(35)

vanish. For the geometrical configurations investigated two iterations are sufficient to satisfy the stress boundary conditions at the hole and those in load direction at the vertical edges. Tab. 2 shows the particular stress deviations to be expanded by the correction functions Fk≥3 . The general form for an Airy stress field Fk describing those Fourier expansions in Eq. (33), (34) is expressed by  r r Fk (r, ϕ) = R2 bσk r ln + cσk r ϕ cos ϕ+ R R r r τ R + dkσr ln sin ϕ + dk rϕ sin ϕ+ R R r N   R 2n  R 2n−2 

+ + Bk,n Ak,n cos 2nϕ+ r r n=1

(39)

Equating the coefficients of the Fourier series in Eq. (33), (34) with those of the Airy stress components in Eq. (37), (39) while taking into account single valued displacements (Bickley, 1928; Timoshenko and Goodier, 1951) ensured by dkσr cσk r

=

1 (1 − ν), 2

(40)

yields bσk r =

N   R 2n+1  R 2n−1 

Ck,n sin(2n + 1) ϕ σ0 , + Dk,n r r n=1

1 σr f , 2 k,0

cσk r = − τ

dk rϕ =

(36)

Ak,n

with the stress components R  R 3   R 2  τ − 2cσk r − dkσr sin ϕ − 2dk rϕ sin ϕ − σrk = bσk r r r r N   R 2n+2

+ −2 n(2n + 1) Ak,n r n=1  R 2n  + (n + 1)(2n − 1) Bk,n cos 2nϕ − r N   R 2n+3

(2n + 1)(n + 1) Ck,n + −2 r n=1  R 2n+1  sin(2n + 1) ϕ σ0 , + n(2n + 3) Dk,n r (37)  R 2   R R 3 τ + dkσr sin ϕ + 2dk rϕ sin ϕ + σϕk = − bσk r r r r N   R 2n+2

n(2n + 1) Ak,n + +2 r n=1  R 2n  (38) cos 2nϕ + + (n − 1)(2n − 1) Bk,n r N   R 2n+3

+ +2 (2n + 1)(n + 1) Ck,n r n=1  R 2n+1  + n(2n − 1) Dk,n sin(2n + 1) ϕ σ0 , r

 1  τrϕ fk,1 + gσk,1r , 2  1  τrϕ fk,1 + dkσr , 2

dkσr =

  1 τrϕ 1 gk,2n + n(2n − 1) Bk,n , =− n(2n + 1) 2

 1 1 − ν cσk r , 2 (41)

 τrϕ  τrϕ σr  1 gk,2n − fk,2n 1  fk,2n+1 , Ck,n = − nDk,n  ,  2 2n − 1 n + 1 2(2n + 1)  1  σr τrϕ + fk,2n+1 . g Dk,n = − 4n k,2n+1 (42) With that all partial functions to build the full solution are determined. Note that the calculus can be extended to the more general case of orthotropic composite laminates using Lekhnitskii complex potential formalism (Lekhnitskii, 1968). Then Hooke’s law represents the plane stress state of an orthotropic material model with effective stiffnesses of the laminate and the governing equations as well as the stress field can be expressed by the complex potentials. The approach to implement finite dimensions by iteratively mirroring and correcting the hole boundary condition stays the same as introduced for quasi-isotropic laminates by means of the Airy stress function. Those results for different geometry ratios w/d and e/d are presented and verified using Finite Element analysis in the next section. Bk,n =

3. Discussion of the stress results The results for the joint configurations with the geometrical dimensions in the range w/d = {3, 20} and e/d = {3, 10} are presented and verified using a Finite Element model implemented in A BAQUS. The mesh consists of CPS8 continuum 7

e

150

−→

w

ξ net section path

Py

∅d

n x,min [−]

100

Fig. 7. Finite Element model for w/d = 3, e/d = 3.

e/d = {3, 5, 10}

50

0 3 4

10 w/d [−]

15

20

−→

Fig. 9. Minimum number of horizontal plates n x,min .

(a) Unmirrored Load introduction.

tios are calculated. For the different bolted joint configurations the present calculus leads to load transfer ratios within the range

(b) Present finite geometry solution.

χσr = 1,

ψτrϕ < 0.005, (45) fulfilling the requirements defined in Eq. (31), (32) and therefore confirming the methodology. Fig. 9 shows n x,min with respect to the geometric properties leading to these load transfer value ranges. The minimum number of horizontal plates n x,min increases with a smaller ratio w/d since the shear stress decay is more limited than in wider joints requiring additional horizontal auxiliary plates for cancellation of τ xy (±w/2, y). Furthermore n x,min increases with the quantity e/d. A small relative hole distance e/d requires the auxiliary plates F2 j being located near the original vertical boundaries at x = ±w/2 ∧ y ≤ e. The shear stresses τ2xyj (±w/2, y) of F 2 j have a reversed sign in comparison to the field F 1 j and have not decayed as much as in joints with larger e/d. Therefore they contribute more significantly to the shear stress elimination by the horizontally aligned auxiliary plates F 1 j . Hence less horizontal plates are required to obtain a given shear stress cancellation ratio χτxy . For further illustration refer to Fig. 5. Let us investigate the circumferential as well as the net section stresses relevant for the crack initiation assessment plotted in Fig. 10. Generally the wider the geometry properties the better the agreement in both stress components is observed. Peaks in σϕ (R, ϕ) at ϕ = {0, π/2, π, 3/2 π} as well σy (x, 0) deviate the more the smaller the dimensions are chosen. Deviations are due to non-zero stress boundary conditions perpendicular to the load direction at the straight edges not covered by the present calculus, particularly τ xy (x, e)  0 and σ x (±w/2, y)  0. This is further confirmed by extracting these non-zero tractions and applying them additionally to the sinusoidal loading in the FE model, which is performed for all configurations except for w/d = 20, e/d = 10 yielding any pronounced stress deviations. The corresponding stress solution then coincides with the results of the present calculus apart from very slight differences (ref. to Fig. 10). They arise since a load transfer value

(c) FE solution. Fig. 8 Force flux of stress vectors ty for bolted joint with w/d = 20, e/d = 10, n x,min = 23.

plane stress elements with 8 nodes. A sinusoidal curve in the radial tractions modelling the bolt contact idealisation is applied at r = R. Convergence in the stresses is reached for 72 elements along the hole boundary. Fig. 7 presents the FE model for w/d = 3, e/d = 3 and introduces the dimensionless coordinate x−R , (43) ξ= w/2 − R

with ξ = 0 at the hole boundary and ξ = 1 at the free straight edge. Extensive Finite Element analyses concerning the pinloaded hole containing a contact model as well as a sinusoidal bolt contact idealisation can be found in Crews et al. (1981). To qualitatively investigate the impact of the virtual auxiliary plates on the load transfer let us exemplarily plot the force flux using the stress vector       0 τ σ τ xy ty = x · (44) = − xy , τ xy σy σy −1 for w/d = 20, e/d = 10 shown in Fig. 8. If superimposing the auxiliary fields the force flux becomes tangent to the straight edges as in the FE solution and thus providing a physical load transfer. For further quantitative approval the load transfer ra-

8

0.98 ≤ χσy ≤ 0.99,

ψσr < 0.005,

2.5

−→

rc,FE = 1.143 mm

FE

1,000

XTL = 845.1 MPa

Airy

σy (ξ, 0)/σ0 [−]

σy (r˜c , 0) [MPa]

−→

rc,Airy = 0.351 mm

500

2 ξrc,FE = 0.381 1.5 1

Airy ξrc,Airy = 0.117

0.5 FE

0

0

0.5

1 r˜c [mm]

1.5 −→

0

2

Fig. 11. Calibration of characterisitc distances using averaged net section stresses at failure of experiment NT2.

w/d [-] 1.75 2

e/d [-] 5.83 4.17

0.4 ξ [−]

0.6 −→

0.8

1

Comparison to experimental data A failure analysis is conducted with the derived characteristic distances. The failure stress for the test setting NT3 is predicted using both analytically and numerically determined stresses. It is obtained by scaling the external load such that the corresponding average net section stresses reach XTL at the characteristic distances. The normalised net section stresses can be found in Fig. 12. The dimensionless characteristic distances are derived by rc . (49) ξrc = w/2 − R Tab. 4 compares the predicted stresses at failure using their deviations normalised to experimental data (Exp). Table 4 Failure stresses for NT3.

Table 3 Test data from Catalanotti and Camanho (2013).

d [mm] 6 6

0.2

Fig. 12. Normalised net section stresses for experiment NT3.

Determination of the characteristic distance Two different methods are introduced to determine the characteristic distance. Let us begin with the approach by Whitney and Nuismer (1974) in which rc is calibrated using the stress field and the experimentally determined failure load for a specific bolted joint configuration. It is thus a non-physical size and generally not applicable to other joint configurations but possibly to a certain extent. For the present study the experiment NT2 performed in Catalanotti and Camanho (2013) is chosen to calibrate rc . This test setting contains the composite plate material Hexcel IM7-8552 with the nominal thickness is t = 3 mm and the layup [90°/0°/ ± 45°]3s providing quasi-isotropic material behaviour. Further test specifications are listed in Tab. 3.

Test ID NT2 NT3

0

TCD-CLB

σF [MPa] 466.2 526.7

TCD-Taylor

FFM Exp Airy FE Airy FE σF [MPa] 540.6 569.1 735.4 641.5 549.5 526.7

[-] 2.6 % 8.1 % 40.6 % 21.8 % 4.3 % -

The quantity σF denotes the failure stress with σF = Py,F /d where Py,F is the external load at failure. The stress field at failure is derived using both FE software and the present calculus. The characteristic distances rc,Airy and rc,FE are those distances where the corresponding average stresses σy (r˜c ) equal XTL = 845.1 MPa. This is further demonstrated in Fig. 11. Large deviations arise since the geometry parameter w/d = 1.75 is very small and non-zero tractions perpendicular to the load direction uncovered by the present calculus strongly affect the net section stresses. In the context of the Theory of Critical Distances by Taylor (2007) (TCD-T) the characteristic distance is derived by strength values of the plain material,  2 2 Kc rc,T = = 1.632 mm, (48) π XTL √ where Kc = 42.8 MPa m (Catalanotti and Camanho, 2013) denotes the fracture toughness.

A good agreement is observed if using the TCD with calibrated characteristic distance involving stresses of the present solution (TCD-CLB-Airy). Since the corresponding failure stress prediction yields a deviation even smaller than that based on FE stresses the calibration approach can be considered as capable to cope with erroneous net section stresses of the present solution for the particular joint configuration. Further experiments with different values for d and w/d should be performed to assess the present method’s capabilities. Results deduced by means of the Finite Fracture Mechanics (FFM) performed in Catalanotti and Camanho (2013) are also shown. The TCDCLB-Airy approach coincidentally leads to a prediction better than the FFM although the latter unifies both stress and energy criterion contrary to the TCD and is therefore the more sophisticated crack initiation model. In the context of the TCD by Taylor poor results are observed. Let us now investigate the predictions for different diameters and analyse how well the hole size effect is modelled by the TCD approaches. 10

5. Failure analysis: modelling of the hole size effect

enabling the calculation of the dimensionless crack length. Back substitution eventually yields the stress at failure.

The failure stress is predicted using both TCD approaches for joint configurations with the same ratios w/d, e/d as in NT3 but with different bolt diameter d. The failure load is determined as in the previous section. Note that the relative characteristic distance ξrc where the stress failure criterion is evaluated decreases the greater d and vice versa. Since there are no experimental failure loads with respect to different diameters accessible to the authors those failure stresses based on the FFM published in Catalanotti and Camanho (2013) are reproduced and used as reference. The FFM has been applied and experimentally confirmed to various problems. Among of them are some concerning an uniaxially loaded open-hole such as Li and Zhang (2006) investigating infinite geometry and isotropic material or Hebel and Becker (2008) for composite with anisotropic behaviour pursuing a numerical approach contrary to Camanho et al. (2012) implementing a closed-form failure model for a finite width geometry. Furthermore the open-hole under combined tension and in-plane bending is numerically solved by Rosendahl et al. (2017). These examples emphasize that one can regard the FFM as an established concept to predict crack initiation in brittle materials. In the context of the FFM (Weißgraeber et al., 2015) the initiation of a crack with the a priori unknown finite length ∆a is assumed which instantaneously and unstably propagates if both stress and energy criteria are fulfilled. This condition is called coupled criterion established by Leguillon (2002) and generally yields an optimisation problem to determine the minimal load and the corresponding crack length yielding its initiation. The bolted joint problem is characterised by a monotonic decrease of the stresses and a monotonic increase of the energy release rate with respect to ∆a. Using the dimensionless crack length ∆ξ = ∆a/(w/2 − R) the coupled criterion then can be simplified to the conditions

∆ξ 1 σy (ξ, 0) dξ = XLT ∆ξ 0 (50)

∆ξ

∆ξ ∧

0

KI2 (ξ) dξ =

0

6. Discussion of the failure analysis results A comparison between the FFM and TCD failure stresses with respect to the hole diameter d is shown in Fig. 13. Be reminded that the abbreviations TCD-CLB-Airy, TCD-CLB-FE refer to the approaches in which the predictions are derived using the TCD with the calibrated characteristic distance rc and stresses based on the present solution or FE respectively. Furthermore TCD-T-FE refers to the TCD by Taylor in which the characteristic distance rc,T and FE stresses are used. The bearing cut-off is defined as σF = 700 MPa similar to Catalanotti and Camanho (2013). Beyond this limit bearing failure occurs and corresponding criteria instead those of the net section crack initiation have to be used. Note that setting the bearing cut-off involves a certain arbitrariness since there exist various definition criteria such as a specific elastic or plastic hole deformation or a particular degradation level of the laminate stiffness noticeable by a non-linearity in the load-strain curve from which one may freely choose. Furthermore the deviations normalised to the FFM prediction are shown in Fig. 14 and 15. A deviation limit of | | ≤ 10 % is regarded as tolerable giving us the diameter range in which the TCD approaches can be used. This is further visualised in Fig. 16. The TCD-T-FE approach leads to acceptable results for d ≥ {12 mm, 13 mm} with w/d = {2, 3} whereas the TCD-T-Airy concept generally yields very poor results. Its deviations lie within the range 15 % ≤ ≤ 40 % and corresponding results are not shown. The TCD-CLB-FE prediction is almost coinciding with the FFM reference beginning with d ≥ 17 mm for w/d = 2 and concerning w/d = 3 in the whole failure stress range below the bearing cut-off. Predictions by the TCD-CLB-Airy fulfil the deviation limit for rather small diameter 1 mm ≤ d ≤ 15 mm in w/d = 2 and can be used within a range where the TCD-CLB-FE approach fails. Beyond this diameter range the TCD-CLB-Airy concept lacks in precision. Contrary for w/d = 3 the TCD-CLB-Airy can be used within 8 mm ≤ d ≤ 354 mm. All predictions using FE stresses involve deviations converging with increasing d to a value smaller than the limit and are thus also applicable for d > 50 mm beyond the range investigated. The TCD-CLB approaches generally lead to smaller deviations for a wider diameter range in comparison to their TCD-T counterpart. Let us further analyse the methodology of the approaches to understand why they lead to different predictions. The TCD is purely based on a stress criterion which is evaluated at a certain characteristic hole distance assumed to be invariant to the hole diameter. Contrary as shown in Fig. 17 for different w/d configurations the FFM additionally incorporating the energy criterion reveals the crack length ∆a and therefore the location where the stress criterion is evaluated in fact being dependant of the hole diameter d. Above a certain hole diameter the curves ∆a(d) take a shape which can be approximated by a constant plateau value. In here the TCD approaches using FE stresses yield accurate predictions if the corresponding characteristic distances lie near or within the plateau. This is true for the

2 KIC (ξ) dξ

where KI is the mode I stress intensity factor of a newly initiated crack. For the R-curve KIC (ξ) a Gompertz function is used modelling a crack length dependency of the fracture toughness. This is essential if the crack length is of the same order as the process zone lpz which particularly occurs for small hole diameters. The net section stresses as well as the stress intensity factor are calculated using FE for joints with different geometry ratios w/d and the unit load σ0 = 1. The results are then approximated using the polynomial fitting functions φ modelling the net section stresses and ψ the energy release rate. Their coefficients can be taken from Catalanotti and Camanho (2013). Scaling these fitting functions with the unknown failure load σF,FFM , inserting then in the coupled criterion Eq. (50) and eliminating σF,FFM leads to $ % % & ' ∆ξ

∆ξ ψ(ξ, w/d)2 dξ 1 0 (51) φ(ξ, w/d) dξ = ∆ξ ' ∆ξ 2 XLT 0 KIC (ξ) dξ 0

11

bearing cut-off

700

TCD-CLB-Airy

σF [MPa]

w/d = {2, 3}

500

TCD-T-FE TCD-CLB-FE

0%

TCD-CLB-Airy FFM

bearing cut-off FFM

−10 %

TCD-T-F

E

400 1

10

20 d [mm]

30 −→

40

10 %

20 d [mm]

TCD-T-FE

30 −→

40

50

11-50 mm w/d = 3

TCD-CLB-Airy

−→

10

Fig. 14. Deviations of predicted failure stresses to FFM for w/d = 2

deviation limit

TCD-T-FE TCD-CLB-FE

0%

deviation limit

1

50

Fig. 13. Predicted bearing stresses at failure

[−]

deviation limit

−→

TCD-CLB-FE

600

[−]

−→

10 %

TCD-CLB-Airy

8-50 mm

TCD-CLB-FE

7-50 mm 12-50 mm

TCD-T-FE bearing cut-off FFM

−10 %

deviation limit

1

10

20 d [mm]

30 −→

40

TCD-CLB-Airy

1

rc,FE

−→ ∆a, rc [mm]

1

w

e Py

∅d

0.4 1

10

20 d [mm]

30 −→

40

20

30

40

50

TCD-CLB-FE approach with rc,FE = 1.143 mm further confirming the calibration methodology which also results in almost identical failure predictions in the plateau. Note that in case ∆a exactly becomes a constant and if the TCD-CLB-FE approach is calibrated to the FFM failure stress in the plateau both predictions are identical in here. This might be performed for w/d = 4 but is eventually neglected since the net section failure prediction is higher than the bearing cut-off within the diameter range d ≤ 50 mm investigated. The bearing cut-off is actually descended for d ≥ 62 mm. The proximity of rc,FE might be also interpreted as follows. The calibrated characteristic distance is obtained using the failure stress of the test and might be therefore considered as an experimentally determined crack length. Since ∆a calculated by means of the FFM lies near to it in the plateau this crack propagation model can be regarded as further confirmed. Contrary, the TCD-T contains rc,T = 1.632 mm significantly differing from the plateau and performs worse in comparison to its calibration counterpart. The quantity r p = rc,T /4 is the effective crack length of a mode I through crack in an infinite plate (Tada et al., 2000). Concerning the present finite bolted joint problem rc,T can be regarded

∆a

0.6

10

Fig. 16. Diameter range yielding an error magnitude | | ≤ 10 %.

1.2

w/d = {1.75, 2, 3, 4}

5-50 mm

TCD-CLB-FE

50

Fig. 15. Deviations of predicted failure stresses to FFM for w/d = 3

0.8

w/d = 2

1-15 mm

50

Fig. 17. Crack length ∆a at failure reproduced using results of Catalanotti and Camanho (2013). rc,Airy = 0.351 mm is out of the range of ∆a and not shown.

12

as fictional length without any physical meaning. Nevertheless the TCD-T-FE stays within the deviation limit for a certain diameter range. Concerning the TCD-T-Airy and TCD-CLB-Airy approaches the corresponding characteristic distances rc,T = 1.632 mm and rc,Airy = 0.351 mm are outside the value range of ∆a. Another impact on the prediction is exerted by the erroneous net section stresses. Depending how these two impacts interact with each other the prediction might be either worse, unchanged or even better compared to the TCD-CLB-FE approach. In particular the TCD-T-Airy yields the most inaccurate prediction of the study whereas the TCD-CLB-Airy is applicable in a wide diameter range for w/d = 3. Concerning w/d = 2 the TCDCLB-Airy can be even applied within 1 mm ≤ d ≤ 15 mm partially lying outside the plateau where its numerical counterpart violates the deviation limit. If using FE stresses there is no doubt that for accurate predictions rather large diameters (d ≥ 5 mm) in the plateau should be chosen to calibrate. However in the TCD-CLB-Airy approach the choice is not obvious due to the interaction of the characteristic distance with the erroneous net section stresses. To analyse this matter further predictions have been derived using failure stresses of different hole diameters for calibration. From the methodology point of view test data should be taken for calibration. To the authors’ knowledge there is no additional experimental data available and results of the FFM are taken assuming providing reliable predictions. Concerning w/d = 2 more precise results are achieved if rc,Airy is calibrated to failure stresses of rather small hole diameters 1 mm ≤ d ≤ 11 mm outside the plateau leading to good correlation nearby them. The results of tiny hole sizes as d = 1 mm are more of theoretical interest since in reality they are hardly manufactured. Any satisfactory result is achieved for large holes even if calibrating to them. Regarding w/d ≥ 3 no net section failure is expected for d < 13.5 mm since the bearing cut-off is ascended by the FFM tension failure prediction. Thus calibration is only performed for d ≥ 13.5 mm yielding acceptable predictions. The extreme cases of failure criteria enabling to model the hole size effect are d → 0 for which the failure stress specialises to the value of the plain material and d → ∞ where predictions coincide with local criteria. The failure stress of a finite hole diameter lies in between and for a lightweight optimal design sufficiently exploiting the material non-local criteria requiring accurate net section stresses and not only stress concentrations are essential. For further validation of the present failure calculus experiments should be conducted since the FFM is just a model of crack propagation mechanisms occurring in real structures with brittle material.

load direction at the straight free edges. Thus the force flux of the finite geometry joint is modelled in a realistic manner. However the net section stresses deviate the smaller the geometry ratios w/d and e/d are chosen due to non-zero tractions normal to the straight free edges of the plate. To further assess the criticality of these deviations a failure analysis is conducted and the calculated predictions are compared to reference values derived by Finite Fracture Mechanics (FFM) in Catalanotti and Camanho (2013). As failure criterion the Theory of Critical Distances (TCD) is chosen using characteristic distances calibrated from the experimental failure load as well as using the effective crack length of a mode I through crack in an infinite domain. The predictions of the present calculus yield deviations to FFM values | | ≤ 10 % for w/d = 2 in the diameter range 1 mm ≤ d ≤ 15 mm and regarding w/d = 3 within 8 mm ≤ d ≤ 354 mm. Furthermore the general concept of the TCD has been shown as capable to model the hole size effect of those defect sizes for which the FFM derived crack length ∆a can be approximated by a plateau value if the characteristic distance has about the same size. References Bickley, W. G., 1928. The Distribution of Stress Round a Circular Hole in a Plate. Philosophical Transactions of the Royal Society of London 227, 383– 415. Camanho, P., Erçin, G., Catalanotti, G., Mahdi, S., Linde, P., 2012. A finite fracture mechanics model for the prediction of the open-hole strength of composite laminates. Composites Part A: Applied Science and Manufacturing 43 (8), 1219 – 1225. Camanho, P., Lambert, M., 2006. A design methodology for mechanically fastened joints in laminated composite materials. Composites Science and Technology 66, 3004–3020. Camanho, P., Matthews, F., 1997. Stress analysis and strength prediction of mechanically fastened joints in frp: a review. Composites Part A: Applied Science and Manufacturing 28 (6), 529 – 547. Catalanotti, G., Camanho, P., 2013. A semi-analytical method to predict nettension failure of mechanically fastened joints in composite laminates. Composites Science and Technology 76, 69–76. Collings, T., 1977. The strength of bolted joints in multi-directional CFRP laminates. Composites 8 (1), 43 – 55. Collings, T., 1982. On the bearing strengths of CFRP laminates. Composites 13 (3), 241 – 252, jointing in fibre-reinforced plastics. Crews, Hong, C., Raju, I., 1981. Stress concentration factors for finite orthotropic laminates with a pin-loaded hole. NASA Technical Paper. De Jong, T., 1976. Spanningen rond een gat in een elastisch orthotrope of isotrope plaat, belast door een pen die zich daarin wrijvingsloos kan bewegen. Technische Hogeschool Delft, Afdeling der Luchtvaart-en Ruimtevaarttechniek, Rapport LR-223. De Jong, T., 1977. Stresses around pin-loaded holes in elastically orthotropic or isotropic plates. Journal of Composite Materials 11, 313–331. De Jong, T., Vuil, H. A., 1981. Stresses around pin-loaded holes in elastically orthotropic plates with arbitrary load direction. Delft University of Technology, Department of Aerospace Engineering, Report LR-333. Echavarría, C., Haller, P., Salenikovich, A., 2007. Analytical study of a pin–loaded hole in elastic orthotropic plates. Composite Structures 79 (1), 107 – 112. Grüber, B., Hufenbach, W., et al., 2007. Stress concentration analysis of fibrereinforced multilayered composites with pin-loaded holes. Composites Science and Technology 67, 1439–1450. Grüber, B., Gude, M., et al., 2018. Calculation method for the determination of stress concentrations in fibre-reinforced multilayered composites due to metallic interference-fit bolt. Journal of Composite Materials 52 (18), 2415– 2429. Hart-Smith, L. J., 1980. Mechanically-fastened joints for advanced composites — phenomenological considerations and simple analyses, 543–574.

7. Conclusion In the present paper the stress field for a symmetric finite bolted joint with quasi-isotropic composite laminate material under in-plane loading is determined by means of the Airy stress function. The bolt contact is idealised by sinusoidal radial tractions along the hole boundary. To obtain the solution use is made of auxiliary functions eliminating non-zero tractions in 13

Hebel, J., Becker, W., 2008. Numerical analysis of brittle crack initiation at stress concentrations in composites. Mechanics of Advanced Materials and Structures 15, 410–420. Kirsch, G., 1898. Die Theorie der Elastizität und die Bedürfnisse der Festigkeitslehre. Zeitschrift des Vereins deutscher Ingenieure 42, 797–807. Knight, R., 1935. The action of a rivet in a plate of finite breadth. Philosophy Magazine 19 (Series 7), 517–540. Kratochvil, J., Becker, W., 2010. Structural analysis of composite bolted joints using the complex potential method. Composite Structures 92 (10), 2512 – 2516. Kretsis, G., Matthews, F., 1985. The strength of bolted joints in glass fibre/epoxy laminates. Composites 16 (2), 92 – 102. Leguillon, D., 2002. Strength or toughness? A criterion for crack onset at a notch. European Journal of Mechanics - A/Solids 21 (1), 61 – 72. Lekhnitskii, S., 1968. Anisotropic plates. Gordon and Breach Science Publishers. Li, J., Zhang, X., 03 2006. A criterion study for non-singular stress concentrations in brittle or quasi-brittle materials. Engineering Fracture Mechanics 73, 505–523. Matthews, F., 1987. Theoretical stress analysis of mechanically fastened joints. Elsevier Applied Science Publishers Ltd, Joining Fibre-Reinforced Plastics„ 65–103. Michell, J. H., 1899. On the direct determination of stress in an elastic solid with application to the theory of plates. Proceedings of the London Mathematical Society 31, 100–124. Ogonowski, J., 1981. Effect of variances and manufacturing tolerances on the design strength and life of mechanically fastened composite joints AFWALTR-81-3041, 3. Technical report, McDonnell Aircraft Company. Rosendahl, P., Weißgraeber, P., Stein, N., Becker, W., 2017. Asymmetric crack onset at open-holes under tensile and in-plane bending loading. International Journal of Solids and Structures 113-114, 10 – 23. Sadd, M. H., 2005. Elasticity: Theory, Applications, and Numerics - Second Edition. Elsevier Butterworth-Heinemann. Tada, H., Paris, P., Irwin, G., 2000. The Stress Analysis of Cracks Handbook. ASME Press. Tan, S. C., 1988. Finite-width correction factors for anisotropic plate containing a central opening. Journal of Composite Materials 22 (11), 1080–1097. Taylor, D., 2007. The Theory of Critical Distances. Elsevier Science. Timoshenko, S., Goodier, J. N., 1951. Theory of Elasticity. McGraw-Hill Book Company, Inc. Waszczak, J. P., Cruse, T. A., 1971. Failure mode and strength prediction of anisotropic bolt bearing specimens. Journal of Composite Materials 5, 421– 425. Weißgraeber, P., Leguillon, D., Becker, W., 2015. A review of finite fracture mechanics: Crack initiation at singular and non-singular stress-raisers. Archive of Applied Mechanics 86, 375–401. Whitney, J. M., Nuismer, R. J., 1974. Stress fracture criteria for laminated composites containing stress concentrations. Journal of Composite Materials 8, 253–265. Xiao, Y., Ishikawa, T., 2005. Bearing strength and failure behavior of bolted composite joints (part i: Experimental investigation). Composites Science and Technology 65 (7), 1022 – 1031. York, J., Wilson, D., Pipes, R., 1982. Analysis of the net tension failure mode in composite bolted joints. Journal of Reinforced Plastics and Composites 1 (2), 141–152.

14