# Accounting for differential swelling in the multi-physics modelling of the diffusive behaviour of polymers

Bun Eang Sar, Sylvain Fréour, Peter Davies, Frédéric Jacquemin

► To cite this version:

Bun Eang Sar, Sylvain Fréour, Peter Davies, Frédéric Jacquemin. Accounting for differential swelling in the multi-physics modelling of the diffusive behaviour of polymers. Journal of Applied Mathematics and Mechanics / Zeitschrift für Angewandte Mathematik und Mechanik, 2014, 94 (6), pp.452-460. 10.1002/zamm.201200272 . hal-01004970

HAL Id: hal-01004970

https://hal.science/hal-01004970

Submitted on 11 Feb 2017

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Accounting for differential swelling in the multi-physics modelling
of the diffusive behaviour of polymers

B. E. Sar<fnref n="1" />, S. Fréour<fnref n="1" />, P. Davies<fnref n="2" />, and F. Jacquemin<fnref n="1" />,*

**Key words** Hygromechanical equation, moisture absorption, thermodynamical approach, unsymmetrical loading of moisture.

The purpose of this paper is to take into account the derivative by moisture content of polymer volume in order to establish a diffusion law within the so-called “thermodynamical approach” for a polymer material which experiences a hygro-mechanical load. In this study, the specific case corresponding to the existence of unsymmetrical hygroscopic boundary conditions was investigated.

# 1 Introduction

Polymers and polymeric composites absorb moisture when exposed to ambient humidity or immersed in liquid. Polymeric matrix composites differ from other materials in the sense that low-molecular weight substances such as water may easily migrate even at room temperature, generating a variation of the material's structure, morphology, and composition. Moreover, many authors have reported that hygro-thermal ageing could induce a loss of the mechanical stiffness and/or strength of organic matrix composites [15, 19]. It is probable that the factors described above will also affect the moisture sorption behavior of polymer matrix composites. In order to predict the time-dependent evolution of the moisture content of composite structures, various models have been developed in the literature. Among them, some are based on the classical Fickian diffusion model [8, 9, 11, 20]. More recently, Fick's model has successfully been combined with scale transition models such as the Eshelby-Kröner self-consistent model for predicting multi-scale distribution of the internal mechanical states during the transient step of the moisture diffusion process experienced by polymer composites [7, 10].

Nevertheless, some experimental data demonstrate that the moisture sorption in composite structures could differ from the typical Fickian uptake [3, 16]. As a consequence, some researchers have developed models in order to reproduce the anomalous sorption curves observed in practice [5, 21]. Among the proposed methods, [17] documented a multi-physics approach to the diffusion mechanisms, compatible with the thermodynamics. The approach is similar to that presented by Larché and Cahn or Aifantis and Gerberich for predicting the diffusion of gases through elastic solids [1, 2, 12]. The multiphysics thermodynamic model proposed by Larché and Cahn was later implemented by Neogi et al. who achieved the successful fitting of experimental results obtained on thin polymer membranes [14]. Nevertheless, in these pioneering works, the differential swelling was treated owing to simplifying assumptions relating the deformation field to the existing penetrant concentration [13].

More recently, other mutliphysics model coupling the mechanical states to mass-transport process were developed in the case that linear viscoelastic solids were considered [4]. An important feature of that formulation, although limited to the one-dimensional case, is that the expressions used for the chemical potential and the stress constitutive equations are thermodynamically consistent, since they come from the equation describing the Helmholtz free energy [4].

In recent works [6,18], other models, focused on the description of anomalous diffusion, were also developed which were compatible with the thermodynamics. Nevertheless, the mathematical formalism presented in both references [6, 18] does not enable the effects on the moisture kinetics induced by the presence of an in-depth heterogeneous profile of the hygro-elastic strain to be accounted for. The present work will present a possible way to address this issue. The developments detailed in this paper will also extend the formalism, so that an unsymmetrical hygroscopic load can be considered, whereas

<footnote n="1">¹ Institut de Recherche en Génie Civil et Mécanique (UMR CNRS 6183), LUNAM Université, Université de Nantes,
Ecole Centrale de Nantes, 37, Bd de l'Université, BP 406, 44602 Saint-Nazaire, France</footnote>

<footnote n="2">² Service Matériaux et Structures, IFREMER - Centre de Brest, B.P. 70, 29280 Plouzané, France</footnote>
only symmetrical cases could be modeled using the historical version of the model [6, 18], as well as according to the original pioneering papers published in this very field of research [1, 2, 12].

## 2 Hygroscopic pressure

Moisture absorption induces swelling strains that actually correspond to the existence of a hygroscopic pressure within the material. The in-depth, time-dependent hygroscopic pressure profile occurring during the transient stage of the diffusion process is determined according to the three following equations: (1) Hygro-elastic Hooke's law, (2) Equilibrium equations, and (3) Compatibility equations.

$$ \varepsilon_{il} = \frac{1 + \nu}{E} \sigma_{il} - \frac{\nu}{E} \delta_{il} \text{tr} \sigma_{il} + \eta C \delta_{il}, \quad (1) $$

$$ \sigma_{il,l} = 0, \quad (2) $$

