Aurélien Doitrand
email: aurelien.doitrand@gmail.com

Dominique Leguillon
email: dominique.leguillon@upmc.fr

Comparison between 2D and 3D applications of the coupled criterion to crack initiation prediction of scarf adhesive joints

Keywords: A. Epoxides, B. Steels, C. Fracture mechanics, C. Finite element stress analysis, Coupled Criterion

Comparison between 2D and 3D applications of the coupled criterion to crack initiation prediction of scarf adhesive joints

Aurélien Doitrand, Dominique Leguillon

Introduction

Adhesively bonded joints are convenient to assemble several kind of materials and are therefore used for the load transfer between structural elements in various industrial domains including aerospace, automotive or civil engineering [START_REF] Adams | Adhesive bonding : science, technology and applications[END_REF][START_REF] Da Silva | Modeling of adhesive bonded joints[END_REF]. In particular, scarf joint [START_REF] Kumara | Tensile failure of adhesively bonded CFRP composite scarf joints[END_REF] is an interesting kind of adhesive bonding that provides a larger adhesive bonding length than typical butt joints [START_REF] Chen | Stress distribution in plane scarf and butt joints[END_REF]. The design of adhesive bondings in structures requires a good knowledge of the mechanical response of the joints, including crack initiation prediction.

Fracture of scarf adhesive joint is often predicted by performing only a stress analysis.

Chen and Cheng [START_REF] Chen | Stress distribution in plane scarf and butt joints[END_REF] studied analytically the 2D stress distributions in a scarf joint as a function of the scarf angle and the elastic mismatch between the adhesive and the adherends.

Chaudhuri and Chiu [START_REF] Chaudhuri | Three-dimensional asymptotic stress field in the vicinity of an adhesively bonded scarf joint interface[END_REF] developed a three-dimensional eigenfunction expansion approach for the prediction of the singular stress field at the scarf joint corner. Suzuki determined through 2D Finite Element (FE) simulations that these singular stress fields vanish for a 52 • scarf angle under tensile loading [START_REF] Suzuki | Adhesive tensile strengths of scarf and butt joints of steel plates : 2nd report, relation between mechanical properties of adhesive and fracture criteria of joints[END_REF][START_REF] Suzuki | Adhesive Tensile Strengths of Scarf and Butt Joints of Steel Plates (Relation Between Adhesive Layer Thicknesses and Adhesive Strengths of Joints) : Solid-Mechanics, Strength of Materials[END_REF]. Through 3D FE simulations, He et al. [14] and Nakano et al. [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF][START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF] found that the singular stress at the interface corner was minimum for a 60 • scarf angle. In several works, scarf joint failure is predicted using a stress criterion in FE simulations. Kimiaeifar et al. [START_REF] Kimiaeifar | Reliability analysis of adhesive bonded scarf joints[END_REF] compared Von Mises, modified Von Mises, and maximum principal stress criterions associated to a probabilistic model and showed that the choice of the stress criterion was significant for scarf failure assessment [START_REF] Kimiaeifar | Reliability analysis of adhesive bonded scarf joints[END_REF]. Adin [START_REF] Adin | The effect of angle on the strain of scarf lap joints subjected to tensile loads[END_REF] studied the effect of scarf angle on joint strain and failure load under tensile loading using FE calculations.

The load at failure was determined when Von Mises stress reached the adhesive strength in the joint. Nakano et al. studied experimentally and numerically the influence of the scarf angle on the joint adhesive failure under tensile loading [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF] and bending [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF]. The failure stress was determined by comparing either Von Mises or maximum principal stress criterion inside the adhesive to the adhesive strength. However, due to the singularity at the adhesive/adherend interface corner and edges, this method seems mesh dependent (the finer the mesh, the higher the stress near the corner). Comparing Von Mises stress criterion to the adhesive strength might be adapted to predict the joint cohesive fracture since the crack initiates within the adhesive. However, in the case of adhesive fracture, the crack initiates at the adhesive/adherend interface. Therefore, this approach is no longer applicable because of the stress singularity at the interface corner, as shown by Afendi et al. [START_REF] Afendi | Strength prediction of epoxy adhesively bonded scarf joints of dissimilar adherends[END_REF]. In this case, the relevant fracture parameters that should be used for adhesive fracture study are those of the interface rather than those of the adhesive [START_REF] Afendi | Strength prediction of epoxy adhesively bonded scarf joints of dissimilar adherends[END_REF][START_REF] Doitrand | 3D application of the coupled criterion to crack initiation prediction in epoxy/aluminum specimens under four point bending[END_REF][START_REF] Leguillon | An attempt to extend the 2D coupled criterion for crack nucleation in brittle materials to the 3D case[END_REF].

From a physical point of view, crack initiation requires two necessary conditions to be fulfilled simultaneously. The first condition states that the energy dissipated by the nucleation of a crack must be higher than the one required to initiate this crack. The second condition compares the tensile stress to the tensile strength prior to crack initiation. Both conditions can be expressed respectively as an energy and a stress criterion. By coupling both criterions, crack initiation can be studied in the framework of Finite Fracture Mechanics (FFM) as proposed by Leguillon [START_REF] Leguillon | Strength or toughness? A criterion for crack onset at a notch[END_REF]. The coupled criterion has already been applied to various crack initiation configurations described in the extensive review proposed by Weißgraeber et al. [START_REF] Weißgraeber | A review of Finite Fracture Mechanics: crack initiation at singular and non-singular stress raisers[END_REF], including fracture initiation of interfaces [START_REF] Mantič | Interface crack onset at a circular cylindrical inclusion under a remote transverse tension. Application of a coupled stress and energy criterion[END_REF][START_REF] Mantič | Crack onset and growth at the fibre-matrix interface under a remote biaxial transverse load. Application of a coupled stress and energy criterion[END_REF][START_REF] Martin | Energetic conditions for interfacial failure in the vicinity of a matrix crack in brittle matrix composites[END_REF][START_REF] Martin | A twofold strength and toughness criterion for the onset of free-edge shear delamination in angle-ply laminates[END_REF][START_REF] Martin | A coupled strength and toughness criterion for the prediction of the open hole tensile strength of a composite plate[END_REF][START_REF] Martin | Initiation of edge debonding: coupled criterion versus cohesive zone model[END_REF][START_REF] Muñoz-Reja | Crack onset and propagation at fibre-matrix elastic interfaces under biaxial loading using finite fracture mechanics[END_REF] and joints [START_REF] Carrère | Comparison between models based on a coupled criterion for the prediction of the failure of adhesively bonded joints[END_REF][START_REF] Moradi | Strength prediction of bonded assemblies using a coupled criterion under elastic assumptions: Effect of material and geometrical parameters[END_REF]. Until now, the coupled criterion has mainly be applied in 2D. Since Leguillon extended the coupled criterion in 3D using matched asymptotic expansions [START_REF] Leguillon | An attempt to extend the 2D coupled criterion for crack nucleation in brittle materials to the 3D case[END_REF], only few 3D applications have been performed [START_REF] Doitrand | Damage onset modeling in woven composites based on a coupled stress and energy criterion[END_REF][START_REF] Doitrand | 3D application of the coupled criterion to crack initiation prediction in epoxy/aluminum specimens under four point bending[END_REF][START_REF] García | A numerical study of transverse cracking in cross-ply laminates by 3D finite fracture mechanics[END_REF]. Compared to the 2D case, where a crack is generally described by its length and/or its orientation, a 3D crack may be described by an infinite number of parameters [START_REF] Doitrand | Damage onset modeling in woven composites based on a coupled stress and energy criterion[END_REF][START_REF] Weißgraeber | A review of Finite Fracture Mechanics: crack initiation at singular and non-singular stress raisers[END_REF], which complexifies crack initiation analysis. Some assumptions can be made in order to reduce the number of parameters describing the crack shape [START_REF] Doitrand | Damage onset modeling in woven composites based on a coupled stress and energy criterion[END_REF][START_REF] García | A numerical study of transverse cracking in cross-ply laminates by 3D finite fracture mechanics[END_REF], hence preventing the real crack initiation shape determination. To overcome this difficulty, Doitrand and Leguillon recently proposed a method for 3D crack shape determination using the coupled criterion and applied it to interface crack initiation prediction in aluminum-epoxy specimens under for point bending [START_REF] Doitrand | 3D application of the coupled criterion to crack initiation prediction in epoxy/aluminum specimens under four point bending[END_REF]. The possible cracks were determined using the stress isocontours, so that the crack topology was described by a single parameter. Then, crack initiation shape was determined as the one fulfilling the coupled criterion that minimized the imposed loading.

