<figure><img src="image_1.png" /></figure>

HAL
open science

# Modeling of Nomex honeycomb cores, linear and nonlinear behaviours

Laurent Gornet, Steven Marguet, Gilles Marckmann

► To cite this version:

Laurent Gornet, Steven Marguet, Gilles Marckmann. Modeling of Nomex honeycomb cores, linear and nonlinear behaviours. Mechanics of Advanced Materials and Structures, 2007, 14 (8), pp.589-601. 10.1080/15376490701675370 . hal-01007067v2

HAL Id: hal-01007067

https://hal.science/hal-01007067v2

Submitted on 14 Dec 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Modeling of Nomex® Honeycomb Cores, Linear and Nonlinear Behaviors

L. Gornet, S. Marguet, and G. Marckmann

Ecole Centrale de Nantes, France

The purpose of this study is to develop tools dedicated to the design of sandwich panels involving Nomex® honeycomb cores. Special attention is paid to the ability to perform full three dimensional calculations up to failure of such structures. In the first part, the determination of effective elastic properties of Nomex® honeycomb cores is carried out thanks to strain based periodic homogenization technique. Using an equivalence in energy between a real honeycomb and a fictitious continuous medium, it becomes possible to evaluate the elastic behavior of Nomex® cores, starting from the knowledge of the behavior of the constitutive paper. Next, relying on experimental observations, the strengths of Nomex® honeycomb cores are evaluated with a linear Euler's buckling analysis. Results are compared with data coming from manufacturers, and give satisfaction. In order to carry out these two first studies, the NidaCore software has been developed using the finite element code Cast3M from CEA. The last part deals with the modeling of the nonlinear compressive response of Nomex® cores. A model based on the thermodynamics of irreversible process is proposed and the identification technique detailed. Good agreement between experimental data and computed values is obtained.

**Keywords** Nomex® honeycomb core, homogenization, failure, nonlinear behavior

## 1. INTRODUCTION

Composite sandwich panels are extensively used in lightweight constructions. There excellent stiffness to weight ratio make them particularly interesting for many applications such as aircraft and aerospace structures or oceanic race sailing boats. In the latter case, typical sandwich panels are made of two carbon fibers and epoxy resin skins, separated with a Nomex® honeycomb core as illustrated in Figure 1. The understanding of the mechanical behavior of the constituents up

to failure is extremely important to ensure safe design of the structures.

Several damage mechanisms contribute to sandwich structure failure. Matrix cracking often appears at first in laminates. Then, the micro-cracks coalesce and initiate delaminations and large scale crack propagation. Interfaces between core and laminates debond or core crushing lead the sandwich to failure [1]. In the context of finite element model, taking into account all these phenomena while keeping an admissible computational cost is not a simple matter. A classical simplifying assumption is to consider brittle linear elastic mechanical behaviors for both the skins and the core. Associated with anisotropic stress failure criteria such as Tsaï-Hill, Tsaï-Wu or Hoffman [2, 3], it then becomes possible to make calculations on the whole structure. As an illustration, Figure 2, presents a mesh of the multihull race sailing boat code B1 which is now well known as Orange II. The related finite element model involves multilayered thick shell elements which allow us to deal with both the skins and the core. Each layer is defined in terms of its orientation, its elastic moduli and a convenient stress failure criteria. With the assumption that no delamination occurs, a first design is then carried out. But in front of hard loadings, which can be wave impacts, slammings or violent swelling of veils, for example, the description of the damage phenomena has to be more precise. It results in more complicated models which involve more parameters can lead to increasing computational costs.

The aim of this study is to propose tools dedicated to the calculation of multihull race sailing boats. The focus is on the mechanical behavior of Nomex® honeycomb core. In the first section, strain based periodic homogenization technique is employed to determine the elastic constants that produce the linear part of the mechanical behavior of the honeycombs. Next, starting from experimental observations, strengths of Nomex® honeycomb core are evaluated thanks to Euler's buckling theory. Comparisons between finite element analysis (performed within Cast3M, the finite element code from CEA) and data coming from manufacturers are presented. In the third section, a unidirectional model dedicated to the nonlinear compressive behavior of Nomex® honeycomb core is proposed. The identification technique is detailed and discussion is made on the results and on the extension of the model in full three dimensions.
<figure><img src="image_3.png" /><figcaption>FIG. 1. Composite sandwich panel used for nautic construction: carbon fibres and epoxy matrix laminated skins separated by a Nomex® honeycomb core.</figcaption></figure>

## 2. ELASTIC PROPERTIES OF NOMEX® HONEYCOMB CORES

This section focuses on the method developed to compute the effective elastic properties of Nomex® honeycomb cores. The motivation for this study comes from the fact that the only available data from manufacturers of Nomex® honeycomb cores [4, 5] are relative to the out-of-plane behavior which involves: the compressive elasticity modulus, the out-of-plane shear modulus and the corresponding strengths. These data are obviously insufficient to perform three dimensional calculations. This paper is in the continuity of previous studies about the plane behavior

of Nomex® honeycomb cores [6, 7] and on other ones, dealing with cellular solids [8–10].

### 2.1. Strain Based Periodic Homogenization

Since modeling the true geometry of Nomex® honeycomb cores on the whole finite element model is unrealistic for reason of computational cost, another approach has to be implemented. Strain based periodic homogenization technique is a powerful tool to model the periodic structures that are honeycombs. The basic idea is to replace a periodic medium of complex geometry by a continuous medium with an equivalent mechanical behavior.

<figure><img src="image_4.png" /><figcaption>FIG. 2. Finite element mesh of code B1 multihull race sailing boat (Multiplast).</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>FIG. 3. Hexagonal and over-expanded unit cells.</figcaption></figure>