$$ \varepsilon_{il,jk} + \varepsilon_{jk,il} - \varepsilon_{jl,ik} - \varepsilon_{ik,jl} = 0, \quad (3) $$

where $\nu$ is Poisson's ratio, $E$ the Young's modulus and $\eta$ the coefficient of moisture expansion of the polymer (CME). $C$ denotes the moisture content (assuming the material to be initially dry) while $\delta_{il}$ stands for the Kronecker's symbol, i.e.

$\delta_{il} = \begin{cases} 1 & (i = l) \\ 0 & (i \neq l) \end{cases}$. For a given set of indices $(i, l)$ in (1)–(2), we use the following replacement rule $j = k = 1, 2, 3$ in Eq. (3), the summation of which yields:

$$ \Delta\varepsilon_{il} + \varepsilon_{kk,il} - (\varepsilon_{ik,lk} + \varepsilon_{lk,ik}) = 0. \quad (4) $$

Accounting for the hygro-elastic Hooke's law (1), the sum $(\varepsilon_{ik,lk} + \varepsilon_{lk,ik})$ appearing above in relation (4) actually satisfies the following equation

$$ \varepsilon_{ik,lk} + \varepsilon_{lk,ik} = \frac{1 + \nu}{E} (\sigma_{ik,lk} + \sigma_{lk,ik}) - \frac{\nu}{E} (\sigma_{kk,lk} \delta_{ik} + \sigma_{kk,ik} \delta_{lk}) + \eta (C_{,lk} \delta_{ik} + C_{,ik} \delta_{lk}). \quad (5) $$

Since $\sigma_{ik,lk} = \sigma_{lk,ik} = 0$; $\sigma_{kk,lk} \delta_{ik} = \sigma_{kk,ik} \delta_{lk} = \sigma_{kk,il}$; $C_{,lk} \delta_{ik} = C_{,ik} \delta_{lk} = C_{,il}$, many terms cancel in Eq. (5) that can be written in the following simplified form:

$$ \varepsilon_{ik,lk} + \varepsilon_{lk,ik} = -\frac{2\nu}{E}\sigma_{kk,il} + 2\eta C_{,il}. \quad (6) $$

Substituting Eq. (6) into Eq. (4) yields

$$ \Delta\varepsilon_{il} + \varepsilon_{kk,il} + \frac{2\nu}{E}\sigma_{kk,il} - 2\eta C_{,il} = 0. \quad (7) $$

Considering the replacement rule $i = 1$ in (7) yields

$$ \Delta\varepsilon_{ll} + \varepsilon_{kk,ll} + \frac{2\nu}{E}\sigma_{kk,ll} - 2\eta C_{,ll} = 0. \quad (8) $$

Actually, the Laplacian of moisture content is written as $C_{,ll} = \Delta C$. Moreover, $\Delta\varepsilon_{ll} = \Delta\varepsilon_{kk} = \varepsilon_{kk,ll}$. As a result, the relation can be simplified as follows

$$ \Delta\varepsilon_{kk} + \frac{\nu}{E}\sigma_{kk,ll} - \eta\Delta C = 0. \quad (9) $$

From Eq. (1), the second derivative of the hygro-elastic strain trace, $\varepsilon_{kk,il}$, featured in Eq. (7), satisfies:

$$ \varepsilon_{kk,il} = \frac{1 - 2\nu}{E}\sigma_{kk,il} + 3\eta C_{,il}. \quad (10) $$

Putting $i = 1$ into Eq. (8) provides the following expression for the Laplacian of the trace of the hygro-elastic strain, $\Delta\varepsilon_{kk}$

$$ \varepsilon_{kk,ll} = \Delta\varepsilon_{kk} = \frac{1 - 2\nu}{E}\Delta\sigma_{kk} + 3\eta\Delta C. \quad (11) $$
Combining (9) to (11) yields

$$ \frac{1-2\nu}{E}\Delta\sigma_{kk} + 3\eta\Delta C + \frac{\nu}{E}\sigma_{kk,ll} - \eta\Delta C = 0 \Leftrightarrow \frac{1-2\nu+\nu}{E}\Delta\sigma_{kk} + 2\eta\Delta C = 0, \quad (12) $$

$$ \Delta\sigma_{kk} = -2\eta \frac{E}{1-\nu} \Delta C. \quad (13) $$

In the present work, the trace of stress tensor $\sigma_{kk}$ is considered to correspond to the following sum of an external mechanical load, $P_{ex}$ and a hygroscopic pressure, $P_{is}$, so that $\sigma_{kk} = -3(P_{ex} + P_{is})$ where $P_{ex}$ is a constant parameter. Thus, $\Delta P_{ex} = 0$ and Eq. (13) can be reduced to

$$ \Delta P_{is} = \frac{2E}{3(1-\nu)} \eta \Delta C = \frac{\alpha}{A_0} \eta \Delta C, \quad (14) $$

where the constants $\alpha$ and $A_0$ have already been defined in previous works [18], as

$$ \frac{\alpha}{A_0} = \frac{2E}{3(1-\nu)}, \quad (15) $$

$$ A_0 = \frac{3\omega_w}{RT\rho_0}, \quad (16) $$