In this paper, 2D and 3D applications of the coupled criterion are proposed that allows determining the crack initiation loading and topology. The studied cases, described in section 2, are crack initiation under tensile loading and bending of scarf adhesive epoxy joints with similar steel adherent [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF][START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. A brief recall of the coupled criterion is presented in section 3.

In section 4, 2D and 3D FE modeling of scarf joint are performed in order to determine crack initiation. In section 5, the adherent/adhesive interface fracture properties (strengths and toughness) are determined by confrontation of the numerical results to experimental data from [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. Then, the proposed method is applied to the blind prediction of crack initiation in scarf joints under bending [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF].

Experimental tests and data

The numerical problem under investigation is adhesive fracture initiation prediction in scarf joints under tensile loading and bending. The studied specimens are composed of SS400 steel adherends bonded by 3M Scotch-Weld1838 epoxy adhesive [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. The specimen dimensions are depicted in Figure 1. Some experimental tests aiming at determining the influence of the scarf joint geometrical parameters on the stress at failure have been performed under tensile loading [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF] and bending [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF]. In these works, Nakano et al. determined the failure stress as a function of the scarf angle (θ ={30,45,52,60,90} • ) for (i) a constant specimen width w=32mm under bending and tensile loading or (ii) a constant adhesive bonding length l=32mm under tensile loading. The studied specimens are h=9mm thick, L=200mm long and the adhesive thickness t is 0.1mm. The varying specimen dimensions in both cases are summarized in Tables 1 and2 For all tested configurations under tensile loading or bending, a crack initiates at the adherend/adhesive interface [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF][START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. Therefore, the relevant fracture parameters that should be involved in crack initiation analysis are the adhesive/adherend interface strength and toughness. It is difficult to identify experimentally these parameters by a direct approach.

Therefore, in the following, the interface fracture parameters are determined by indirect confrontation to the experimental data in the case of tensile loading [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF]. Then, a blind application to crack initiation prediction under bending [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF] is performed using the obtained interface fracture parameters.

The coupled criterion for crack initiation prediction

The coupled stress and energy criterion

Two separate conditions must be fulfilled in order to ensure crack initiation [START_REF] Leguillon | Strength or toughness? A criterion for crack onset at a notch[END_REF]. On the one hand, stress must be high enough to damage the material and on the other hand, sufficient energy is required in order to open the crack. Both conditions can be expressed mathematically as a stress and an energy criterion. The energy condition states that the incremental energy release rate G inc (i.e. the ratio between the change in potential energy W due to crack initiation and the crack area S) must be higher than the material critical energy release rate G c (Equation ( 1)). Under the assumption of elasticity, the potential energy of both the damaged and the undamaged material, and hence the incremental energy release rate, are proportional to the square applied load σ 0 (Equation ( 1)). In Equation ( 1), A only depends on the specimen geometry and of the material mechanical properties.

G inc = ∆W S = W (0) -W (S) S = A(S) • σ 2 0 ≥ G c (1)
If the crack shape can be described by a set of parameters d = (d 1 , .., d n ), A(S) can be determined by generating response surfaces of the potential energy change, ∆W , to the varying parameters d by means of FE calculations with the corresponding cracked configuration [START_REF] Doitrand | Damage onset modeling in woven composites based on a coupled stress and energy criterion[END_REF].

The stress condition states that at crack initiation, the crack surface must be overloaded.

It means that a stress criterion f must be fulfilled on the whole crack surface prior to crack nucleation. Under the assumption of elasticity, the stress is proportional to the applied load

Energy criterion fulfilled

Stress criterion fulfilled

Coupled criterion fulfilled

S S * σ * σ 0

Stress criterion

Energy criterion

Figure 2: Imposed loading at crack initiation as a function of the crack area obtained using stress (in blue) and energy (in red) criterions. Among all admissible initation configurations (zone in purple), crack initiation configuration is determined as the one that minimizes the imposed loading (i.e., for S = S * ). σ 0 (Equation ( 2)). Therefore, only one elastic calculation is needed to compute the stress condition.

f (σ) = k(S) • σ 0 ≥ 1 (2)
The crack configuration at damage initiation (i.e. d * and S * ) can be determined as the configuration for which both the stress and the energy conditions are satisfied that minimizes the imposed loading. If A and k are monotonic function of S, such a configuration is determined by combining and solving Equations ( 1) and (2) (cf. Equation (3)).

A(S) k(S) 2 = G c σ 2 c ( 3
)
The solution of the coupled criterion in this case is shown in Figure 2, which illustrates the imposed loading necessary to fulfill the stress (in blue) and the energy (in red) criterions as a function of the crack area S. The purple area represents all the configurations for which the coupled criterion is fulfilled, i.e., all cracks that are likely to initiate at a certain imposed loading. Among them, the crack initiation configuration is determined as the one that minimizes the imposed loading, which corresponds to a crack of area S * and a loading level σ * .

The case of non monotonic stress evolution

Under tensile loading, it has been shown that the stresses do not exhibit a monotonic evolution along the adhesive/adherend interface of scarf joints [START_REF] Carrère | 3D models of specimens with a scarf joint to test the adhesive and cohesive multi-axial behavior of adhesives[END_REF][START_REF] Davies | Influence of adhesive bond line thickness on joint strength[END_REF][START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. The different configurations that can thus be encountered when solving the coupled criterion are shown in Figure 3.

If the solution S * of Equation ( 3) is small enough so that σ * < σ min (Figure 3a), or high enough so that the σ * > σ max (Figure 3b), the crack initiation area and load are S * and σ * , similarly to a monotonic stress evolution (cf. Figure 2). However, the case where σ min < σ * < σ max leads to alternative crack initiation configurations (Figure 3c andd),

depending on the value of S * . If S * > S 1max (where S 1max is the crack area corresponding to the local maximum of the stress criterion, cf. Figure 3c), solving Equation (3) leads to an initiation stress σ * < σ max . However, as stated previously, the stress criterion must be fulfilled on the whole surface of the crack prior to crack initiation, which implies that a crack of area S * may initiate only for a crack initiation loading σ * ≥ σ max . Therefore, the crack initiation area and loading are S 2max and σ max rather than S * and σ * (Figure 3c). This configuration has been encountered by Leguillon et al. [START_REF] Leguillon | Application of the coupled stress-energy criterion to predict the fracture behaviour of layered ceramics designed with internal compressive stresses[END_REF] in fracture initiation prediction of layered ceramics. The last possible configuration is obtained if σ min < σ * < σ max and S * < S 1max (Figure 3d). It is clear that a crack of area S * (1 st crack in Figure 3d) is an admissible initiation crack. However, contrary to the other configurations, another crack of maximum area S 2 -S 1 (2 nd crack in Figure 3d) may also initiate since the stress criterion is also fulfilled in this zone. In this case, it is necessary to study the competition (using the energy criterion, see e.g. [START_REF] Martin | Competition between deflection and penetration at an interface in the vicinity of a main crack[END_REF]) between the nucleation of both cracks in order to determine which crack (and possibly both) initiates.

The choice of the stress criterion

Kimiaeifar et al. [START_REF] Kimiaeifar | Reliability analysis of adhesive bonded scarf joints[END_REF] demonstrated that the choice of the stress criterion is of primary importance for scarf joint failure prediction. Interface failure is classicaly determined using a stress criterion combining the tensile and shear stress to strength ratios [START_REF] Carrère | Comparison between models based on a coupled criterion for the prediction of the failure of adhesively bonded joints[END_REF][START_REF] García | Mixed-mode crack initiation at a v-notch in presence of an adhesive joint[END_REF][START_REF] Liao | Effect of adhesive thickness, adhesive type and scarf angle on the mechanical properties of scarf adhesive joints[END_REF][START_REF] Martin | Initiation of edge debonding: coupled criterion versus cohesive zone model[END_REF]. In the studied case, it can be shown by studying an idealized model with an infinitesimally small interface that such a criterion may not be adapted to failure stress prediction.

Stress criterion

Energy criterion The adhesive thickness (t = 0.1mm) is really small compared to the specimen dimensions (cf. Section 2). In a first approximation, a simple scarf joint specimen model with an infinitesimally small interface under uniform tensile loading may be studied. The interface normal (σ nn ) and shear (τ ns ) stresses can be expressed as a function of the imposed stress σ 0 and of the scarf angle θ (Equation ( 4)).

σ * S * S * σ * σ * σ * (a) (b) (c) (d
   σ nn = σ 0 • sin 2 θ τ ns = σ 0 • cosθ • sinθ (4)
Therefore, a criterion (f 1 ) combining shear and normal stress to strength ratios may be deduced from Equation (4) (cf. Equation [START_REF] Carrère | 3D models of specimens with a scarf joint to test the adhesive and cohesive multi-axial behavior of adhesives[END_REF], where σ c and τ c are the interface normal and shear strengths).

f 1 (θ, σ 0 ) = ( σ nn σ c ) 2 + ( τ ns τ c ) 2 = σ 0 σ c • sin 2 θ • (sin 2 θ + ( σ c τ c ) 2 • cos 2 θ) (5)
For a given scarf angle θ, the failure stress σ failure (θ) such that f 1 (θ, σ failure (θ)) = 1 can be determined as a function of the normal and shear strengths. Figure 4a shows the obtained failure stress evolution as a function of the scarf angle for several shear to normal strength ratios. It is clear that, whatever the tensile to shear strength ratio, this criterion does not allow reproducing the experimentally observed failure stress trend that increases from

45

• to 60 • scarf angles and then decreases for 90 • scarf angle [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. This trend can, however, qualitatively be well reproduced by employing a stress criterion consisting of the maximum between normal or shear stress to strength ratios (Figure 4b). Coupled with the energy criterion, this stress criterion will therefore be used in the following for scarf joint crack initiation determination.

Scarf joint and crack initiation modeling

The application of the coupled criterion to the studied case requires the modeling of scarf joint specimens and boundary conditions representing bending and tensile loading. The full modeling of four point bending or tensile test as illustrated in [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF] and [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF] would be nonlinear, due to the contact between the spans and the specimens. Talmon l'Armée et al. [START_REF] Talmon L'armée | Nonlinear crack opening integral : Mode mixity for finite cracks[END_REF] highlighted the possiblity to apply the coupled criterion in a nonlinear framework, however triggering much more complexity than in a linear elastic case. Therefore, the choosen boundary conditions, depicted in Figure 5a-b in the 2D case, ensures a linear elastic framework for the coupled criterion application, which means that the stress and the energy are respectively proportional to the imposed loading and the square imposed loading. For all the scarf joint configurations indicated in Tables 1 and2 been generated and extruded in order to obtain 3D meshes, the meshes being refined in the vicinity of the adhesive/adherend corner (Figure 5c-d). The material mechanical properties used in the FE calculations are taken from [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF] : E adherend = 209 GPa, ν adherend = 0.29 for steel and E adhesive = 3.34 GPa, ν adhesive = 0.38 for epoxy.

The computation of the energy criterion requires the calculation of the potential energy of both the undamaged and damaged configuration, which requires the generation of the corresponding cracked meshes. 2D cracked mesh generation is straightforward since the adherend/adhesive interface crack is described by only one parameters, i.e., its length.

Therefore cracked meshes are obtained by doubling the nodes along the adhesive/adherend interface to obtain a cracked mesh with a given crack length. 3D cracked mesh generation is more complex, since it requires the determination of the 3D crack shape. The same method as proposed in [START_REF] Doitrand | 3D application of the coupled criterion to crack initiation prediction in epoxy/aluminum specimens under four point bending[END_REF] is used herein, which consists in determining the possible crack shape based on the stress criterion isocontours. The advantages of this method is that the 3D crack topology is described by a single parameter, and that the stress condition of the coupled criterion (i.e. the stress criterion must be fulfilled on the whole surface of the crack prior to crack initiation) is strictly fulfilled. For each studied configuration, a first calculation has been carried out in order to determine the 2D (under plane strain assumption) and 3D stress fields at the adhesive/adherend interface and hence compute the stress condition of the coupled criterion (cf. Section 3.1) and the stress criterion isocontours in order to determine the possible crack initiation shapes. Figure 6a shows the interface normal stress field under bending for a 15 • scarf joint specimen whose dimensions are given in section 2 and Table 1. Given a value of the stress criterion isocontours (depicted in Figure 6b), an automated approach has been developped in order to generate a 2D mesh of the interface whose topology includes the crack shape based on this isocontour value (Figure 6c). The 3D scarf joint specimen mesh is then obtained by extrusion and shearing of this 2D mesh. The nodes located on the crack surface are finally unbuttoned in order to obtain the cracked mesh and compute the potential energy release (cf. Section 3.1). It can be noted that the same mesh topology is used to compute the potential energy release, the unbutonned crack surface nodes being the only difference between the undamaged and the crack mesh, hence reducing the numerical errors. For each specimen configuration (cf. Tables 2 and1), 15 stress isocontour values have been chosen and the corresponding meshes have been generated in order to compute W (S).

The FE calculations have been performed using Zébulon FE solver [36]. In previous studies [START_REF] Carrère | 3D models of specimens with a scarf joint to test the adhesive and cohesive multi-axial behavior of adhesives[END_REF][START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF], it has been shown that the stress singularity at the adherend/adhesive interface corner (and therefore crack initiation) depends on the scarf angle.

In the present work, the interface corners for which the stresses are singular have been determined using calculations on undamaged configurations. In 2D (respectively in 3D), the stress field is singular at two (respectively four) adherend/adhesive corners, symmetrically located with respect to the joint center. In 2D, whether one or two corner cracks initiate can be determined by computing the incremental energy release rate for both configurations (Figure 7), the stress criterion being fulfilled for the same loading level for both cracks due to the stress field symmetry. For small enough cracks, similar incremental energy release rates associated with one or two crack initiation are obtained, since both cracks are far enough so that they do not interact with each other. For larger crack, however, both cracks interact and one crack is more likely to initiate than two, which might be due to a shielding effect. Similar results are obtained in 3D, leading to the conclusion than a single corner crack is more likely to initiate. The same phenomenon was also observed for 3D crack initiation in bimaterial specimen under four point bending [START_REF] Doitrand | 3D application of the coupled criterion to crack initiation prediction in epoxy/aluminum specimens under four point bending[END_REF]. Therefore, in the following, all the presented results concern the initiation of a single corner crack.

In section 3.2, it has been shown that for a non monotonic stress evolution, the competition between a crack starting from the corner or not (both crack lengths (in 2D) or shapes (in 3D) being determined based on stress isovalues) must be studied. This can be determined by computing the incremental energy release rate for the following configurations : (i) a single corner crack (depicted as the 1 st crack in Figure 3), (ii) a single non-corner crack (illustrated as the 2 nd crack in Figure 3) and (iii) both cracks. Similar incremental energy release rates associated to either a single non-corner crack or both cracks are obtained, which are much larger than the one associated to a single corner crack initiation. Therefore, either a single non-corner crack or both cracks are more likely to initiate than a single corner crack. This may be explain by the stress criterion local maximum (cf. Figure 3) that is attained for a really small crack length (in 2D) or surface (in 3D) for which the energy release rate is small.

Numerical and experimental crack initiation results

2D and 3D applications of the coupled criterion to crack initiation prediction in scarf joints require the knowledge of the adhesive/adherend interface critical energy release rate and normal and shear strengths. These parameters are difficult to identify experimentally in a direct maneer, but can be determined by indirect confrontation to experimental results, which is done using the experimental results under tensile loading from [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF]. Then, crack initiation is predicted in the case of bending and confronted to the experimental results from [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF].

Crack initiation under tensile loading

The coupled criterion presented in section 3.1 is applied to crack initiation determination in scarf joints under tensile loading. By indirect confrontation to the experimental results taken from [START_REF] Nakano | Three-dimensional FEM stress analysis and strength prediction of scarf adhesive joints with similar adherends subjected to static tensileloadings[END_REF], the interface fracture parameters are determined in both 2D and 3D cases. This is done by minimizing the difference between experimental crack initiation stresses and numerical predictions in a least square sense. for a constant specimen width ((Figure 8(a andc)) or bonding length (Figure 8(b andd))

using the fracture parameters obtained by indirect confrontation to the experimental results.

Even if the fracture parameters have been determined so as to best fit the experimental results, a satisfying agreement between numerical and experimental results is obtained since the experimental trend is well reproduced numerically. 3D simulations lead to a slighty better quantitative agreement with the experimental data compared to 2D results.

Crack initiation prediction under bending

The coupled criterion is now applied to crack initiation prediction in scarf joint specimen under bending, using the fracture parameters determined in the previous section. Figure 9 shows crack initiation bending moment as a function of the scarf angle obtained experimentally [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF] and numerically through 2D (Figure 9a) and 3D (Figure 9b) simulations. The stress evolution observed experimentally is qualitatively well reproduced using both 2D and 3D approach. Similarly to the tensile loading case, 3D simulations lead to a better quantitative agreement with experimental data then the 2D approach.

Conclusion

The proposed application of the coupled criterion allows adhesive fracture initiation pre- Figure 9: Bending moment at crack initiation as a function of the scarf angle θ obtained experimentally [START_REF] Nakano | FEM stress analysis and strength prediction of scarf adhesive joints under static bending moments[END_REF] and numerically using (a) 2D and 3D applications of the coupled criterion (CC).

a function of the scarf angle. This trend is well repoduced using a stress criterion consisting of the maximum shear or tensile stress to strength ratio. 3D crack shapes are determined based on the stress criterion isocontours. Due to the non monotonic stress evolution at the adherend/adhesive interface under tensile loading, a non corner crack is more likely to initiate than a crack starting from the corner.

Adhesive/adherend toughness and tensile and shear strengths are determined by indirect confrontation to experimental results under tensile loading in the case of a constant width or adhesive bonding length. While similar values of tensile strength and toughness are obtained through 2D or 3D simulations, a lower shear strength is obtained using 3D simulations. The blind application of the coupled criterion to crack initiation in scarf joint specimens under bending leads to a good qualitative and quantitative agreement with experimental crack initiation stress evolution. Both under tensile loading and bending, 3D results give a slight better quantitative agreement with experimental data.

Figure 1 :

 1 Figure 1: The scarf joint specimen consisting of two steel adherends bonded by an epoxy layer adhesive.

Figure 3 :

 3 Figure 3: Possible coupled criterion solutions for a non-monotonic stress evolution. The purple zone indicates the configurations fulfilling the coupled criterion, among which the crack initiation configuration area S * and loading σ * can be determined (a-c). Depending on the energy criterion, a second crack is also likely to initiate (d).

Figure 4 :

 4 Figure 4: Failure stress as a function of scarf angle obtained using a stress criterion f consisting of (a) a combination of or (b) the maximum between tensile and shear stress to strength ratios for several shear to normal strength ratios.

Figure 5 :

 5 Figure 5: (a) 2D FE model boundary conditions for (a) tensile loading and (b) bending. (c) 2D and (d) 3D meshes of a specimen with a 30 • scarf angle and zoom on the adherend/adhesive interface corner.

Figure 6 :

 6 Figure 6: 15 • scarf joint specimen adherend/adhesive interface normal stress (a) field and (b) isovalues. (c) 2D mesh whose topology follows the crack shape determined using the stress isocontour value indicated in (b).

Figure 7 :

 7 Figure 7: 2D Incremental energy release rate as a function of the normalized cracked length for (a) 15 • and (b) 45 • scarf joints.

o= 15 . 2 Figure 8 :

 1528 Figure 8: Crack initiation stress σ failure as a function of the scarf angle θ obtained experimentally [31] and numerically using (a-b) 2D and (c-d) 3D applications of the coupled criterion (CC) in the case of (a and c) a constant specimen width w and (b and d) a constant specimen bonding length l.

  diction in scarf joints. Crack initiation topology and critical load under tensile loading and bending are determined through 2D and 3D FE simulations of undamaged and cracked scarf joint specimens. Unlike most approach involving a single stress criterion (generally Von Mises or Maximum principal stress criterion) for scarf joint crack initiation prediction, both an energy criterion and a stress criterion are employed. A criterion consisting of the combination of normal and shear stress to strength ratios, classicaly used for interface crack initiation prediction does not allow reproducing the experimentally observed crack initiation stress as

Table 1 :

 1 . Adhesive bonding length l as a function of the scarf angle θ in the case of a constant specimen width

	Scarf angle θ ( • )	30	45	52	60	90
	Adhesive bonding length l (mm) 64 45.25 40.6 36.95 32.0
	w=32mm.
	Scarf angle θ ( • )	30	45	52	60	90
	Specimen width w (mm) 16.0 22.6 25.2 27.7 32.0

Table 2 :

 2 Specimen width w as a function of the scarf angle θ in the case of a constant adhesive bonding

	length l=32mm.


Claudiu Herteliu
email: claudiu.herteliu@gmail.com

Peter Richmond
email: peterrichmond@ymail.com

Bertrand M Roehner
email: roehner@lpthe.jussieu.fr

Mr Je

Ap Jl

Se De Oc

Ja No

Fe De Mr

Ja Ap Fe

Se De

Oc Ja

Deciphering the fluctuations of high frequency birth rates

Keywords: birth, fluctuations, holidays, religion, mass mortality, heat waves

published or not. The documents may come

Deciphering the fluctuations of high frequency birth rates

Claudiu Herteliu, Peter Richmond, Bertrand

Introduction

That the determinants of conceptions and births are still not well understood is shown fairly clearly by the fact that measures taken to boost birth rates by governments of countries experiencing shrinking populations, e.g. Germany, Japan, Singapore or South Korea has proved largely ineffective. The fertility rates of these countries have been below 1.50 (that is to say more than 30% lower than the minimal reproduction level of 2.1) since 1995. Despite incentives introduced by the respective governments their fertility rate curves remain flat or decreasing. In 2015 the fertility rates in Singapore and South Korea were as low as 1.25.

Here, we do not focus on such long-term changes but on short-and medium-term changes. Nevertheless, some effects studied here may also have long-term implications. For example, consider the effect on births of religious interdicts. In the Middle Age in western Europe only about 70 days were left annually for permitted sexual relations. Had they been strictly followed, such rules would have drastically restricted fertility.

Seemingly random fluctuations Fig. 1a,b represents daily birth numbers in France in 1969 and 2015. This fast succession of upward and downward spikes resembles a graphical representation of white noise except that instead of being stationary the signal displays also some mediumterm fluctuations. From 1969 to 2015 there was a marked increase in the standard deviation that will be explained in a short moment.

Flawed tentative interpretations

To make sense of such a jerky graph seems a real challenge. Usually, the standard approach is to perform a moving window average which will eliminate high frequency changes. The remaining ups and downs are then interpreted in terms of specific human activities, e.g. vacation, harvesting in rural societies and so on which are assumed either to favor or reduce conception rates. Examples of such interpretations can be found in [START_REF] Houdaille | Mouvement saisonnier des conceptions en France de 1740 à 1829 [Seasonal pattern of conceptions in France from 1740 to 1829[END_REF][START_REF] Houdaille | Le mouvement saisonnier des naissances dans la France rurale de 1640 à 1669[END_REF]. However, there is a fundamental flaw in this approach for one does not know the effect that vacationing and harvesting have on the sexual behavior of people. Any assumption made in this respect will be just wishful thinking. Take harvesting for example. In former rural societies it was a time of hard labor but also of celebration especially when the harvest was good. Which effect dominates?

Fact-driven methodology

We shall avoid such wishful thinking1 . In fact, we take the problem by the other end. The blue labels give calendar months whereas the red labels in italic indicate the corresponding conception months. Thus a conception in April would result in a birth in January of the following year. The same notation is used throughout the present paper. For the purpose of making relative variations more readable the series was normalized by making its average equal to 100. The substantial increase in dispersion has a simple explanation given below in the "Days off" section. Source: Répartition quotidienne des naissances vivantes, France métropolitaine. Institut National de la Statistique et des Études Économiques (INSEE). [Daily distribution of live births. Overseas territories are excluded. National Institute for Statistics and Economic Studies]. Available on Internet under the title: "Les naissances en 2016. Tableaux de séries longues".

We start by identifying remarkable patterns. For example, discarding leap years, the average of daily births reveals huge dips. Then we connect such patterns with sharply defined events, for instance public holidays. In this way we can build a kind of dictionary of the connections between specific patterns and the events which triggered them.

Parallel with the deciphering of Egyptian hieroglyphs

The graphs of the time series of births can be seen as pictographs the meaning of which we try to discover. Pictographs were also used in the earliest forms of writing. More specifically, the challenge of the present research can be described through a parallel with the situation experienced in 1822 by Jean-Franc ¸ois Champollion (1790-1832) after he was able to decipher the hieroglyphs displayed on the Rosetta Stone2 . It gave him some clues in the form of phonetic signs and ideograms but it took several decades until a fairly comprehensive understanding of the language could be "The rules of the sociological method" [START_REF] Durkheim | Les règles de la méthode sociologique[END_REF]. He suggested to investigate social behavior in the same way as one studies natural phenomena that is to say from "outside". Nowadays, his advice seems to be largely forgotten or overlooked. achieved. Our present situation is similar in the sense that we are able to explain a number of patterns which gives us some clues but many unsolved questions remain. For instance, later on we describe a new pattern that appeared around 2000 and for which we can at present only offer a tentative explanation. The deciphering of the Egyptian hieroglyphs was a long-term collective endeavor. Likewise we hope that the deciphering of the patterns of vital events (not only births but also deaths, marriages and divorces) will attract the attention of other researchers.

Outline of the paper

The paper proceeds through the following steps.

(1) First, we identify daily spikes attributable to holidays during which medically induced births are drastically reduced. In this way, all national holidays which occur on fixed days (such as 4 July in the US or 14 July in France) will become apparent. This is a fairly easy step which only requires an averaging procedure over a sufficiently large number of years.

(2) Secondly, we focus on adverse conditions which create troughs nine months later. This part builds on results obtained in a recent paper (Richmond et al. 2018a).

(3) The two previous parts give interpretations of dips and narrow troughs but not of wide troughs nor of peaks. Here we show that broad troughs can result from religious interdicts.

"Days off" effect Fig. 2a shows the curve that one obtains by summing up 35 annual curves (restricted to non-leap years) such as those represented in Fig. 1a,b.

Reduction in the dispersion

Not surprisingly, the dispersion is reduced. Whereas the curve of Fig. 1a moves within 100 ± 10 the curve of Fig. 2 moves within 100 ± 5.

If successive years were independent their average would have a standard deviation 6 times (i.e.

√ 35) smaller than the initial series. The fact that the dispersion is only divided by 2 shows that the successive annual series are correlated. Let us consider this point more closely. The formula which gives the standard deviation of an average of n pair-wise correlated random variables

Y n = (1/n)ΣX i is: σ(Y n ) = [σ/ √ n] 1 + (n -1)
r where σ is the standard deviation of the X i and r their pair-wise correlation (Richmond et al. 2018b). The special case r = 1 has a clear intuitive interpretation. A correlation of 1 means that all X i are identical: X i ≡ X, that is to say represented by the same time-series. Then Y n = X and σ(Y n ) = σ which means that in this case averaging does not at all reduce the standard deviation. If, on the contrary, r = 0 one is in the well-known case of the average of independent variables which results in:

σ(Y n ) = σ/ √ n. Now, in the present case one has: √ 1 + 34r = √ n/[σ/σ(Y n )]
= 6/2 = 3 which leads to r = 0.24. In words, the series are moderately correlated but nevertheless the fact that there are as many as 35 series results in a sizable effect on the standard deviation of the average.

Identification of the cause of the downward spikes

However, there are downward pointing spikes which reach well beyond this interval. From their positions it is easy to recognize that they correspond to fixed public holidays. The reduction in the number of births on non work days relates to medical induction of labor where drugs are used to induce labor in advance of it beginning naturally. Obviously such interventions will rather be planned on working days. In 2014 some 23% of American pregnant women experienced labor induction up from a percentage of 10% in 1990. Cesarean delivery has also increased: from 15% in 1996 to some 22% in 2011. Note that the figures for induced labor and Cesarean delivery should not be added together since the former may lead to the latter. 3 . The proportion of induced delivery is about the same in France (Enquête nationale périnatale [Perinatal national survey] 2010).

For the period of 35 years following 1969 (excluding leap years) the average reduction in birth numbers for public holidays is about 12%.

Effect of Saturdays and Sundays

On Fig. 2a one cannot clearly see the effect of mobile public holidays. Examples are Easter or "Labor Day" in the US which is on the first Monday of September. For the same reason, the effect of Saturdays and Sundays does not appear. However, the weekends appear in Fig. 2b. Between Fig 2a and2b there is only a small difference in design but it makes a drastic difference in shape. In Fig. 2a, because the annual series start at 1 January which can be any day of the week the weekends are positioned fairly randomly with the result that in successive years the summation mixes weekends and working days. However in Fig. 2b the annual series start at the first Monday of January with the result that the weekends are at the same locations in successive years; thus, the summation singles them out in a clear way. By direct examination of the data, one finds the following reductions in birth numbers for 1969 and 2016 (Table 1). Notes: The data are for France. The reduction indicated in the table refers to the difference between Saturday or Sunday and the average of adjacent working days. Sources: Same as for Fig. 1. The topic of weekend birth rates is also studied and discussed in [START_REF] Lerchl | Where are the Sunday babies? Observation on a marked decline on weekend births in Germany[END_REF], [START_REF] Lerchl | Where are the Sunday babies? II. Declining weekend birth rates in Switzerland[END_REF] and [START_REF] Lerchl | Where are the Sunday babies? III. Caesarean sections, decreased weekend births, and midwife involvement in Germany[END_REF].

The more induced Cesarean deliveries the deeper the dips in birth numbers on Saturday and Sunday.

Removal of the "Days off" effect

In order to estimate the global effect of days off on birth numbers, we replaced these data by the average of the births observed on adjacent working days. The resulting shape (see Fig. 3) is now much less jerky than in Fig. 1a. The removal of the days off was performed in 3 steps which resulted in the following reductions of the standard deviation: Initial: σ 1 = 6.1, No public holidays: σ 2 = 6.0, No Sundays: σ 3 = 5.4, No Saturdays: σ 4 = 5.3

Effect of adverse living conditions on conceptions

The removal of the "Days off" effect has a dramatic impact on the shape of the curve of daily births, but it is a fairly obvious effect. Less trivial is the effect resulting from adverse living conditions which can be revealed through a surge in mortality. This effect is less obvious for two reasons.

• It occurs at time of conception that is 9 months before the births.

• From conception to birth there are several successive steps which involve social as well as biological phenomena. It is through a comparative examination of several case studies that the key role of social as opposed to biological factors was revealed (Richmond et al. 2018a,b).

In this section, this effect will be examined using two cases.

(1) The impact of influenza outbreaks on conceptions.

(2) The impact of heat waves on conceptions.

Effect of influenza on January births

In Richmond et al. (2018a) the authors used data for influenza outbreaks of exceptional gravity such as in 1889 or 1918. Here we consider "ordinary" outbreaks of the kind that occur every year either in December or (most often) in January. Naturally, the effect is less spectacular but this analysis will show that even small outbreaks have a sizable effect on conceptions.

Our strategy is to select years which have the highest January death rates. In this respect it is important to observe that there is a high degree of correlation between January death rates either for influenza or for all causes of deaths. For instance, in the United States for 18 years from 1999 to 2016 the correlation is as high as 0.82. This is a remarkable result for influenza and pneumonia (code J09-18 of the International Classification of Diseases) deaths in fact represent only a small fraction (about 3%) of the total deaths of January. Despite this small fraction influenza nevertheless controls the January deaths because, apart from influenza changes, the total January deaths are very stable having a coefficient of variation of only about 5%. Thanks to the high correlation, we can select the appropriate years by using either the monthly influenza deaths or, in the countries where such data are not available, the monthly total deaths..

The "death-nobirth" effect

Fig. 4a shows that the death peak of January 1949 (which is the largest in the whole time period) has an amplitude peak/average of 1.8; this can be read easily on the graph because the curve was normalized so that the average corresponds to an index value of 100.

In order to facilitate the comparison between the death and birth curves, the latter was moved 9 months to the left (which in fact transforms it into a conception curve) and inverted through a change in its sign. The scale on the left-hand side shows that the amplitude of the birth curve is only 1.1. As shown in Richmond et al. (2018a), this is a general property: whether the mortality is due to influenza, food scarcity or earthquakes, the amplitude of the birth trough is always smaller than the amplitude of the death peak.

For the sake of brevity, the present effect will be called the death-nobirth effect. Fig. 4a The effect on births of an influenza outbreak in France in January 1949. There are influenza outbreaks every year but the one of January 1949 was particularly severe which is why it was selected. The birth series has been moved 9 months to the left so that it represents conceptions rather than births; moreover it has been inverted so that the birth trough appears as a peak. This kind of representation was introduced in

Relationship between the amplitudes of death peaks and birth troughs

In the previous subsection we considered a single, fairly spectacular, case. However in order to find out the relationship between the death peaks and the ensuing birth troughs we need to select and analyze a whole set of cases. First of all we explain how these cases should be selected.

For low January death rates, the impact on births will be masked by the background noise. In other words, in such cases the effect may still exist but it will not be visible.

In Fig. 4b we selected the years in the interval 1949-1974 for which the January deaths were higher than the annual average. Why limit ourselves to the interval 1949-1974? At the end of this section we show (in Fig. 7b) that toward the end of the century a new effect appeared which created secondary troughs some 2 or 3 months after October. Although the influenza troughs are still present the secondary troughs complicate the picture. By limiting ourselves to 1947-1974 we avoid this difficulty. A b ∼ A α d , α = 0.12 ± 0.04 The fact that similar exponents are obtained in fairly different conditions shows that these relationships are fairly robust. In what respect are the conditions different? First, the countries are not the same. Secondly, as already observed, the US data are for exceptional massive mortality numbers whereas the French data are for "normal" influenza mortality. Thirdly, the US data are for different states in the same year (1918) whereas the French data are for the whole country in different years.

Effect of heat waves on death rates

It seems common sense that heat waves may affect fragile persons particularly babies and elderly people. Naturally, one would expect this effect to be more pronounced in countries (such as France or Switzerland considered below) where air conditioning is not common. Although there have been many papers which document this effect only few give age-specific results. Among the latter, two fairly recent papers are particularly noteworthy in this respect, namely [START_REF] Rey | Vagues de chaleur, fluctuations ordinaires des températures et mortalité en France depuis 1971 [Heat waves, temperature fluctuations and mortality in France from 1971 to 2003[END_REF] for France and [START_REF] Vicedo-Cabrera | Excess mortality during the warm summer of 2015 in Switzerland[END_REF] for Switzerland.

For the present investigation it is particularly important to know which age groups are most affected since, if only babies and elderly are affected, there will be no reduction in births 9 months later. Fortunately, the paper by [START_REF] Rey | Vagues de chaleur, fluctuations ordinaires des températures et mortalité en France depuis 1971 [Heat waves, temperature fluctuations and mortality in France from 1971 to 2003[END_REF] has a table which gives all the information that we need (see Table 2). Table 2 shows that, not surprisingly, young people under 35 are under-represented among heat wave victims while elderly people over 75 are over-represented. However, the same table shows that for influenza the proportions are almost the same. Actually, we are not interested in the persons who die but rather in those who are affected and survive. Thus, the real question is: "Are influenza survivors more affected than heat wave survivors?" By studying the effect on births, we can answer this question at least for the specific meaning of "being affected in terms of conception capacity". In the next subsection we will see that influenza affects surviving people more seriously than heat waves do.

How conceptions are affected by heat waves

Whereas many papers have considered the effect of heat waves on death rates, to our best knowledge only one author has investigated the effects on birth rates nine months later. In a pioneering study Arnaud Régnier-Loilier (2010a) analyzed the 1975,1976,1983,1990,2001,2003. The influenza data are for the United States and are an average for 1999-2016. For a same population there are about 2.6 times more influenza deaths than heat wave deaths. Moreover, contrary to influenza outbreaks which basically occur every year, heat waves take place on average only every 5 years. Sources: Heat waves: Rey et al. (2007, Table 1). Influenza: CDC-Wonder database.

effect on conceptions and births of 3 of the 6 heat waves that occurred in France, namely: 1976France, namely: , 1983France, namely: and 2003. . Here the percentages of deaths with respect to annual deaths are of the order of a few percent that is about 10 times smaller than for influenza deaths. This makes the analysis more tricky. In order to remove the seasonal fluctuations Régnier-Loilier divided the births of the year under consideration by the average of the series of the two adjacent years. This leads to the identification of birth troughs whose width is about 40 days.

Such shocks belong to the same death-nobirth category as those for famine, influenza or earthquakes considered above.. That is why it makes sense to plot the results in the graph of Fig. 4b. We can see that for the same percentage of deaths one gets smaller conception troughs. This provides an answer to the question raised above. It turns out that for the same number of fatalities the heat wave survivors are less affected than the influenza survivors, intuitively a rather natural conclusion. This statement is a description in words of the fact that the data points for the heat waves are above the data points of influenza. As the y-axis has a negative scale, "above" means a smaller birth troughs.

We conclude the discussion of the effect of heat waves on conceptions by explaining in more detail the calculation procedure used in Régnier-Loilier (2010a).

Identification of the birth troughs due to heat waves in France

Between 1975 and 2005 there were 6 heat waves but the one of 2003 was by far the most severe. In this case there is absolutely no doubt about its effect on conceptions; this effect appears clearly in Fig. 7b.

In the other 5 cases the heat wave was much weaker. As a result, the birth troughs do not appear on the monthly data. Régnier-Loilier (2010a) used daily birth data. Due to the large daily fluctuations (illustrated at the beginning of this paper), the expected birth trough are hidden by the background noise. To make them visible the following operations were performed.

• First, each annual birth series b was "renormalized" through division by the average of the two adjacent birth numbers, e.g. b = (1977) was replaced by b ′ = 1977/ [(1976 + 1978)/2]. This is a fairly common way to smooth the series and reduce variability. However that is not what is observed here: in the three cases of 1977,1984,2004 the coefficient of variation (CV) of the initial series is smaller than the CV of the renormalized series; for the three cases the ratio is on average 0.97,

• As the series are still fairly bumpy further smoothing is performed via a moving average with a window width of 21 days. These operations lead to series which have more ups and downs than the initial series but which also display birth troughs nearly where expected. In short, the effect may well exist but due to its small magnitude it is at the limit of being detectable.

It is natural to test if, as in Fig. 4b, there is a correlation between death peaks and birth troughs. As such a test would be meaningless with only 3 cases, we repeated the calculation done in Régnier-Loilier (2010a) and extended it to the three other heat waves. However, it turns out that no correlation could be detected. One should perhaps not be surprised for a sample of only 6 cases can hardly provide a clear answer unless the correlation is fairly high which is not expected here.

Removal of the death-nobirth effect for daily births

To what extent are the fluctuations of births shown in Fig. 3 reduced when the deathnobirth effect is removed, that is to say when the dips due to excess-deaths 9 months earlier are compensated for. This compensation process is made by using the linear regression equation for the data of Fig. 4b. The standard deviation is reduced from the previous value of σ 4 = 5.3 to σ 5 = 4.7 (a reduction of 11%). This is a fairly small reduction and Fig. 5 explains why it is not larger. In order to affect the births of 1969, the deaths must occur in the time interval shown in magenta in Fig. 5.5 ; within the two years 1968-1969 this defines a specific time interval. In this interval there are only small sub-intervals (marked in red) which have excess-deaths. Therefore the curve of the births during 1969 is only slightly changed. Despite this small impact the standard deviation of births is nevertheless reduced by 11% showing that the correction works in the expected way, namely by pushing up the troughs. The modified birth curve is omitted because visually the changes are rather subdued.

Removal of the death-nobirth effect for monthly births

Graphically the effect of the removal of the death-nobirth effect becomes clearer at the level of monthly births, as shown in Fig. 6. It can be seen that all deep

Months Monthly births (in % with respect to mean)

Monthly births, compensated series Monthly births, initial series Fig. 6 Removal of the death-nobirth effect, France 1950-1954. For the compensation we used the regression coefficients a, b obtained in Fig. 4b; then the compensation was done as follows: ∆b = a∆d + b, where ∆d, ∆b are the excess deaths and excess births respectively; then b c = b i + ∆b where b i , b c are the initial and compensated birth series respectively. Note that this compensation was done only when there were excess deaths, that is to say when the death series was above its average. In the graph, the compensated series was translated upward by 2% to avoid a superposition with the initial series in all intervals in which there is no compensation. Sources: Same website as for Fig. 4.

troughs occurring around September are replaced by fairly shallow troughs. This corresponds to a reduction of the standard deviation from 5.11 for the initial series to 2.90 for the compensated series.

The previous procedure was quite successful but now we must try to explain unexpected birth troughs which appeared in recent decades (see Fig. 7b).

Beyond the death-nobirth effect

Before coming to the "new" birth troughs we show in Fig. 7a a decade where they are well accounted for by the death-nobirth effect. In Fig. 7b the birth troughs follow a different pattern which displays the following features.

• The death peaks are still in January (indicated by the vertical lines in red) with the only exception of the heat wave of the summer of 2003. Here, in most of the cases the time lags between death spikes and birth troughs are longer than 9 months. However the death spike associated with the heat wave of August 2003 is followed by a birth trough 9 months later (i.e. in May 2004). In addition, we observe a time lag close to 9 months between birth peaks and birth troughs (see the U shaped line in cyan) for which an explanation is suggested in the text. Source: Same website as for Fig. 4.

• Except in 1997 and 1999 he birth troughs are no longer in October (indicated by the vertical lines in magenta) but rather in December.

• Between the birth peaks and the birth troughs there is nearly a 9 month interval.

• The death-nobirth effect did not really disappear but the October dips are shadowed by the presence of lower dips. What can be the origin of the new troughs? By comparison with Fig. 7a the births dips appear more disorderly; for instance often they display double dips. Similarly the birth peaks appear also more disorderly in the sense of a single spike there are often double spikes. In other words, there seems to be a connection between birth spikes and birth dips. This observation leads to a possible explanation.

It is self evident that for a couple, conception at time t excludes conception during the following 9 months. If all couples (of child bearing age) would conceive in October, then during the subsequent 9 months there would not be a single conception. As a result, one would see a huge birth peak in October + 9 months=July, followed by a 9-month interval without any birth. In the real world where couples are not highly synchronized one may expect a conception peak to be followed by a conception trough during the following nine months. For the sake of brevity this effect will be referred to as the single conception effect. Needless to say, as it is a basic biological effect the single conception effect will always be present. However, it will be clearly visible only when the birth peaks are high and fairly narrow. With broad and rather low birth peaks the effect will be "diluted" and hardly visible.

If our assumptions are correct the birth pattern will be a superposition of two effects: the death-nobirth effect and the single conception effect but the weight of each effect may change in the course of time. When the January influenza death spike becomes lower as a result of more persons being vaccinated, then the weight of the deathnobirth effect will dwindle and, as a result, the single conception effect will become more apparent.

We conjecture this is what happens here. From the 1950s to the 1990s the average amplitude of the January death spikes decreases from 1.5 to 1.15 and as a result the death-nobirth effect is eclipsed by the single conception effect. This explanation should be tested more thoroughly in a subsequent paper.

Influence of the religious environment on conceptions

Some religious precepts rule out sexual intercourse during special periods. For instance, during the month of Ramadan, sexual activity is not allowed during the day but is permitted during the night. The Orthodox religion advises against sexual activity during Lent which in Orthodox countries is the 7-week period (49 days) preceding Easter.

In present day Catholic and protestant religions there are no specific precepts on this point. However, although the scriptures do not say anything in this respect in past centuries the Catholic Church prohibited sexual activity not only during Lent (46 days before Easter) but also in many other periods, e.g. 40 days before Christmas, 40 days after giving birth, or on all the nights from Saturday to Sunday. Conversely, on some specific days, sexual intercourse is encouraged. This is the case for Saturday (or more precisely Shabbat which starts already on Friday night) in the Jewish religion.

With respect to birth rates this leads to the following questions.

(1) Is there empirical evidence showing that some religious precepts affect the pattern of births?

(2) Is it possible to document the impact of religious rules over the course of time and across countries? Empirical evidence: the Fire Horse years in Japan An easy way to show the impact of beliefs on births may be found in Japan. No monthly or daily data are requested, annual birth data are sufficient.

Number of births (millions)

Japan China Fig. 8 Annual number of births in China and Japan. The sharp reduction of births in Japan in 1966 can be related to the fact that in the Chinese astrological calendar it was a "Fire Horse" year. In Japan there is a kind of prejudice with respect to girls born in a Fire Horse year. The fact that this belief exists in Japan but not in China suggests that it has little to do with Buddhism but is rather a local tradition. Sources: Japan: Historical Statistics of Japan. China: Ren kou (1988, p. 268) The number of people (both males and females) born in 1966 is much smaller than in 1965 or in 1967. The difference is of the order of 25%. It turns out that 1966 was a Hinoeuma year, which means a Fire Horse year in the Chinese calendar6 . According to popular belief in Japan girls born in that year grow up to be "Fire Horse women" who are reputed to be headstrong and to bring bad luck to their families. In 1966, as a baby's sex could not be reliably identified before birth, there was a sharp fall in birth rate, partly due to a fall in conceptions and partly to an increase in abortions. According to the Chinese calendar, Fire Horse years occur every 60 years; thus, the three previous ones were in 1906, 1846 and 1786. Are the effects in those years similar to the one in 1966? For 1786 there are no data but for the three subsequent years the answers are given in Table 3. Fall in annual births (%) -11% -11% -24%

Notes: The percentages refer to the differences between the year under consideration and the mean of the 9 other years in the same decade. Was there was a similar effect in China? In spite of the fact that 1966 has less births than adjacent years, there is no sudden drop, which means that the trough around 1966 may be due to other reasons. Note also that in Japan 1966 had a larger male/female birth ratio than adjacent years; the jump from 1.053 to 1.076 may seem small but is significant because in normal years the male/female ratio is very stable.

In China the male/female ratio does not show any special behavior in 19667 .

Empirical evidence: the case of Lent in Romania

There are two good reasons for examining the case of Romania.

• As mentioned above, in contrast with the Catholic or Protestant religions, the Orthodox religion has clearly defined rules regarding sexual relations during Lent. In the census of 2011, 86% of the Romanian population identified as part of the Eastern Orthodox Church.

• Because the time of Lent changes from year to year it is a much better "marker" than a fixed religious day. Although a single observation may display a reduction in births it will not reveal the reason of the reduction. In contrast, if in successive years the birth reductions follow closely the Lent time-intervals, then there can be little doubt that Lent is the causal factor.

The effect of religious rules on conception during Lent was investigated in [START_REF] Herteliu | Effect of religious rules on time of conception in Romania from 1905 to 2001[END_REF]. For such a study one of the main difficulties is the fact that strictly speaking monthly birth data are not detailed enough to follow the shifts of Lent. One needs weekly or daily data. Fortunately, the Romanian Bureau of the Census asked all persons surveyed in the censuses of 1992 and 2002 for their birth day dates.

As some elderly persons may not remember their birthday with absolute accuracy, in order to limit this uncertainty the data set was limited to the persons born after 1905; thus, those were less than 87 year old in 1992. Because the size of a cohort decreases rapidly after 65 the dataset would include only small samples of persons born in the years immediately following 1905; that may result in fairly large statistical fluctuations.

For the whole period from 1930 to 2000 the reduction of conceptions during Lent was found to be on average equal to: 14.1% ± 1%. In contrast there was almost no reduction during Advent, the 40-day "Nativity Fast" that precedes Christmas.

Comparative perspective across countries

From the above discussion one would expect a smaller reduction effect in Catholic or Protestant countries than in Orthodox countries How can we test this conjecture?

As daily or weekly birth data are not be available for the early time periods that we wish to consider, we made the following simplified argument which will allow us to use monthly data.

(1) Firstly, our comparison will concern the time-interval 1931-1935 because this allows us to use the international monthly data given in [START_REF] Bunle | Le mouvement naturel de la population dans le monde de 1906 à 1936[END_REF].

(2) Secondly, we observe that in western countries Easter is on average on 8 April. This statement means that if we identify Easter with its daily location from 1 January on, then the average of these numbers over 30 (non-leap) years will be 98 days which corresponds to 8 April. Note that this is for Catholic and Protestant countries. For Orthodox Easter the average date is somewhat later, namely on 21 April. As in western countries Lent lasts 40 days, on average western Lent will last from 1 March to 8 April which means that if we identify Lent with March we will miss only 8/40 = 20%. For Orthodox Easter the situation is less favorable. As Lent lasts 49 days it will be from 1 March to 21 April and will therefore be divided more equally (4/7 against 3/7) between March and April. By identifying Lent with March we will miss 21/49 = 43%.

There is an additional factor which makes the approximation more acceptable. If the religious injunction is followed fairly closely one would expect a rebound of conceptions after Easter. Such a rebound can indeed be observed; for instance in Bulgaria conceptions surged from 672 in March to 983 (+46%) in April. In western countries a surge can also be observed but which is much smaller. For instance in France from March to April conceptions rise from 939 to 1042 (+11%).

In other words, in Orthodox countries any reduction between 1 and 21 April may be offset by a surge in conceptions in the rest of the month. Thus, even in Orthodox countries, March appears to be a better month to consider than April.

Under the previous assumptions the conception reduction due to Lent can be represented by the ratio8 : R defined as: (conceptions in March)/(average of the conceptions in February and April) By adding 9 months the conception ratio is then translated into a birth ratio. Source: Bunle (1954, p. 92-93).

Comparative perspective across time

One would expect a decrease of the Lent effect in the course of time. The data given in [START_REF] Herteliu | Effect of religious rules on time of conception in Romania from 1905 to 2001[END_REF] allow us to follow the reduction from about 1930 to 2000.

Needless to say, one would not be surprised to see a fairly sharp fall of the impact of religion during the 40 years of Communist regime. However, in order to assess the rapidity of the fall we should compare Romania with a non-Communist country. If we use again the approximation of Lent by the month of March, we can use French monthly birth data which cover a broad time interval ranging from 1650 to 2000. The comparison presented in Fig. 9 shows that in Romania the fall was about 6 times faster than in France.

Conclusion

In this paper we have given interpretations of several kinds of birth rate troughs.

• For daily rates the main factor is the difference between working days and days off. The results for France are based on monthly birth data whereas those for Romania are based on daily data which is why they have a larger dispersion. Actually, prior to the 19th century, the records give the dates of baptisms. In order to accept them as proxies of birth dates one must assume that baptism followed birth fairly closely. Overall the fall of the Romanian curve is 6 times steeper than the fall of the French curve. Sources: [START_REF] Herteliu | Effect of religious rules on time of conception in Romania from 1905 to 2001[END_REF], Houdaile (1979Houdaile ( , 1985)), Régnier-Loilier (2010b), [START_REF] Dupâquier | Le mouvement saisonnier des naissances en France de 1853 à 1973[END_REF].

• For monthly rates any mass mortality event brings about a birth rate trough 9 months later because surviving people (and especially young people who can have children) suffer adverse living conditions. This effect was seen to be much stronger for influenza than for heat waves.

• The influence on birth rates of the religious environment can be seen statistically during special events such as the Fire Horse year or the Lent period.

The hallmark of the rules that we have studied is that, to varying degrees, they are valid in all countries for which we were able to find reliable data. Although in present-day France the influence of religious beliefs is small the analysis of historical data permitted to measure their strength in former centuries. In short, this analysis provides an objective way to assess the beliefs of people and how they change in the course of time.

In Fig. 7b we have seen that when two effects are in competition the resulting birth pattern becomes fairly intricate. Note that in this figure the death peaks remained unchanged; any change in the death pattern in the course of time would make the interpretation even more difficult. Similarly, in the language parallel proposed in the introduction, if the script or the meaning of a hieroglyph changes in the course of time this will affect the meaning of all sentences in which it occurs.

We plan to extend this investigation further when new data become available. For instance, using monthly birth data for China it would be interesting to test if there are conception surges in the "Golden week" (1-7 October) or during the "Spring festival", that is to say the week centered on the Chinese New Year. It can also be observed that so far we have no interpretation of birth rate peaks.

Because the question of birth seasonality has been studied for almost two centuries it gives an excellent opportunity to examine how the methodology has changed in the course of time. The discussion given in Appendix A suggests that by renouncing the problem-oriented approach based on comparative analysis, present-day researchers put themselves in an unpromising condition.

Appendix A: A historical view of former studies It would be an impossible and useless task to try to summarize the immense literature on this topic. However, as our study brought us in contact with both historical and present-day papers it may be of interest to characterize the prevailing trend during the previous two centuries.

Monthly birth data series became available in European countries in the early 19th century and in the United States about one century later (due to the fact that it took some time to set up registration procedures in all states). Explaining the monthly distribution of births is a question which attracted the attention of researchers very early.

A pioneer study published in 1831

One of the first studies was published in 1831 by Louis René Villermé (1782Villermé ( -1863)), a French physician and a pioneer of social epidemiology. In several respects his article is typical of those published at that time and it is interesting to note how it differs from present-day publications.

(1) It is a long article of one hundred pages.

(2) It is a comparative study which uses data not only for France but as well from many other countries: Great Britain, Germany (or rather Bavaria, Prussia, Wurtemberg for at that time Germany was not unified), Italy, the Netherlands, Russia, Sweden. It also uses data from various time periods going back for some places (e.g. Paris) to the 17th century. This kind of comparative methodology remained in use during the entire 19th century. Durkheim's renowned study of suicide (1897) used exactly the same methodology which, needless to say, was borrowed from physics where experimenters also perform as many observations as possible of any new phenomenon.

(3) Interestingly, 30 pages of the paper are statistical tables which contain all the data used by the author. Several of these data sets are still useful, as for instance the monthly birth numbers for Belgium and the Netherlands from 1815 to 1826. The thanks expressed in the notes to the tables reveal that the author had personal contacts with several foreign scholars from whom he received many of his data.

(4) In line with other scholars of the 19th century (e.g. Alfred Espinas, 1878) Villermé does not hesitate to include in his study observations about animal species. This is of course in accordance with the perspective of comparative analysis. He observes for instance that for foxes, hares and wolves conception takes place between December and February which is fairly different from humans for whom the conception peak is in spring. Although nowadays animal "models" are also investigated it is not common to refer to them in demographic studies.

(5) Basically the paper states most of the general rules we currently know about. In other words, progress seems to have been fairly limited.

In summary, Villermé used a problem-oriented approach which leads him to collect such data that are the most appropriate for solving the questions under consideration.

Main characteristics of studies published in the past 50 years

The previous approach is in marked contrast with more recent studies. As illustrations we briefly discuss the following papers: (1) MacFarlane ( 1974), ( 2) [START_REF] Seiver | Trend and variation in the seasonality of US fertility, 1947-1976[END_REF], (3) [START_REF] Lam | Seasonality of births in human populations[END_REF][START_REF] Lam | Global patterns of seasonal variation in human fer-tility[END_REF], ( 4) [START_REF] Mcquillan | Does religion influence fertility?[END_REF] (5) [START_REF] Huber | Perinatal winter conditions after later reproductive performance in Romanian women: intra and intergenerational effects[END_REF],.

(1,2,3) and ( 4) are fairly broad papers whereas (5) focuses on a very specific point in a single country, namely Romania. In a problem-oriented research it is impossible to limit oneself to just one country for then it is impossible to know if the answer obtained in one case has any broad validity. Nonetheless. the one-case approach has become the norm nowadays.

Whereas (2) is also a single country investigation, (1,3) include data for several countries. However, they are not problem-oriented. As a matter of fact, they avoid raising any question. For instance, [START_REF] Lam | Global patterns of seasonal variation in human fer-tility[END_REF] show that between the monthly birth patterns of Britain and the US state of Georgia there is a 6 month time lag. For two regions that are fairly similar in many respects this is rather unexpected; however the authors do not try to explain it. They do not even mention that this observation would warrant a further investigation.

At first sight, (4) seems the most promising paper for its title has a strong connection with our investigation. Unfortunately, it is a review paper which contains almost no data whatsoever. It is a collection of qualitative statements which provide little real understanding. How can one discuss fertlity, basically a quantitative notion, without providing tables or graphs with data?

In summary, because they are not problem-oriented, the previous papers offer very little in terms of new understanding. They do not even give an "observational understanding" that is to say rules of broad validity which, thanks to their predictive power, allow further progress. Naturally, we do not claim that the previous papers represent the whole literature on this topic, but at least they show a trend.

Fig

  Fig. 1a,b Daily births in France in 1969 and 2015. The blue labels give calendar months whereas the red labels in italic indicate the corresponding conception months. Thus a conception in April would result in a birth in January of the following year. The same notation is used throughout the present paper. For the purpose of making relative variations more readable the series was normalized by making its average equal to 100. The substantial increase in dispersion has a simple explanation given below in the "Days off" section. Source: Répartition quotidienne des naissances vivantes, France métropolitaine. Institut National de la Statistique et des Études Économiques (INSEE). [Daily distribution of live births. Overseas territories are excluded. National Institute for Statistics and Economic Studies]. Available on Internet under the title: "Les naissances en 2016. Tableaux de séries longues".

  Fig. 2b Average of daily births in France over 35 years starting with first Monday of January. The data cover 1969-2014 (leap years were excluded). The downward spikes correspond to Saturdays and Sundays. It is for the sake of clarity that the graph was limited to 200 days. The rest of the 52 weeks would not show any new effect. Source: Same as for Fig.1.

Fig. 3

 3 Fig. 3 Daily births inFrance in 1969 after anomalous non-working day data were smoothed out. Data for days off were replaced by the average of the births on the two adjacent working days. This resulted in a reduction of the standard deviation from 6.12 in Fig.1to 5.29 here. Source: Same as for Fig.1.

  Fig.4bThe effect on births of influenza outbreaks in France between 1948 and 1974. Relationship between monthly excess deaths in January (with respect to the annual average and reduction in births in October of the same year. The coefficient of correlation is -0.89 and the equation of the regression line is given in the text. The three squares in red correspond to heat waves in France (see the text). Sources: Births: Same website as for Fig.1. Deaths: "Nombre moyen de décès par jour selon le mois depuis 1946 pour la France métropolitaine", INSEE. [Monthly numbers of deaths divided by the length of each month from 1946 to 2016]. Available on Internet under the title: "Les décès en 2016. Tableaux de séries longues". Effect of heat waves:

Fig. 5

 5 Fig. 5 Removal of the death-nobirth effect, for daily births, France 1968-1969. This figure explains why in 1969 the removal of the death-nobirth effect has only a small incidence. Sources: Same website as for Fig. 4.

Fig. 7a

 7a Fig. 7a Illustration of the death-nobirth effect. The data are for France and correspond to monthly deaths and births divided by the number of days of each month. All red lines of the death graph indicate January and all magenta lines of the birth graph indicate October (i.e. January+9months). It can be seen that the birth troughs follow the death spikes with the expected time lag of 9 months. However, as shown in the next figure, in the decades after 1990 the birth pattern becomes more complicated. Source: Same website as for Fig. 4.

  Fig. 7bApparition of a new effect which produces birth troughs not related to death spikes. Here, in most of the cases the time lags between death spikes and birth troughs are longer than 9 months. However the death spike associated with the heat wave of August 2003 is followed by a birth trough 9 months later (i.e. in May 2004). In addition, we observe a time lag close to 9 months between birth peaks and birth troughs (see the U shaped line in cyan) for which an explanation is suggested in the text. Source: Same website as for Fig.4.

Fig. 9

 9 Fig.9Reduction of conceptions during Lent in France and Romania. The results for France are based on monthly birth data whereas those for Romania are based on daily data which is why they have a larger dispersion. Actually, prior to the 19th century, the records give the dates of baptisms. In order to accept them as proxies of birth dates one must assume that baptism followed birth fairly closely. Overall the fall of the Romanian curve is 6 times steeper than the fall of the French curve. Sources:[START_REF] Herteliu | Effect of religious rules on time of conception in Romania from 1905 to 2001[END_REF],Houdaile (1979Houdaile ( , 1985)),Régnier-Loilier (2010b),[START_REF] Dupâquier | Le mouvement saisonnier des naissances en France de 1853 à 1973[END_REF].

Table 1 : Saturday and Sunday reductions in birth numbers with respect to working days.

 1

	Year Saturday Sunday
	1969 -1.2% -6.8% 2016 -15% -19%

Table 2 : Age distribution of the victims of heat waves and influenza respectively.

 2 Notes: The heat wave data are for France and are an average of 6 heat waves, namely:

	All ages	All ages	Age	Age	Age
	0 -100	0 -100 (per million pop)	0 -35 35 -75 75 -100

Table 3 : Birth reductions in the last three Fire Horse years in Japan

 3

	1846	1906	1966

  The high reduction in 1966 is due to two factors. (i) Contrary to what happened in 1846 and 1906, in 1966 the male-female sex ratio at birth increased only by 1.3% instead of 20% in 1846 and 4.3% in 1906; this means almost no infanticide in 1966. (ii) In 1946 there was a sharp fall in births which generated a small cohort; then, in 1966 this cohort reached the age of 20 at which point it started to have children. Source: Roehner (2007, p. 42).

Table 4 : Conception reduction during Lent, 1931-1935

 4 Notes: R equal to 1 would mean that March which represents Lent has same conceptions as the two adjacent months that is to say no reduction due to Lent. The monthly numbers given in the table are an index which sums up to 1,000 for one year. Bulgaria and Romania are Orthodox countries, France and Spain are Catholic countries, Norway and Sweden are Protestant countries.

	Conceptions in: Feb. (f ) Mar. (m) Apr. (a) R = m/[(f + a)/2]
	Bulgaria	891	688	983	0.73
	Romania	911	686	919	0.75
	France	985	930	1042	0.95
	Spain	911	947	1096	0.94
	Norway	857	899	992	0.97
	Sweden	904	941	979	1.00

Incidentally, this is the advice already given by Emile Durkheim, one of the founding fathers of sociology, in his book

See: https://en.wikipedia.org/wiki/Jean-Fran%C3%A7ois Champollion.

The sources are: "Quick facts about labor induction", August 2016; Centers for Disease Control: Recent declines in induction of labor by gestational age, June 2014.

The error bars are for a confidence level of 0.95.

More precisely we have taken a pregnancy duration of 267 days as explained inRichmond et al. (2018a).

Because the Chinese New Year occurs in late January, the Fire Horse Year does not exactly coincide with 1966; in fact, it started on 21 January 1966 and ended on 8 February 1967.

Although a Hong Kong -Canadian movie produced in 2005 and entitled "Eve and the Fire Horse" suggests that in 1975 the "Fire Horse" belief was still present in China, this view was not confirmed by discussions with Chinese colleagues; as a matter of fact, they had no knowledge of it.

Because of the "surge effect" the conceptions of April will be higher than they would be without such an effect. As a result, R will be lower than it should be.