The first step of the approach is to define a Representative Volume Element (RVE) which can, by series of translations, reproduce the whole discontinuous geometry [11]. In the case of this study, two geometries presented in Figure 3 are taken into account: hexagonal and over-expanded rectangular unit cells of honeycombs.

By convention, and in all that follows, the $e_1$ (strong) and $e_2$ (weak) axis are in plane axis whereas $e_3$ is the out-of-plane axis. O is the origin of basis $\Re = (O, e_1, e_2, e_3)$. Starting from the center of a unit cell, two orthogonal plans of symmetry can be defined: $(O, e_1, e_3)$ and $(O, e_2, e_3)$. As a result, Nomex® honeycomb cores exhibit an orthotropic mechanical behavior.

The link between the micro scale of the RVE and the macro scale of the equivalent continuous medium is done thanks to Hill Mandel's theorem. It postulates the equality between the mean elastic energy on the RVE and the energy in the equivalent medium, which reads:

$$ \underline{\underline{\Sigma}} : \underline{\underline{E}} = \frac{1}{V} \int_V \underline{\underline{\sigma(x)}} : \underline{\underline{\epsilon(x)}} dV \quad (1) $$

where:

* $\underline{\underline{\sigma(x)}}$: is the stress at point $(x)$ in the RVE;

* $\underline{\underline{\epsilon(x)}}$: is the strain measure at the same point;

* $\underline{\underline{V}}$: the volume of the RVE ($dV = dx_1dx_2dx_3$ is the elementary volume);

* $\underline{\underline{\Sigma}}$: is the homogeneous macro stress of the equivalent continuous medium;

* $\underline{\underline{E}}$: is the homogeneous macro strain loading.

By introducing the mechanical behaviors at the micro scale (known behavior of Nomex® paper) and at the macro scale (unknown behavior of the equivalent medium), Eq. (1) deals to:

$$ \left\{ \begin{array}{l} \underline{\underline{\sigma(x)}} = \underline{\underline{K}}_{\underline{\underline{micro}}} : \underline{\underline{\epsilon(x)}} = \underline{\underline{K}}_{\underline{\underline{micro}}} : \underline{\underline{\epsilon(x)}} \\ \\ \underline{\underline{\Sigma}} = \underline{\underline{K}}_{\underline{\underline{MACRO}}} : \underline{\underline{E}} \\ \\ \underline{\underline{E}} : \underline{\underline{K}}_{\underline{\underline{MACRO}}} : \underline{\underline{E}} = \frac{1}{V} \int_V \underline{\underline{\epsilon(x)}} : \underline{\underline{K}}_{\underline{\underline{micro}}} : \underline{\underline{\epsilon(x)}} dV \end{array} \right. \quad (2) $$

with:

* $\underline{\underline{K}}_{\underline{\underline{micro}}}$: elasticity tensor representative of the mechanical behavior of Nomex® paper;

* $\underline{\underline{K}}_{\underline{\underline{MACRO}}}$: elasticity tensor of the mechanical behavior of equivalent medium.

A macroscopic strain field $\underline{\underline{E}}_{\underline{\underline{imposed}}}$ is then imposed at the micro level of the RVE by prescribing a displacements field: $\underline{\underline{u(x)}} = \underline{\underline{E}}_{\underline{\underline{imposed}}}, \forall x \in V$. On the boundary of the RVE, periodicity conditions also have to be taken into account. If $\partial V$ is the boundary of the RVE, it can be decomposed into $\partial V_\sigma$ which is the boundary with imposed stress and $\partial V_\epsilon$ the boundary with imposed strain. As the study is strain based, no stress is imposed on $\partial V_\sigma$, but on $\partial V_\epsilon$, periodicity conditions are introduced by writing that the field of displacements is the sum of two terms:

* the field of displacements due to the imposed macro strain field $\underline{\underline{E}}_{\underline{\underline{imposed}}}$;

* a periodic field of displacement that links the opposite points on the boundary of the RVE.

If A and B are two opposite points on the boundary of the RVE, periodicity relations are:

$$ \left\{ \begin{array}{l} \underline{\underline{u(A)}} = \underline{\underline{E}}_{\underline{\underline{imposed}}} \underline{\underline{x}}_A + \underline{\underline{u_{periodicity}(x_A)}} \\ \\ \underline{\underline{u(B)}} = \underline{\underline{E}}_{\underline{\underline{imposed}}} \underline{\underline{x}}_B + \underline{\underline{u_{periodicity}(x_B)}} \\ \\ \underline{\underline{u(x_A)}} - \underline{\underline{u(x_B)}} = \underline{\underline{E}}_{\underline{\underline{imposed}}} (\underline{\underline{x}}_A - \underline{\underline{x}}_B) \end{array} \right. \quad (3) $$

since $\underline{\underline{u_{periodicity}(x_A)}} = \underline{\underline{u_{periodicity}(x_B)}}$. To summarize the previous considerations, the elementary mechanical problem to solve is now introduced within Eq. (4):

$$ \left\{ \begin{array}{ll} \underline{\underline{div(\sigma(x))}} = \underline{\underline{0}} & \text{in } V \quad (\text{equilibrium}) \\ \underline{\underline{\sigma(x)}} = \underline{\underline{K}}_{\underline{\underline{micro}}} : \underline{\underline{\epsilon(x)}} & \text{in } V \quad (\text{behavior}) \\ \underline{\underline{\epsilon(x)}} = \underline{\underline{grad(\underline{u(x)})}} & \text{in } V \quad (\text{kinematics}) \\ \underline{\underline{\sigma(x)}} = \underline{\underline{0}} & \text{on } \partial V_\sigma \quad (\text{boundary conditions}) \\ \underline{\underline{u(x)}} = \underline{\underline{E}}_{\underline{\underline{imposed}}} \underline{\underline{x}} & \text{on } \partial V_\epsilon \quad (\text{boundary conditions}) \\ + \underline{\underline{u_{periodicity}(x)}} & \end{array} \right. \quad (4) $$

It has to be noted that, in finite element analysis, the last equation of Eq. (4) is introduced with linear relations between displacement degrees of freedom of the opposite nodes located on $\partial V_\epsilon$. Last step is the choice of the macroscopic strain field $\underline{\underline{E}}_{\underline{\underline{imposed}}}$ by reminding us that due to symmetries, Nomex® honeycomb cores present an orthotropic mechanical behavior, the
elasticity tensor is:

$$
\begin{equation}
\left\{
\begin{aligned}
\underline{E} &= \underline{K}^{-1}_{\equiv MACRO} : \underline{\Sigma} \equiv \text{or, with engineering notations:} \\
\left[
\begin{array}{c}
E_{11} \\
E_{22} \\
E_{33} \\
\sqrt{2}E_{12} \\
\sqrt{2}E_{23} \\
\sqrt{2}E_{31}
\end{array}
\right]
&=&
\left[
\begin{array}{ccccccc}
\frac{1}{E_1} & \frac{-\nu_{12}}{E_1} & \frac{-\nu_{13}}{E_1} & 0 & 0 & 0 \\
\frac{-\nu_{21}}{E_2} & \frac{1}{E_2} & \frac{-\nu_{23}}{E_2} & 0 & 0 & 0 \\
\frac{-\nu_{31}}{E_3} & \frac{-\nu_{32}}{E_3} & \frac{1}{E_3} & 0 & 0 & 0 \\
0 & 0 & 0 & \frac{1}{2G_{12}} & 0 & 0 \\
0 & 0 & 0 & 0 & \frac{1}{2G_{23}} & 0 \\
0 & 0 & 0 & 0 & 0 & \frac{1}{2G_{31}}
\end{array}
\right] \\
& & \times &
\left[
\begin{array}{c}
\Sigma_{11} \\
\Sigma_{22} \\
\Sigma_{33} \\
\sqrt{2}\Sigma_{12} \\
\sqrt{2}\Sigma_{23} \\
\sqrt{2}\Sigma_{31}
\end{array}
\right]
\end{aligned}
\right.
\tag{5}
\end{equation}
$$

Nine macroscopic loadings are then applied to the RVE in order to compute the macroscopic behavior of the equivalent continuous medium. They are presented in Figures 4, 5 and 6 for over-expanded rectangular cells.

The stored elastic energies of these finite element models are then evaluated and, by using Hill Mandell's theorem (Eq. 2), the characteristics of the equivalent continuous medium are computed.

## 2.2. Results, Effective Elastic Properties

The mechanical properties of numerous geometries of Nomex® cores have been derived from the method presented above. The microscopic mechanical behavior of Nomex® paper has been evaluated by rescaling a given out-of-plane shear modulus determined by periodic homogenization technique with the

value given by manufacturers [4, 5]. By this way, with only one set of data from the manufacturer, it is possible to estimate the whole three dimensional mechanical behavior of a Nomex® honeycomb core. Table 1 presents the elastic properties for three kinds of over-expanded rectangular honeycombs of cell size 4.8 mm, sheet thickness 51 µm and densities: 29 kg.m⁻³, 48 kg.m⁻³ and 64 kg.m⁻³. Very good agreement is found with available data concerning the out-of-plane properties.

Next, the evolutions of the elastic moduli in function of the core density are plotted on Figure 7 for two types of honeycombs which are: 4.8 mm over-expanded rectangular and hexagonal cells. An increasing density of the core appears to have significant effect on the stiffness of the structure.

With this first study, relying on the strain based periodic homogenization technique, it becomes possible to determine the full three dimensional mechanical behavior of any Nomex® honeycomb core, any in the sense of any geometry, any cell size, any density of Nomex® sheet and any thickness. This aspect is particularly interesting since available catalogues only provide the out-of-plane moduli of these materials for a given thickness of 12.7 mm.

A pragmatic approach to evaluating the strengths of honeycomb cores is presented in next section.

# 3. STRENGTHS OF NOMEX® HONEYCOMB CORES

## 3.1. Experimental Observations

In the approach presented here, the determination of the strengths of Nomex® honeycomb cores relies on the assumption that the buckling phenomenon drives the failure of the structural materials. To justify this strong hypothesis, several tests are performed on out-of-plane compressive and shear specimens which are presented on Figure 8 with a photograph of the mean of test.

Nomex® honeycomb cores used during this experimental campaign are made of over-expanded rectangular cells of size 6.4 mm, density 64 kg.m⁻³ and thickness 40 mm. Both the compressive and simple shear specimens use TR-50 carbon fibres and epoxy resin skins of stacking sequence: [0 90]s for a

<figure><img src="image_6.png" /><figcaption>FIG. 4. Over-expanded rectangular unit cell submitted to simple stretches: $E_{11}$, $E_{22}$, and $E_{33}$.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>FIG. 5. Over-expanded rectangular unit cell submitted to simple shear: E₁₂, E₂₃, and E₃₁.</figcaption></figure>

total thickness of 1.28 mm. These skins are stiff enough to enable the characterization of the mechanical behavior of Nomex® cores only. Their main advantage is to impose boundary conditions similar to the ones observed in real constructions between the core and the skins. For a compressive specimen, the block of Nomex® has a 50 mm length for a width of 50 mm and a thickness of 40 mm. In the case of simple shear specimen, four blocks of dimensions: 70 by 50 by 40 mm are used. This choice ensures a sufficient number of cells in the loading direction. The tests are carried out on a MTS hydraulic tensile testing machine. The displacement is measured by the magnetic captor that drives the test, and the force is obtained with a 25 kN cell force. In the case of compressive tests the imposed velocity is 0.04 mm.s⁻¹ equivalent to a strain rate of 0.001 s⁻¹. For simple shear tests, the imposed velocity is 0.16 mm.s⁻¹ which also leads to an average strain rate of 0.001 s⁻¹.

Figure (9) plots the crushing stress in function of crushing strain response of two over-expanded rectangular Nomex® honeycomb cores.

The photograph in Figure 9 has been taken just after the appearance of a global mode of instability. It shows the formation of a macroscopic buckling ply on the whole block of

Nomex® core and is directly related to the brutal loss of stiffness visible on the response stress curve. The link between the failure of Nomex® honeycomb core and the buckling of its walls is clearly demonstrated. Figure 10 does the same as Figure 9 but for simple shear specimens. W (or *e*₁) and L (or *e*₂) stress versus strain curves are plotted on the graph. After a slightly nonlinear reversible response, a subit loss of stiffness is observed. For W direction, this decrease in stress appears in a range of strain evolving from approximately 0.01 and 0.012 m.m⁻¹. The photograph on the upper left draws the lower right block of Nomex® just before the brittle failure (presented on the other photograph). At this stage, no macro crack is visible and moreover, a global bending mode appears on the whole block. This leads us to think that the loss of stiffness is not related to a brittle failure phenomenon but to the apparition of an instability which is the first buckling mode under simple shear loadings. The conclusion is that, for compressive loadings, in simple shear the failure of the structural material is due to the buckling phenomenon that is prior to the brittle failure of the block.

In the next part, and starting from this remark, a method based on Euler's buckling analysis is presented. Critical loads

<figure><img src="image_8.png" /><figcaption>FIG. 6. Over-expanded rectangular unit cell submitted to biaxial stretches: E₁₁ & E₂₂, E₂₂ & E₃₃ and E₃₃ & E₁₁.</figcaption></figure>
<table><caption>TABLE 1
Computed elastic properties of Nomex® honeycomb cores</caption><thead><tr><th rowspan="3"></th><th colspan="4">ECA-R 4.8 mm 29 kg.m<sup>-3</sup> 51 µm</th><th colspan="4">ECA-R 4.8 mm 48 kg.m<sup>-3</sup> 51 µm</th></tr><tr><th colspan="2">FE analysis</th><th colspan="2">Euro-composite</th><th colspan="2">FE analysis</th><th colspan="2">Euro-composite</th></tr><tr><th>Minimal</th><th>Typical</th><th>Minimal</th><th>Typical</th><th>Minimal</th><th>Typical</th><th>Minimal</th><th>Typical</th></tr></thead><tbody><tr><td>ν<sub>12</sub></td><td colspan="2">0.2606</td><td colspan="2"></td><td colspan="2">0.2606</td><td colspan="2"></td></tr><tr><td>ν<sub>23</sub></td><td colspan="2">0.0219</td><td colspan="2"></td><td colspan="2">0.0219</td><td colspan="2"></td></tr><tr><td>ν<sub>31</sub></td><td colspan="2">0.224</td><td colspan="2"></td><td colspan="2">0.224</td><td colspan="2"></td></tr><tr><td>E<sub>1</sub></td><td>44.5</td><td>72.5</td><td colspan="2"></td><td>101.1</td><td>128.7</td><td colspan="2"></td></tr><tr><td>E<sub>2</sub></td><td>4.4</td><td>7.1</td><td colspan="2"></td><td>9.9</td><td>12.6</td><td colspan="2"></td></tr><tr><td>E<sub>3</sub></td><td>79.4</td><td>129.5</td><td colspan="2"></td><td>180.5</td><td>229.8</td><td colspan="2"></td></tr><tr><td>G<sub>12</sub></td><td>0.5</td><td>0.9</td><td colspan="2"></td><td>1.2</td><td>1.6</td><td colspan="2"></td></tr><tr><td>G<sub>23</sub></td><td><strong>8.6</strong></td><td>14.1</td><td><strong>9</strong></td><td>14</td><td><strong>19.6</strong></td><td>25</td><td><strong>18</strong></td><td>24</td></tr><tr><td>G<sub>31</sub></td><td><strong>14.6</strong></td><td>23.9</td><td><strong>14</strong></td><td>24</td><td><strong>33.3</strong></td><td>42.4</td><td><strong>36</strong></td><td>44</td></tr></tbody></table>

<figure><img src="image_9.png" /><figcaption>FIG. 7. Evolution of elastic properties with the density.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>FIG. 8. Compressive and simple shear specimens.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>FIG. 10. Over-expanded rectangular honeycomb cores submitted to simple shear tests.</figcaption></figure>

obtained by homogenization technique, it is possible to determine the ultimate stress $\Sigma^{c,i}$ corresponding to the ith mode of buckling:

$$ \Sigma^{c,i} = \underset{\equiv MACRO}{K} : \underset{\equiv MACRO}{E}^{c,i} \qquad (11) $$

where:

* $\Sigma^{c,i}$ has only one non null value $\Sigma^{c,i}$ in adequation with the imposed strain field $\underset{\equiv MACRO}{E}^{c,i}$.

The failure stress limits are evaluated on the first buckling mode (the one of lowest energy).

### 3.3. Results, Stress Failure Criterion

To determine the full three dimensional stress failure limits of Nomex® honeycomb cores, the method presented below has been applied with the following loadings:

* $\underset{\equiv MACRO}{E}^{unit}$ such as $E_{ij} = 1$ if $i \neq j$, 0 elsewhere, with $(i; j) = \{1; 2; 3\}^2$ (pure shear loadings);

* $\underset{\equiv MACRO}{E}^{unit}$ such as $\underset{\equiv MACRO}{E}^{unit} = K^{-1} : \Sigma^{unit}$ with $\Sigma_{ij} = 1$ if $i = j$, 0 elsewhere and $(i; j) = \{1; 2; 3\}^2$ (equivalent to uniaxial pure compressive tests in direction $\mathbf{e}_{i=j}$).

The deformed shapes of the RVE submitted to out-of-plane shearing are plotted in Figure 11. A comparison between computed ultimate stress values and data available from catalogue of manufacturer Euro-composite [5] is given in Table 2 for 4.8 mm cell size, 29 kg.m⁻³ density and 51 µm of thickness paper honeycombs.

Satisfactory correlation is founded on the out-of-plane ultimate shear stress values. A more important difference appears in compressive ultimated stress, difference which might be due to

the boundary conditions on the skins of the RVE that can't prevent rotations of nodes and results in a lowest energy buckling mode (under evaluation of ultimate stress).

### 4. NONLINEAR COMPRESSIVE BEHAVIOR

As mentioned in introduction, hard loadings often occur during trans-oceanic races, inducing numerous abandons due to structural failure. To safely design the boats in front of these loadings, it appears necessary to understand the physics of the failure process and to be able to reproduce them. It means understanding how skins and core react not only in their linear domain but also in their nonlinear part. For example, when a wave impacts the hull of a boat, the energy is transmitted to the inside skin

<figure><img src="image_12.png" /><figcaption>FIG. 11. Deformed shapes of over-expanded rectangular RVE submitted to simple out-of-plane shear</figcaption></figure>
<table><caption>TABLE 2
Ultimate stress of over-expanded rectangular honeycomb cores.</caption><thead><tr><th>ECA-R</th><th>FE (min)</th><th>Tests (min)</th><th>FE (typ)</th><th>Tests (typ)</th></tr></thead><tbody><tr><td>σ<sup>c</sup><sub>11</sub></td><td>0.185</td><td>—</td><td>0.31</td><td>—</td></tr><tr><td>σ<sup>c</sup><sub>22</sub></td><td>0.282</td><td>—</td><td>0.48</td><td>—</td></tr><tr><td>σ<sup>c</sup><sub>33</sub></td><td><strong>0.39</strong></td><td>0.6</td><td><strong>0.62</strong></td><td>0.85</td></tr><tr><td>σ<sup>c</sup><sub>12</sub></td><td>0.72</td><td>—</td><td>1.2</td><td>—</td></tr><tr><td>σ<sup>c</sup><sub>23</sub></td><td><strong>0.25</strong></td><td>0.31</td><td><strong>0.415</strong></td><td>0.42</td></tr><tr><td>σ<sup>c</sup><sub>31</sub></td><td><strong>0.33</strong></td><td>0.32</td><td><strong>0.55</strong></td><td>0.44</td></tr></tbody></table>

through the outside skin first, and then through the Nomex® honeycomb core (or through the foam, depending on the localization of the impact). Predicting the response of the core under out-of-plane loadings is then essential to make the calculation of the sandwich panel. The design of shock absorbers is another application that make intensive use of nonlinear models to enable the evaluation of the absorbed energy during the impact.

## 4.1. Compressive Tests on Sandwich Specimens

In this section, a campaign of compressive tests on sandwiches made of Nomex® honeycomb core and carbon fibers and epoxy resin composite skins is presented. The mechanical behavior of two kinds of Nomex® honeycomb cores from Hexcel is investigated. The first one, an over-expanded rectangular Nomex® honeycomb core has the following characteristics: 50 mm length, 50 mm width and 40 mm thick; cell size: 6.4 mm; density: 64 kg.m⁻³. The second one, an hexagonal Nomex® honeycomb core is: 50 mm length, 40 mm width and 25 mm thick; for a cell size of 3.175 mm and a density of 80 kg.m⁻³.

Two TR-50 carbon fibers and epoxy resin composite skins of thickness 0.68 mm and stacked at [45–45] are stacked on each side of the core to complete the specimen. The tests are carried out on the same tensile testing machine as presented before. The specimens are placed between two plates of compression and crushed to a maximum strain of: 0.085 m.m⁻¹ with an imposed velocity of 0.04 mm.s⁻¹ for the first type of honeycombs and 0.08 m.m⁻¹ with an imposed velocity of 0.025 mm.s⁻¹ for the second type of honeycombs. This corresponds to a compressive quasi-static strain rate of 0.001 s⁻¹ in both cases. Moreover, cyclic loadings and unloadings are imposed to the first type of specimens.

The results are summarized below for compressive tests on hexagonal honeycombs. Figure 12 presents the response of the hexagonal Nomex® honeycomb cores. The fully crushed specimen is also visible on the graph. For compressive tests on over-expanded rectangular honeycombs, results are available in Figure 9.

Four major crushing steps appear from the stress versus strain curves:

* First step: the linear growth of the stress versus strain— it corresponds to the reversible elastic response of the Nomex® honeycomb cores and occurs at low level of strain (less than 0.05 m.m⁻¹);

* Second step: the brutal decrease of stress—this phenomenon is linked to the apparition of a huge macroscopic buckling ply orthogonal to the direction of the load as clearly visible on Figure 9, here, failure stress is around 3.6 MPa;

* Third step: the slightly rising plateau—it is due to the successive formation of folds on the walls;

<figure><img src="image_13.png" /><figcaption>FIG. 12. Hexagonal honeycomb core submitted to compressive test.</figcaption></figure>
* fourth step: the increasing stiffness of the specimen response—it traduces the fact that there is no longer formation of buckling ply, only densification occurs at this point.

The last observation coming from the stress versus strain curves is that the tests are very repetitive which give more confidence in the data. In what follows a model dedicated to the crushing behavior of Nomex® honeycomb cores is presented. This model tries to represent the response of the strucutral material that is Nomex® honeycomb core in the sense of an homogenized macroscopic behavior.

## 4.2. Modeling of the Compressive Behavior

The unidirectional model presented here relies on the thermodynamics of irreversible process [13]. The Helmholtz's volumic free energy potential is chosen of the following form:

$$ \begin{cases} 2\psi(\epsilon_e; d; p) = \epsilon_e E (1-d) \epsilon_e & \text{if } d < d_c \\ 2\psi(\epsilon_e; d; p) = \epsilon_e E \epsilon_e + 2h(p) & \text{elsewhere} \end{cases} \quad (12) $$

where the following variables are respectively:

* $\psi$: the Helmholtz's volumic free energy [J.m$^{-3}$];

* $\epsilon_e$: the recoverable part of the total strain, i.e., elastic strain [m.m$^{-1}$];

* $E$: the elasticity modulus of the undamaged Nomex® honeycomb core [Pa];

* $d$: a damage variable [m$^2$.m$^{-2}$] (damage area for initial undamage unit area);

* $d_c$: a critical damage that traduces the initiation of the plateau regime;

* $h(p)$: a function of the an accumulated plastic strain $p$ associated with a nonlinear isotropic hardening [J.m$^{-3}$].

The state laws can now be derived from the Helmholtz's volumic free energy leading to the expressions of $\sigma$, the Cauchy's stress, $Y$ the thermodynamic force associated to $d$ and $R$ the stress linked to nonlinear isotropic hardening. If $d < d_c$ then:

$$ \begin{cases} \sigma = \frac{\partial \psi}{\partial \epsilon_e} = E(1-d)\epsilon_e & [\text{Pa}] \\ Y = -\frac{\partial \psi}{\partial d} = \frac{1}{2}\epsilon_e E \epsilon_e & [\text{J.m}^{-3}] \end{cases} \quad (13) $$

Else, if $d \ge d_c$ then:

$$ \begin{cases} \sigma = \frac{\partial \psi}{\partial \epsilon_e} = E(\epsilon_t - \epsilon_p) & [\text{Pa}] \\ R = \frac{\partial \psi}{\partial p} = \left(\frac{\partial h(p)}{\partial p}\right) = (Q_1 + Q_2 \langle p - p^y \rangle_+^\beta) & [\text{Pa}] \end{cases} \quad (14) $$

with:

* $\epsilon_t$: the total strain [m.m$^{-1}$];

* $\epsilon_p$: the plastic strain which stands for the inelastic strains observed in experiments [m.m$^{-1}$];

* $p^y$: the cumulated plastic strains at the end of the plateau regime [m.m$^{-1}$];

* $Q_1, Q_2$ and $\beta$ material constants ([Pa] for $Q_1$ and $Q_2$, no unit for $\beta$).

Once the state laws are obtained, a dissipation potential has to be defined to determine the evolution laws. This potential is related to all the dissipative phenomena, which means, in the case of this study, to damage and to plasticity with nonlinear isotropic hardening:

$$ \phi(\sigma; Y; R; (d)) = \phi_p(\sigma; Y; R; (d)) + \phi_d(\sigma; Y; R; (d)) \quad [\text{J.m}^{-3}] \quad (15) $$

where:

* $\phi_p(\sigma; Y; R; (d)) = f(\sigma; Y; R; (d)) = |\frac{\sigma}{1-d}| - \sigma^y - R$: is the plastic dissipation potential;

* $\phi_d(\sigma; Y; R; (d))$: is the damage dissipation potential;

with:

* $\sigma^y$: the plastic yield stress [Pa].

Here it can be noted that the model is constructed under the assumption of associated plasticity since the dissipation potential is taken to equal the plastic flow surface $f$ which delimits the recorvable states from the irreversible ones. The coupling between plasticity and damage is made as proposed by P. Ladevèze and E. Le Dantec [14] with the introduction of effective stress $|\frac{\sigma}{1-d}|$ which is the stress that might be applied on an undamaged specimen to obtain the same state of strain as that obtained by imposing $\sigma$ on the damaged specimen. In what concerns the damage dissipation potential, it is supposed to exist, to be convex in the space of stress variables, to be positive and to take a zero value at its origin ($\phi_d(0; 0; 0; 0) = 0$).

By derivation, the evolution laws are then obtained:

$$ \begin{cases} \dot{d} = \frac{1}{\tau}(1 - e^{-\alpha(\frac{\sqrt{Y} - \sqrt{Y_c}}{\sqrt{Y_c}} - d)_+}) & \text{if } d < d_c, 0 \text{ elsewhere} \\ \dot{e}_p = \dot{\lambda}^p \frac{\partial f}{\partial \sigma} & \text{if } d \ge d_c, 0 \text{ elsewhere} \\ \dot{p} = -\dot{\lambda}^p \frac{\partial f}{\partial R} & \text{if } d \ge d_c, 0 \text{ elsewhere} \end{cases} \quad (16) $$

with:

* $\dot{\lambda}^p$: the plastic multiplier that must check the consistency condition $\dot{\lambda}^p \dot{f} = 0$:

* $f < 0$: elastic domain;

* $f = 0$ and $\dot{\lambda}^p > 0$: plastic loading;

* $f = 0$ and $\dot{\lambda}^p = 0$: plastic unloading;

* $Y_o$: the threshold that drives the initiation of damage [J.m$^{-3}$];
<figure><img src="image_14.png" /><figcaption>FIG. 13. Compressive response of Nomex® honeycomb cores.</figcaption></figure>

* $Y_c$: a critical energy release rate [J.m⁻³];

* $\alpha$ and $\tau$: material parameters (no unit for $\alpha$, [s.m².m⁻²] for $\tau$).

The interest of using a delayed damage model has been highlighted in many works such as those of P. Ladevèze [15] and A. Suffis [16]. It introduces an internal length that results in the apparition of an upper bound for the damage rate $\dot{d}$ which is equal to $1/\tau$. As a result, the localization phenomenon of damage into a single element (or spatial discretization dependence of the results in finite element analysis) is avoided.

Equation (17) gives the expression of the plastic multiplier:

$$ \dot{\lambda}^p = \frac{f}{\partial_{\sigma} f E \partial_{\sigma} f + \partial_{R} f \partial_{p} R \partial_{R} f} \quad (17) $$

The last point is the distinction between tensile and compressive behavior. It is done by a test on the sign of the stress at time $t_n$. If $\sigma(t_n)$ is positive, the tensile behavior (brittle elastic on stress failure criterion [2]) is integrated, if not, the compressive behavior presented above is used. Figure 13 illustrates the different steps.

The delayed damage model associated with linear elasticity is used until the damage variable $d$ has reached $d_c$, a critical damage level. Then, at this stage ($d \ge d_c$), the following variables are updated:

* the damage variable is locked: $d = d_c$;

* the yield stress is evaluated: $\sigma^y = E(1-dc)\epsilon_e$ (with $\epsilon_e = \epsilon_t$ at this point);

* inelastic strains are computed after elastic unloading: $\epsilon_p = d_c\epsilon_t$;

* and elastic strains are updated: $\epsilon_e = (1-d_c)\epsilon_t$.

The stress is then computed: $\sigma = E\epsilon_e = E(1-d_c)\epsilon_t$ with the mechanical behavior relation (it is equal to $\sigma^y$ at the moment when instability appears).

Next, a usual return mapping plasticity algorithm [17] is used to update variables. Figure 13 illustrates the two major steps of the model.

### 4.3. Identification and Results

To determine the values of the material parameters involved in the model, an optimization procedure coupled with classical empirical rules is developed. Hence, the elasticity modulus $E$ is measured as usual on the compressive stress versus compressive strain curve as the director coefficient in the elastic step of the response. Then the critical damage $d_c$ that initiates the plasticity regime, comes from the apparent elasticity modulus at the starting point of the regime plateau (see Figure 13). The yield cumulated plastic strain $p^y$ is defined as the inelastic strain that appears after elastic unloading just before the densification step.

For the other coefficients of the model:

* $\alpha$, $\tau$, $Y_o$ and $Y_c$: for the damage model;

* $Q_1$, $Q_2$ and $\beta$: for the nonlinear isotropic behavior;

an optimization procedure based on a direct search method: the pattern search algorithm, is employed [18]. The only restriction is the choice of the value of $\tau$ (which is still an open question) that is the upper bound of the damage rate. Its value is fixed at $\tau = 10^{-5}$ s⁻¹ as usually found in literature [16]. The pattern search algorithm is presented below:

1. initialize: $k = 0$; a direction of search: $\eta_k$ in the space of parameters; a range of search: $r_k$; the initial guess of parameters: $\Omega_k$;

2. starting from the initial set of parameters $\Omega_0 = [\alpha; Y_o; Y_c; Q_1; Q_2; \beta]_0$, compute its cost function $\Upsilon_0$ (defined later) to estimate the quality of this set;

3. next, make a translation in direction $\eta_k$ of range $r_k$ in the space of the parameters to find a new candidate $\Omega_{k+1}$ and evaluate its associated cost function $\Upsilon_{k+1}$;
4. make a test:

a) if $\Upsilon_{k+1} \le \Upsilon_k$ then, the direction of search is a promising direction, and so increase the range of search $r_{k+1} > r_k$ without changing the direction of search $\underline{\eta}_{k+1} = \underline{\eta}_k$; go to step 3;

b) if $\Upsilon_{k+1} > \Upsilon_k$ then, the direction of search is not a suitable direction, and so change $\underline{\eta}_{k+1} \ne \underline{\eta}_k$ (with an arbitrary empirical rule to test in all the directions of the parameters space); go to step 3 without changing the range of change $r_{k+1} = r_k$ until a better candidate is has been found; if such a candidate can't be found, decrease the range of search $r_{k+1} < r_k$ and go to step 3.

This algorithm, quite simple, has several interesting advantages:

* its ability to deal with numerous parameters;

* its robustness (if an incompatible set of parameters is met leading to an unreal cost function for example, the algorithm changes its direction of search);

* its capacity to check a wide area in the space of parameters.

However, it has some non negligible drawbacks such as a starting point dependency and a tendency to converge on a local minimum. In this study, the drawbacks are not a limiting point since even if the solution does not correspond to the global minimum, it has no real effect on the ability to perform calculations with the model and on the capacity to match experimental data. Many starting points were checked (initial guess $\Omega_0$), leading for an important part of them, to similar results in terms of set of parameters. To finish with the presentation of the method, the cost function is defined as an error in the sense of least squares

<table><caption>TABLE 3
Parameters values for the model.</caption><tr><td>E</td><td>128 MPa</td><td>Q<sub>1</sub></td><td>2.42 J.mm<sup>-3</sup></td><td>α</td><td>10</td></tr><tr><td>d<sub>c</sub></td><td>0.87 m<sup>2</sup>·m<sup>-2</sup></td><td>Q<sub>2</sub></td><td>450 J.mm<sup>-3</sup></td><td>τ</td><td>10<sup>-5</sup> s<sup>-1</sup></td></tr><tr><td></td><td></td><td>β</td><td>0.54</td><td>Y<sub>o</sub></td><td>0.06 J.mm<sup>-3</sup></td></tr><tr><td></td><td></td><td>p<sup>y</sup></td><td>0.53 m.m<sup>-1</sup></td><td>Y<sub>c</sub></td><td>0.01 J.mm<sup>-3</sup></td></tr></table>

between experimental data and computed values:

$$ \Upsilon (\Omega_k) = \sum_{i} (\text{data}_i - \text{computedvalues}_i (\Omega_k))^2 \quad (18) $$

at each point i of the curves. After the step of optimization all the parameters values are available. They are summarized in Table 3 for the over expanded rectangular honeycombs:

With these parameters, the crushing stress versus crushing strain curve computed is compared to experimental data in Figure 14:

From a general point of view, the computed values are very similar to the experimental data that were used to make the identification. The four steps of the failure of Nomex® honeycomb cores: elasticity, apparition of the instability, plateau and densification are reproduced with a satisfactory accuracy. However, at low level of strain, just after the apparition of the instability the unloading is not well represented, it is stiffer with the model than in experiment. This can be explained by the instantaneous change of mechanical behavior in the model, change of mechanical behavior which might be softer for real. Another small drawback of the model appears at very high level of strains. Experiments exhibit an unloading stiffer than the initial elasticity modulus. This is due to the honeycomb sheets which are under

<figure><img src="image_15.png" /><figcaption>FIG. 14. Comparison between computed and experimental crushing curves.</figcaption></figure>
compressive loading between the two skins of the sandwich. The model is not able to deal with this phenomenon.

As a conclusion for this first unidirectional model dedicated to the nonlinear behavior of Nomex® honeycomb cores, it can be said that the obtained results are quite promising. The next step is to understand the coupling between out-of-plane compression and shear in order to propose a full three dimensional nonlinear behavior for Nomex® honeycomb cores.

## CONCLUSION

In this work, the first point was to develop software dedicated to the determination of the effective elastic properties of any Nomex® honeycomb core (any geometry of cell, any cell size, any density, any thickness). This tool, relying on the strain based periodic homogenization technique, gave satisfactory results compared with manufacturer's data. The next point was the determination of the strengths of Nomex® honeycomb core. In agreement with experimental observations, the linear Euler's buckling analysis was introduced to evaluate the critical macroscopic strain field that leads to failure. Thanks to the homogenized mechanical behavior, the ultimate stress was computed, and, once again, good results have been obtained from data coming from manufacturers. The last point was a prospect on the nonlinear behavior of Nomex® honeycomb cores submitted to out-of-plane compressive loadings. A model based on the thermodynamics of irreversible effects, including damage and plasticity was developed, and an identification technique proposed. Very good adequation was obtained between experimental data and computed values.

## REFERENCES

1. Gornet, L., Marckmann, G., and Ollier, G., "Interactions modèles expériences sur des âmes nids d'abeilles Nomex®: Application au design

d'un voilier multicoque de course océanique," *Revue des composites et des matériaux avancés*, a-b, Editions Hermès (2006).

2. Tsaï, S. W., and Wu, E. M., "A general theory of strength for anisotropic materials," *Int. J. Composites and Materials*, **5**, 58-80 (1971).

3. Hoffman, O., "The brittle strength of orthotropic materials," *Int. J. Composites and Materials*, **1**, a-b (1967).

4. HexWebTB Honeycomb Attributes and properties, 21, Hexcel Composites (1999).

5. Euro-composite, Nomex® cores, http://www.euro-composites.com

6. Gellatley, R. A., and Clark, B. W., "The shear modulus of foil honeycomb cores," *Aircraft Engineering*, **30**, 294-302 (1958).

7. Gibson, L. J., and Ashby, M. F., "Cellular solids structures and properties," a-b, Pergamon Press, Oxford (1988).

8. Hohe, J., and Becker, W., "Effective elastic properties of triangular grid structures," *Composite Structures*, **45**, 131-145 (1999).

9. Aminanda, Y., Castanié, B., Barrau, J.-J., and Thevenet, P., "Experimental analysis and modeling of the crushing of honeycomb cores," *Applied Composite Materials*, **12**, 213-227 (2005).

10. Grediac, M., "A finite element study of the transverse shear in honeycomb cores," *Int. J. Solids Structures*, **30**, 1777-1788 (1993).

11. Suquet, P., "Plasticité et homogénéisation," *Ph.D thesis., Université Paris, 6* (1982).

12. Bathe, K.-J., "Finite element procedures," a-b, Prentice-Hall (1996).

13. Lemaître, J., and Chaboche, J. L., "Mechanics of solid materials," a-b, Cambridge University Press (1990).

14. Ladevèze, P., and Le Dantec, E., "Damage modelling of the elementary ply for laminated composites," *Composites Science and Technology*, **43**(3), 257-267 (1992).

15. Ladevèze, P., Allix, O., Deü, J.-F., and Lévêque, D., "A mesomodel for localization and damage computation in laminates," *Computer Methods in Applied Mechanics Engineering*, **183**, 105-122 (2000).

16. Suffis, A., Lubrecht, T. A. A., and Combescure, A., "Damage model with delay effect: Analytical and numerical studies of the evolution of the characteristic damage length," *Int. J. Solids and Structures*, **40**(13-14), 3463-3476 (2003).

17. Simo, J. C., and Hugues, T. J. R., "Computational Inelasticity," a-b, Springer, New York (1998).

18. Lewis, R. M., Torczon, V., and Trosset, M. W., "Direct search methods: then and now," *Int. J. Computational and Applied Mathematics*, **124**, 191-207 (2000).