where $\rho_0$ is the density of polymer resin at free strain state, whereas $\omega_w$ stands for the molar mass of water, $T$ the temperature and $R$ the ideal gas constant.

We consider a plate whose lateral dimensions are large compared to the thickness. As a consequence, the diffusion is considered to occur along the direction $x$, only. The unidirectional solution of Eq. (14) satisfies the following general form

$$ P_{is}(x,t) = \frac{\alpha}{A_0} \eta C(x,t) + k_1(t)x + k_2(t). \quad (17) $$

The constants $k_1(t)$ and $k_2(t)$ are deduced from the equilibrium conditions, in which $L$ stands for the thickness of the sample

$$ \left\{ \begin{array}{l} \displaystyle \int_0^L P_{is}(x,t)dx = 0, \\ \\ \displaystyle \int_0^L P_{is}(x,t)xdx = 0. \end{array} \right. \quad (18) $$

The solutions satisfying the system of Eqs. (18) are

$$ k_1(t) = \frac{6}{L^3} \frac{\alpha}{A_0} \eta \left( L^2 \overline{C(t)} - 2I \right), \quad (19) $$

$$ k_2(t) = \frac{2}{L^2} \frac{\alpha}{A_0} \eta \left( 3I - 2L^2 \overline{C(t)} \right), \quad (20) $$

where

$$ \overline{C(t)} = \frac{1}{L} \int_{0}^{L} C(x,t) dx, \quad (21) $$

$$ I = \int_{0}^{L} C(x,t) x dx. \quad (22) $$

Introducing (19)-(20) in the general expression (17) for the internal pressure yields

$$ P_{is}(x,t) = \frac{\alpha}{A_0} \eta \left( C(x,t) - 4\overline{C(t)} \right) + \frac{6}{L^3} \frac{\alpha}{A_0} \eta x \left( L^2 \overline{C(t)} - 2I \right) + \frac{6}{L^2} \frac{\alpha}{A_0} \eta I. \quad (23) $$

In the case that a symmetrical hygroscopic load is applied on the boundaries of the structure, $C(x,t) = C(L-x,t)$. As a result, the integration $I$ is equal to $\frac{L^2 \overline{C(t)}}{2}$. Hence, the corresponding hygroscopic pressure is given by the simplified form $P_{is}(x,t) = \frac{\alpha}{A_0} \eta \left( C(x,t) - \overline{C(t)} \right)$.
### 3 Chemical potential

The chemical potential of water $\tilde{\mu}_w$ is defined as the partial derivative of free energy of Helmholtz, $F = F_0 + n f_w(C) + V_\epsilon W$, with respect to the amount of water $n_w$. Where $F_0$ is the free energy of the dry stress-free polymer, $f_w(C)$ is the variation of the free energy per mole of dry polymer, due to the addition of water when the polymer is free to swell, $n$ and $V_\epsilon$ are respectively the amount of polymer and its volume at any stage, whereas $W$ denotes the hydro-elastic strain energy [6]

$$ \tilde{\mu}_w(C) = \frac{\partial F}{\partial n_w} = \frac{\partial F}{\partial C} \frac{\partial C}{\partial n_w}. \qquad (24) $$

The moisture content in the polymer is calculated through $\frac{\partial C}{\partial n_w} = \frac{\omega_w}{n\omega} \frac{\partial c}{\partial n_e} = \frac{\partial}{\partial n_e} \left( \frac{n_e \omega_e}{n_p \omega_p} \right) = \frac{1}{n_p} \frac{\omega_e}{\omega_p}$, here $\omega_w, \omega$ stand for respectively the molar mass of water and polymer.

The hydro-elastic strain energy written as a function of both the bulk modulus $k$ and shear modulus $G$, is defined by

$$ W = \frac{1}{2} \sigma : \varepsilon^{el} = \frac{k}{2} (\operatorname{tr}\varepsilon - 3\eta C)^2 + Ge : e, \qquad (25) $$

where $\varepsilon^{el}$ is the elastic strain, $\varepsilon$ being the total strain, whereas $e$ is the deviatoric strain tensor.

Introducing $f'_w(C) = \frac{\partial[f_w(C)]}{\partial C}$, one obtains the following expression for the derivative of the Helmholtz free energy with respect to the moisture content

$$ \frac{\partial F}{\partial C} = \frac{\partial F_0}{\partial C} + n f'_w(C) + W \frac{\partial V_\epsilon}{\partial C} + V_\epsilon \frac{\partial W}{\partial C}. \qquad (26) $$

During the moisture diffusion process, we take into account the evolution of the volume occupied by the polymer, and the resulting variation of its density, through:

$$ \frac{V_{\epsilon}}{V_0} = \frac{\rho_0}{\rho_{\epsilon}} = \operatorname{tr}\epsilon + 1, \qquad (27) $$

where $V_\epsilon, V_0, \rho_0, \rho_\epsilon$ stand respectively for the polymer volume and its density at present (strained) and initial (strain-free) states

$$ \frac{\partial V_{\epsilon}}{\partial C} = \frac{\partial V_{\epsilon}}{\partial \operatorname{tr}\epsilon} \frac{\partial \operatorname{tr}\epsilon}{\partial C} = V_0 \frac{\partial \operatorname{tr}\epsilon}{\partial C} = \frac{n\omega}{\rho_0} \frac{\partial \operatorname{tr}\epsilon}{\partial C}. \qquad (28) $$

Let us consider (27) as well as the equation $e : e = 0$ (which comes from the specific case, considered here, of a macroscopically isotropic polymer submitted to a hydrostatic pressure) in the expression of the hydro-elastic strain energy (25). As a result, the partial derivative of Helmholtz free energy $F$ with respect to the moisture content (26) transforms as follows

$$ \frac{\partial F}{\partial C} = n f'_w(C) + V_\epsilon k (\operatorname{tr}\epsilon - 3\eta C) \frac{\partial (\operatorname{tr}\epsilon - 3\eta C)}{\partial C} + \frac{k}{2} (\operatorname{tr}\epsilon - 3\eta C)^2 \frac{n\omega}{\rho_0} \frac{\partial \operatorname{tr}\epsilon}{\partial C}. \qquad (29) $$

Accounting for Eq. (29), the chemical potential (24) eventually satisfies

$$ \tilde{\mu}_w(C, \operatorname{tr}\epsilon) = \frac{\omega_w}{\omega} f'_w(C) + k \frac{\omega_w}{\rho_0} (\operatorname{tr}\epsilon - 3\eta C) \left( \frac{\partial \operatorname{tr}\epsilon}{\partial C} - 3\eta \right) (\operatorname{tr}\epsilon + 1) + \frac{\omega_w}{\rho_0} \frac{k}{2} (\operatorname{tr}\epsilon - 3\eta C)^2 \frac{\partial \operatorname{tr}\epsilon}{\partial C}. \qquad (30) $$

Besides, the trace of the strain tensor can be expressed as a function of the total pressure as follows

$$ \operatorname{tr}\epsilon = \frac{\sigma_{kk}}{3k} + 3\eta C = -\frac{P}{k} + 3\eta C = -\frac{P_{\text{ex}} + P_{\text{is}}}{k} + 3\eta C. \qquad (31) $$

Taking into account the expression (21) obtained for the internal pressure induced notably by the differential swelling, the derivative of relation (31) with respect to the moisture content satisfies

$$ \frac{\partial \operatorname{tr}\epsilon}{\partial C} = \eta \frac{3A_0 k - \alpha}{A_0 k}. \qquad (32) $$

According to Eq. (32) and considering that $\frac{\omega_w}{\omega} f'_w(C) = \mu_0 + RT \ln \frac{C}{C_0}$, the relation previously obtained for the expression of the chemical potential (30) can be developed as follows

$$ \tilde{\mu}_w(C, \operatorname{tr}\epsilon) = \mu_0 + RT \ln \frac{C}{C_0} - \frac{3\eta\omega_w k}{\rho_0} (\operatorname{tr}\epsilon - 3\eta C) (\operatorname{tr}\epsilon + 1) \\ + \frac{\eta\omega_w}{A_0\rho_0} (3A_0k - \alpha) (\operatorname{tr}\epsilon - 3\eta C) (\operatorname{tr}\epsilon + 1) + \frac{\eta\omega_w}{A_0\rho_0} \frac{3A_0k - \alpha}{2} (\operatorname{tr}\epsilon - 3\eta C)^2. \qquad (33) $$
## 4 Equation of model

Generally, the diffusion equation is deduced from the conservation of mass equation [13] in which the first derivative of moisture with respect to time, $\dot{C}$, relates to the diffusion flux of moisture, $J_i$, as follows

$$ \dot{C} + J_{i,i} = 0. \qquad (34) $$

In the present work, the diffusion flux of moisture was written in term of chemical potential $\tilde{\mu}_w$ as proposed, for instance, in [14]

$$ J_i = -\frac{DC}{RT} \text{grad} \tilde{\mu}_w, \qquad (35) $$

where D is the diffusion coefficient in [mm²/s], R is the gas constant in [J/(mol.K)], and T the absolute temperature [K].

We obtain the constitutive equation by using the mass conservation equation (34) in which the chemical potential of water has been written as a function of both the trace of the strains and the moisture content

$$ \dot{C} = \frac{D}{RT} \text{div} \left[ C \left( \text{grad} \tilde{\mu}_w (C, \text{tr}\varepsilon) \right) \right]. \qquad (36) $$

In vector calculus, the gradient operator satisfies the following linear property

$$ \text{grad} \tilde{\mu}_w (C, \text{tr}\varepsilon) = \frac{\partial \tilde{\mu}_w}{\partial C} \text{grad} C + \frac{\partial \tilde{\mu}_w}{\partial \text{tr}\varepsilon} \text{grad} \text{tr}\varepsilon. \qquad (37) $$

Introducing the development (37) in (36) yields

$$ \dot{C} = \frac{D}{RT} \text{div} \left[ C \left( \frac{\partial \tilde{\mu}_w}{\partial C} \text{grad} C + \frac{\partial \tilde{\mu}_w}{\partial \text{tr}\varepsilon} \text{grad} \text{tr}\varepsilon \right) \right]. \qquad (38) $$

The partial derivatives of the chemical potential by either the moisture content C or the strain trace trε can respectively be written as

$$ \frac{\partial \tilde{\mu}_w (C, \text{tr}\varepsilon)}{\partial C} = \frac{RT}{C} + \frac{9\eta^2 k \omega_w k}{\rho_0} (\text{tr}\varepsilon + 1) - \frac{3\eta^2 \omega_w}{A_0 \rho_0} (3A_0 k - \alpha) (2\text{tr}\varepsilon - 3\eta C + 1), \qquad (39) $$

$$ \frac{\partial \tilde{\mu}_w (C, \text{tr}\varepsilon)}{\partial \text{tr}\varepsilon} = -\frac{3\eta\omega_w k}{\rho_0} (2\text{tr}\varepsilon - 3\eta C + 1) + \frac{\eta\omega_w}{A_0\rho_0} (3A_0 k - \alpha) (3\text{tr}\varepsilon - 6\eta C + 1). \qquad (40) $$

Combining Eqs. (39) and (40) with the diffusion equation (38) leads to

$$ \begin{aligned} \dot{C} = \frac{D}{RT} \text{div} \left[ C \left\{ \left[ \frac{RT}{C} + \frac{9\eta^2 k \omega_w}{\rho_0} (\text{tr}\varepsilon + 1) \right] \text{grad} C - \frac{3\eta\omega_w k}{\rho_0} (2\text{tr}\varepsilon - 3\eta C + 1) \text{grad} \text{tr}\varepsilon \right\} \right] \\ + \frac{D}{RT} \text{div} \left[ C \left\{ -\frac{3\eta^2 \omega_w}{A_0 \rho_0} (3A_0 k - \alpha) (2\text{tr}\varepsilon - 3\eta C + 1) \text{grad} C \right\} \right] \\ + \frac{D}{RT} \text{div} \left[ C \left\{ \frac{\eta\omega_w}{A_0 \rho_0} (3A_0 k - \alpha) (3\text{tr}\varepsilon - 6\eta C + 1) \text{grad} \text{tr}\varepsilon \right\} \right]. \end{aligned} \qquad (41) $$

Further simplifications applied to the previous form (41) enable us to write

$$ \dot{C} = \frac{D}{RT} \text{div} \left[ C \left\{ \left[ \frac{RT}{C} + \frac{9\eta^2 k \omega_w}{\rho_0} (\text{tr}\varepsilon + 1) \right] \text{grad} C - \frac{3\eta\omega_w k}{\rho_0} (2\text{tr}\varepsilon - 3\eta C + 1) \text{grad} \text{tr}\varepsilon \right\} \right] + g, \qquad (42) $$

where

$$ \begin{aligned} g = \frac{D}{RT} \text{div} \left[ C \left\{ -\frac{3\eta^2 \omega_w}{A_0 \rho_0} (3A_0 k - \alpha) (2\text{tr}\varepsilon - 3\eta C + 1) \text{grad} C \right\} \right] \\ + \frac{D}{RT} \text{div} \left[ C \left\{ \frac{\eta\omega_w}{A_0\rho_0} (3A_0 k - \alpha) (3\text{tr}\varepsilon - 6\eta C + 1) \text{grad} \text{tr}\varepsilon \right\} \right]. \end{aligned} \qquad (43) $$
Equation (43) yields the following developed expression

$$g = D \frac{\eta (3A_0 k - \alpha)}{3} \left\{ \begin{aligned} & -3\eta (2 \operatorname{tr}\varepsilon - 3\eta C + 1) C \Delta C - 3\eta (2 \operatorname{tr}\varepsilon - 6\eta C + 1) \overrightarrow{\operatorname{grad}} C \overrightarrow{\operatorname{grad}} C \\ & + (3 \operatorname{tr}\varepsilon - 6\eta C + 1) C \Delta (\operatorname{tr}\varepsilon) + (3 \operatorname{tr}\varepsilon - 18\eta C + 1) \overrightarrow{\operatorname{grad}} C \overrightarrow{\operatorname{grad}} \operatorname{tr}\varepsilon + 3C \overrightarrow{\operatorname{grad}} \operatorname{tr}\varepsilon \overrightarrow{\operatorname{grad}} \operatorname{tr}\varepsilon \end{aligned} \right\}. \quad (44)$$

Finally, the factor g could be written as follows

$$g = D\xi \left[ z_1 \frac{\partial^2 C}{\partial x^2} + z_2 \left(\frac{\partial C}{\partial x}\right)^2 + z_3 \frac{\partial C}{\partial x} + 3C \left\{ \frac{6}{L^3} \frac{\alpha}{A_0 k} \eta \left( L^2 \overline{C(t)} - 2I \right)^2 \right\} \right], \quad (45)$$

where

$$z_1 = -3\eta (2 \operatorname{tr}\varepsilon - 3\eta C + 1) C + (3 \operatorname{tr}\varepsilon - 6\eta C + 1) C \frac{3A_0 k - \alpha}{A_0 k} \eta, \quad (46)$$

$$z_2 = -3\eta (2 \operatorname{tr}\varepsilon - 6\eta C + 1) + (3 \operatorname{tr}\varepsilon - 18\eta C + 1) \frac{3A_0 k - \alpha}{A_0 k} \eta + 3C \left( \frac{3A_0 k - \alpha}{A_0 k} \eta \right)^2, \quad (47)$$

$$z_3 = -\frac{6}{L^3} \frac{\alpha}{A_0 k} \eta \left( L^2 \overline{C(t)} - 2I \right) \left( 3 \operatorname{tr}\varepsilon - 18\eta C + 1 \right) - \frac{36C}{L^3} \frac{\alpha}{A_0 k} \left( L^2 \overline{C(t)} - 2I \right) \frac{3A_0 k - \alpha}{A_0 k} \eta^2, \quad (48)$$

$$\xi = \frac{3A_0 k - \alpha}{3} \eta. \quad (49)$$

Using the same method, the first term of the right hand side of Eq. (42) was developed, and then simplified. The resulting time-dependent diffusive behavior for a polymer plate subjected to an unsymmetrical humid ambient load is given by

$$\dot{C} = D \left[ \begin{aligned} & (1 + V_1 \eta^2 C + V_2 \eta^3 C^2) \frac{\partial^2 C}{\partial x^2} + \eta^2 (V_3 + V_4 C) \left(\frac{\partial C}{\partial x}\right)^2 \\ & - \frac{6}{L^3} \frac{\alpha}{A_0 k} \eta \left( L^2 \overline{C(t)} - 2I \right) (V_5 + V_6 \eta^2 C) \frac{\partial C}{\partial x} - \frac{72}{L^6} \frac{\alpha^2}{A_0 k} \eta^3 \left( L^2 \overline{C(t)} - 2I \right)^2 C \end{aligned} \right] + g, \quad (50)$$

where

$$V_1 = -3A_0 k \operatorname{tr}\varepsilon + 2\alpha \operatorname{tr}\varepsilon + \alpha, \quad V_2 = 9A_0 k - 3\alpha, \quad V_3 = -3A_0 k \operatorname{tr}\varepsilon + 2\alpha \operatorname{tr}\varepsilon + \alpha,$$

$$V_4 = \eta V_2 - \frac{2\eta\alpha^2}{A_0 k}, \quad V_5 = 2\eta A_0 \operatorname{tr}\varepsilon + A_0 \eta, \quad V_6 = 3A_0 - 4\frac{\alpha}{k}.$$

Significant simplifications of Eq. (50) can be made when the polymer structure is subjected to symmetrical moisture conditions. This requires that the equation $L^2 \overline{C(t)} - 2I = 0$ be satisfied.

The resulting behavior law then respects the following form

$$\dot{C} = D \left[ (1 + V_1 \eta^2 C + V_2 \eta^3 C^2 + \xi z_1) \frac{\partial^2 C}{\partial x^2} + \eta^2 (V_3 + V_4 C + \xi z_2) \left(\frac{\partial C}{\partial x}\right)^2 \right]. \quad (51)$$

## 5 Boundary conditions

The boundary condition is obtained by equating the chemical potential of water in humid air,

$$\hat{\mu}_w = \hat{\mu}_0 + RT \ln \frac{p_w}{p_0} \quad (\text{where})$$

$\hat{\mu}_0$ is the chemical potential of water in humid air at the reference pressure

$p_0$, the partial pressure of water being $p_w$), with the generalized chemical potential of the polymer, $\tilde{\mu}_w$, the expression for which is given by Eq. (33) in the present work. This statement yields the following equation
$$ \tilde{\mu}_w(C, \text{tr}\varepsilon) = \mu_0 + RT \ln \frac{C}{C_0} - \frac{3\eta\omega_w k}{\rho_0} (\text{tr}\varepsilon - 3\eta C) (\text{tr}\varepsilon + 1) \\ + \frac{\eta\omega_w}{A_0\rho_0} (3A_0k - \alpha) (\text{tr}\varepsilon - 3\eta C) (\text{tr}\varepsilon + 1) + \frac{\eta\omega_w}{A_0\rho_0} \frac{3A_0k - \alpha}{2} (\text{tr}\varepsilon - 3\eta C)^2. \quad (52) $$

The boundary condition is obviously only satisfied at the specific positions $x_b$ denoting the boundaries between the ambient fluid and the polymer. The equalization between the chemical potential of water in humid air and the generalized chemical potential of water in the system leads to the following moisture conditions at the structure boundaries

$$ C(x_b, t) = \frac{p_w}{p_0} C_0 \exp \left[ \frac{\hat{\mu}_0 - \mu_0}{RT} + \eta k A_0 (\text{tr}\varepsilon - 3\eta C) (\text{tr}\varepsilon + 1) \right. \\ \left. - \frac{\eta}{3} (3A_0k - \alpha) (\text{tr}\varepsilon - 3\eta C) (\text{tr}\varepsilon + 1) - \frac{\eta}{6} (3A_0k - \alpha) (\text{tr}\varepsilon - 3\eta C)^2 \right]. \quad (53) $$

Equation (53) could also be written as a function of the total pressure $P$ instead of $\text{tr}\varepsilon$ owing to their relation as expressed by (31). One can then write:

$$ C(x_b, t) = \frac{C_0}{p_0} p_w \exp \left( \frac{\hat{\mu}_0 - \mu_0}{RT} \right) \exp \left[ \eta k A_0 \left( -\frac{P}{k} \right) \left( -\frac{P}{k} + 3\eta C + 1 \right) \right. \\ \left. - \frac{\eta}{3} (3A_0k - \alpha) \left( -\frac{P}{k} \right) \left( -\frac{P}{k} + 3\eta C + 1 \right) - \frac{\eta}{6} (3A_0k - \alpha) \left( -\frac{P}{k} \right)^2 \right]. \quad (54) $$

Introducing Henry's law, $S = \frac{C_0}{p_0} \exp(\frac{\hat{\mu}_0 - \mu_0}{RT})$ into Eq. (54), the boundary condition for the moisture content becomes

$$ C(x_{b,t}) = S p_w \exp \left[ \left( \frac{\eta A_0}{k} - \frac{\eta}{3k^2} (3A_0k - \alpha) \right) (P^2 - 3\eta C k P - k P) - \frac{\eta}{6k^2} (3A_0k - \alpha) P^2 \right]. \quad (55) $$

# 6 Numerical results

The numerical simulations correspond to a 4 mm thick plate made of epoxy resin whose Young modulus is 3.65 GPa and Poisson's ratio is 0.36. The polymer plate is subjected to moisture diffusion while experiencing a hydrostatic pressure load.

We would like to simulate the moisture absorption within the above described material in the cases when a whether symmetrical or an unsymmetrical moisture condition takes place at the opposite edges of the plate.

## 6.1 Symmetrical moisture load

The opposite surfaces of the plate are assumed to be submitted to the same relative humidity, which correspond to a reference moisture content level $C_0$ of 5% (in the case that the multiphysics effects are neglected). The mathematical equation governing the diffusion corresponds to Eq. (51), whereas the boundary condition is obtained owing to expression (55).

Figure 1 shows the time-dependent evolution of the macroscopic (average) moisture content, as a function of the CME: $\eta = 0$; $\eta = 0.6$ or $\eta = 1$, at an imposed pressure of 1 MPa. Increasing CME reduces the maximum moisture absorption capacity of the polymer as indicated by the evolution of the average moisture content in the steady state.

According to Fig. 2, the moisture uptake in the polymer resin decreases with the reduction of the relative humidity on the second side of the plate. Non-linearities, similar to those observed on Fig. 1, occur in the presently considered cases, also. The previously so-called “delay time” before the establishment of a Fickian-like diffusion process can be observed, as well. The dependence of the apparent diffusion coefficient with time increases with the deviation of the environmental conditions applied to the opposite boundaries of the plate. Besides, the maximum moisture content attained in permanent regime clearly does not vary linearly with the boundary condition applied to the second surface of the plate.

According to Fig. 1, the multiphysics model predicts a fickian diffusion process in the case that the coefficient of moisture absorption of the polymer is assumed to be null. Discrepancies from the typical Fickian kinetics occur is predicted by the model when $\eta \neq 0$. In particular, the apparent moisture diffusion coefficient of the polymer plate (i.e. the slope of the curves drawn on Fig. 1) varies at the beginning of the diffusion process (i.e. the slope of the curves are not independent from the ratio $\frac{\sqrt{t}}{e}$ anymore when at the initial stage, when $t$ tends towards 0). Thus, a sort of delay time is predicted, during which
the instantaneous moisture diffusion coefficient varies. This short period of time is followed by a pseudo-fickian diffusion regime with a constant apparent diffusion coefficient.

These discrepancies significantly increase with the coefficient of moisture expansion. Eventually, the coefficient of moisture absorption affects both the transient and permanent stages of the diffusion process predicted according to the multiphysics model.

## 6.2 Unsymmetrical moisture load

Let us consider the case when the opposite surfaces of the plate are submitted to different relative humidity levels. The environmental conditions correspond to a reference moisture content level $C_0$ of 5% on the first side of the plate, whereas it is either equal to 0%, 2.5% or 5% on the second side.

The moisture diffusion process is computed through Eq. (50), assuming the polymer to present a typical coefficient of moisture expansion: $\eta = 0.6$. The obtained results for the volume average of the moisture content are shown as a function of the classical $\sqrt{t}$ ratio on Fig. 2 below.

<figure><img src="image_2.png" /><figcaption>**Fig. 1** Effect of CME on moisture sorption ($P_{ex}$ = 1 MPa).</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Fig. 2** Average moisture content predicted in a polymer plate submitted to unsymmetrical humid conditions. Cs2 stands for the moisture content reference level assumed to be applied on the second side of the plate.</figcaption></figure>

## 7 Conclusions

This work is focused on developing an enhanced version of the model describing the diffusion of moisture in polymers based on the so-called thermodynamical approach first introduced by Derrien and Gilormini [6], and then improved by Sar et al. [18]. For the first time, in contrast to both those references, the present paper handles the differential swelling experienced by the polymer during the moisture diffusion process. The effects induced by the through-thickness differential swelling on the time-dependent diffusion are properly taken into account in the mathematical development of the model, through additional terms involving partial derivatives of the volume strain by the moisture content. Obviously, the resulting multi-physics kinetics law changes by comparison with the original (simplified) version of the model. The expressions satisfied by the boundary conditions for the moisture have been determined both for the cases where the material is considered as subjected to symmetrical moisture loads, and in cases when heterogeneous humid conditions are experienced by the polymer structure.

Some preliminary results obtained through computations demonstrate that the developed model enables to predict anomalous (i.e. non Fickian) moisture uptakes. The anomalies of diffusion do mostly take place at both the very beginning of the diffusion process and the permanent regime. Non-linearities of the weight gain are thus predicted when the moisture sorption starts. After a short time, these non-linearities vanish, so that a pseudo-fickian moisture uptake follows. This pseudo-fickian regime corresponds to an instantaneous moisture diffusion coefficient independent from the $\sqrt{t}$ ratio.
At the end of the process, the permanent regime is characterized by a maximum moisture absorption capacity, the value of which depends on materials properties such as the coefficients of moisture expansion, as an example.

Further work will be dedicated to a thorough investigation of this new version of the model through extensive numerical tests.

Future work will focus on further enhancements, such as accounting for reversible plasticization effects experienced by the polymer during the moisture diffusion process, (i.e. the reduction of the material stiffness induced by the presence of water).

## References

[1] E. C. Aifantis and W. W. Gerberich, Gaseous diffusion in a stressed-thermoelastic solid. Part II: The thermomechanical formulation, Acta Mech. 28, 1-24 (1977).

[2] E. C. Aifantis and W. W. Gerberich, Gaseous diffusion in a stressed-thermoelastic solid. Part II: Thermodynamic structure and transport theory, Acta Mech. 28, 25-47 (1977).

[3] L.-W. Cai and Y. Weistman, Non-Fickian moisture diffusion in polymeric composite, J. Compos. Mater. 28, 130-154 (1994).

[4] R. G. Carbonell and G. C. Sarti, Coupled deformation and mass-transport processes in solid polymers, Ind. Eng. Chem. Res. 29, 1194-1204 (1990).

[5] H. G. Carter and K. G. Kibler, Langmuir-type model for anomalous moisture diffusion in composite resins, J. Compos. Mater. 12, 118-131 (1978).

[6] K. Derrien and P. Gilormini, The effect of moisture-induced swelling on the absorption capacity of transversely isotropic elastic polymer-matrix composites, Int. J. Solids Struct. 46, 1547-1553 (2009).

[7] S. Fréour, F. Jacquemin, and R. Guillén, On an analytical self-consistent model for internal stress prediction in fiber-reinforced composites submitted to hygroelastic load, J. Reinf. Plast. Compos. 24, 1365-1377 (2005).

[8] M. Gigliotti, F. Jacquemin, J. Molimard, and A. Vautrin, Transient and cyclical hygrothermoelastic stress in laminated composite plates: Modelling and experimental assessment, Mech. Mater. 39, 729-745 (2007).

[9] F. Jacquemin and A. Vautrin, Modelling of the moisture concentration field due to cyclical hygrothermal conditions in thick laminated pipes, Eur. J. Mech. A, Solids 21, 845-855 (2002).

[10] F. Jacquemin, S. Fréour, and R. Guillén, Prediction of local hygroscopic stresses for composite structures – analytical and numerical micromechanical approaches, Compos. Sci. Technol. 69, 17-21 (2009).

[11] J. Jedidi, F. Jacquemin, and A. Vautrin, Design of accelerated hygrothermal cycles on polymer matrix composites in the case of a supersonic aircraft, Compos. Struct. 68, 429-437 (2005).

[12] F. C. Larché and J. W. Cahn, A linear theory of thermochemical equilibrium of solids under stress, Acta Metall. 21, 1051-1063 (1973).

[13] F. C. Larché and J. W. Cahn, The effect of self-stress on diffusion in solid, Acta Metall. 30, 1835-1845 (1982).

[14] P. Neogi, K. Myungsoo, and Y. Yang, Diffusion in solids under strain, with emphasis on polymer membranes, AIChE J. 32, 1146-1157 (1986).

[15] B. P. Patel, M. Ganapathi, and D. P. Makhecha, Hygrothermal effects on the structural behaviour of thick composite laminates using higher-order theory, Compos. Struct. 56, 25-34 (2002).

[16] R. M. V. G. K. Rao, H. V. Shylaja Kumari, and K. Suresh Raju, Moisture behavior of T300-914C laminates, J. Reinf. Plast. Compos. 14(5), 513-522 (1995).

[17] S. Roy, Modeling of anomalous diffusion in polymer composites: A finite element approach, J. Compos. Mater. 33(14), 1318-1343 (1999).

[18] B. E. Sar, S. Fréour, P. Davies, and F. Jacquemin, Coupling moisture diffusion and internal mechanical states in polymers – A thermodynamical approach, Eur. J. Mech. 36, 38-43 (2012).

[19] R. Selzer and K. Friedrich, Mechanical properties and failure behavior of carbon fiber-reinforced polymer composites under the influence of moisture, Compos. A, Appl. Sci. Manuf. 28, 595-604 (1997).

[20] C. H. Shen and G. S. Springer, Moisture absorption and desorption of composite materials, J. Compos. Mater. 10, 2-20 (1976).

[21] I. Verpoest and G. S. Springer, Moisture absorption characteristics of aramid-epoxy composites, J. Reinf. Plast. Compos. 7(1), 2-22 (1988).