On shape stability of panel paintings

exposed to humidity variations

A numerical study between Science &
Art
Part 1: Modelling isothermal moisture movement
S. Reijnen and A.J.M. Jorissen
Eindhoven University of Technology, the Netherlands

A substantial part of museums collection is painted on wood, so called panel paintings. In an

attempt to protect this part of our cultural heritage against degradation, museums apply strict
humidity requirements mainly based upon practical experience: Relative Humidity RH = 55%
5.0%. Today scientist and conservators ask themselves: Are these strict requirements really
needed and what are the consequences of changing these strict requirements?. Finite
Element Analysis could provide an answer and possible lead to a more permanent solution.
Due to its hygroscopic behaviour, wood is sensitive to variations in humidity. As a
consequence of changing environmental conditions and the hygro-expansional behaviour of
wood it tends to deform. When these deformations exceed the elastic limit, it could lead to
permanent deformations. Simulation of this behaviour is quite complex, due to wood being a
heterogeneous, hygroscopic, cellular and anisotropic material. Without numerical simulation,
it is almost impossible to predict the deformation when exposed to moisture variations. This
paper discusses the development of a constitutive model for simulations, the analysis
procedure to solve a moisture diffusion problem and the influence of material properties and
coating layers on the shape stability of sawn wood.

Introduction
It is important to preserve our cultural heritage against degradation. Different materials
demand different conservation methodologies. Perception, based on years of experience,
plays an important role in making the decision on which method of conservation is to be
applied. History learns that this approach does not (always) deliver the desired results.
HERON Vol. 58 (2013) No. 2/3

119

Despite the high value of this experience and empirical approach, it is not (always)
sufficient. Thorough knowledge of the more or less changing environmental factors,
mainly characterized by temperature, relative humidity, light and even more the effect of
these changes on the internal stresses and deformations in the wooden panel, as well as the
response of the paint layers, could lead to a more permanent solution.

A panel painting is painted on wood. Not many people are aware that famous paintings,
painted by great masters, are painted on wood. The Mona Lisa by Leonardo da Vinci is a
good example. The Mona Lisa is painted on a wooden panel made from poplar and is
almost 500 years old. Wood, and the better known canvas, varies widely in material
properties.

In order not to expose the artefacts to big changes in temperature and relative humidity,
museums, invest in climatic control systems. Strict requirements, completely based on
empirical evidence, should protect the art against degradation. Besides the cultural
responsibility, there is also a financial stimulation to review the current way of thinking
and to try to provide a scientific basis on acceptable climate variations.

Wood is sensitive to fluctuations in relative humidity. Absorption from and desorption of

moisture to the immediate environment is unfortunately not without conflict. Absorption
and desorption of moisture give rise to swelling and shrinking in radial and tangential
directions. The shrinkage and swelling values in tangential direction are roughly twice the
values in radial direction, as indicated in figure 3, resulting in unwanted shape
deformation and strain development. Deformations exceeding the elastic range can lead to
permanent deformation and/or cracking as shown in figure 1.

The objectives of this research are: (1) to examine the possibility of ABAQUS standard to
model moisture flow and shape stability using the procedure of heat conduction instead of
mass diffusion within a multi-physical environment, see also Mirianon, Fortino and
Toratti (2008), (2) to examine the influence of changing environmental conditions and
changing material properties on shape stability and the related tress field, see also
Ormarsson, Dahlblom and Petersson (1997), (3) to examine the influence of a gesso coating
layer (mixture of hide glue, gypsum or sometimes ground chalk and water for smoothing a
panel surface) on shape stability.

120

The content of this paper is limited to the development of a constitutive model for
simulations, to solve moisture diffusion and the influence of gesso on shape stability of
sawn wood.

Modelling isothermal moisture movement in wood, using ABAQUS

transient heat conduction
When a wooden object is exposed to varying environmental conditions, especially
fluctuations in relative humidity, the wood moisture content also varies due to the
hygroscopic behaviour of wood. Especially fluctuations in relative humidity below the
fibre saturation point result in deformation in the form of shrinking or swelling. If these
deformations do not exceed the elastic range, theoretically there is no problem. Beyond this
elastic range, the so called plastic range, the material no longer returns into the initial state
after removing of the load. This is called plastic deformation. Jakiela, Bratasz and
Kozlowski (2008) developed a numerical model describing moisture movement due to
changing environmental conditions to calculate the related stress field in the elastic range.
This model has been applied to lime wood cylinders. Schellen and Schijndel (2011) verified
the work done by Jakiela, Bratasz and Kozlowski (2008) with help of a numerical model
developed in COMSOL, a finite element package designed to solve building physics
problems.

To validate the outcome by Jakiela et al. and Schellen et al., an ABAQUS CAE transient heat
conduction analysis is performed using a transient heat conduction analysis procedure
instead of mass diffusion analysis. One of the possible great advantages of using ABAQUS
CAE (Complete ABAQUS Environment) heat conduction analysis is the ease of modelling.
Within the CAE, it is possible to perform a heat conduction analysis and apply the
outcome as a predefined field to a static stress/strain analysis. This is a so called
sequentially coupled thermal stress analysis. A sequentially coupled thermal stress
analysis can be used when stress/displacement is dependent on a temperature field and
there is no inverse dependency.

ABAQUS uses Fouriers law of heat conduction to analyse heat transfer, see equation (1),
(2) and equation (3) to perform a mass diffusion analysis. Equation (3) is an extension of
Ficks law, see equation (4) and (5). The difference can be found in the fact that the
equations used by ABAQUS allow for a non-uniform solubility (the ability of a liquid
121

Figure 1: The boy from AlFayum;

2nd century; Encaustic panel painting; The figure shows local cracking in the painting

(solute) to dissolve into a solid (solvent)) of the solute through the solvent and for mass
diffusion driven by gradients of temperature and pressure. Ficks first law of mass
diffusion is a linear equation (4). The extended law used by ABAQUS, equation (3),
becomes non-linear since the diffusion coefficient, the Soret factor s and the pressure
stress factor p depend on the concentration. Due to the analogy between Fouriers
equations and Ficks equations, mass diffusion can be modelled using a heat transfer
analysis and vice versa.

qx = k

122

T
Fouriers law (steady state)
x

(1)

T
2T
=
Fouriers law (transient)
t
x 2

(2)

J x = sD + s (ln(T T Z )) + P
x
x
x

(3)


Ficks first law (steady state)
x

(4)

2
Ficks second law (transient)
=D
t
x 2

(5)

J x = D

qx ... heat flux [W/m]

k thermal conductivity [W/(mK)]

T temperature [K]

c p
thermal conductivity (k) [W/(mK)]]

thermal diffusivity [m/s] =

density [kgm]

c p specific heat capacity [J/(kgK)]

J x diffusion flux [mol/(ms)]

D diffusion coefficient [m/s]

concentration [mol/m]

D(c, T, f) .. diffusivity [m/s]

s(T, f) solubility [ppm]
K s (C , T , f ) .... Soret factor, providing diffusion because of a temperature gradient [-]

T ... temperature [K]

T Z .. value of the absolute zero on the temperature scale used [-]

K p (C , T , f ) .... pressure stress factor, providing diffusion, driven by the gradient of the
pressure stress [-]
C ... concentration of the diffusing material [kg/m]
f . other predefined field variables (potential) [-]
Although mass diffusion can be modelled using a transient heat transfer analysis, heat
transfer analysis and mass diffusion analysis are not the same. For example, heat transfer
analysis based on Fouriers law can only use a temperature gradient as the driving force
behind the diffusive process. Within mass diffusion analysis other driving forces
(ABAQUS calls this chemical potentials) like pressure, temperature and concentration can
control the diffusive process. ABAQUS has developed a strong capability over a long time
to solve multi-physics problems. The advantage of ABAQUS multi-physics is the ease with
which multi-physics problems can be solved. It has the ability to utilise the same model,
element library, material data and load history. A single-physics analysis can easily be
123

extended to a multi-physics analysis, without the need for additional tools, interfaces or
simulation methodology.
2.1

description

What is happening in a lime wood cylinder exposed to relative humidity fluctuations? A

reproduction of the research done by Jakiela, Bratasz and Kozlowski (2008).
2.1.1 Geometry (figure 2)

= 0.13 m lime wood

Figure 2: Geometry of the lime wood model

2.1.2 Plane strain situation

Strain directed perpendicular to x-y plane equals zero

normal strain z = 0
shear strain xz = 0
shear strain yz = 0
The wood reaction to moisture parallel to the grain is neglected, the surface can be
regarded as fully closed.
2.1.3 Thermal boundary conditions (isothermal)

t<0

124

TS = T0

t>0

q = h (T0 TS )

t=

T = T0

T0 ... 20 [C]

T .. 20 [C]
TS .. surface temperature [C]

q ... heat flux at surface [W/m]

h heat transfer coefficient [W/(mK)] h =
Q .. heat flow [W]

Q
A T

A .. surface area [m]

2.1.4 Hygric boundary conditions

t<0

uS = u0 = 14%

t>0

g = (u0 uS )

t=

u = 6%

u0 = 14% RH0 = 70%

u = 6% RH = 30%

g . moisture flux at surface [kg/(ms)]

. moisture transfer coefficient [kg/(ms)]
2.1.5 Mechanical boundary conditions

The model surface is constraint (type: coupling) to a reference point in space. The
constraint degree of freedom is UR3 (vertical movement). By coupling the surface to a
reference point in space, rigid body rotation cannot take place and the model is free to
move in radial direction.
2.1.6 Mechanical stress and strain (generalized Hookes law)

E
x
x

yx
y =
Ex

xy

xy
Ey
1
Ey
0

x x

x

0 y + y + y w

0
xy

0
1

Gxy

(6)

125

 x , y ... normal strain components [-]

xy .... shear strain component [-]

Ex , Ey .. moduli of elasticity [N/m]

xy , yx ... Poissons ratios [-]

Gxy .... shear modulus [N/m]
x , y ... normal stress components [N/m]
x , y ... linear thermal expansion coefficient [1/K]
.. temperature increment [K]
w .. moisture content increment [kg/m]

x , y ... linear relative deformation due to changing moisture content [1/(kgm)]

2.1.7 Material properties

For such a complicate material as wood, one should not expect ideal elastic behaviour, as
described by Hookes law. The stress-strain diagram is therefore not the same as for an
ideal elastic body. Ideally behaving materials seem hard to be found in nature, but a
material like rubber approaches this behaviour quite well. Table 1 shows the modulus of
elasticity of lime wood used to model the lime wood cylinder. Notice that the modulus of
elasticity in tangential direction is approximately half the value of the modulus of elasticity
in the radial direction.
It is known that there is a difference between tangential and radial shrinkage and that
characteristic cracks can develop when wood dries. There is of course more than one
reason for this behaviour. Different wood species have different properties, and show
different shrinkage and swelling behaviour. One thing is always the same for all species:
wood shrinks the most in the tangential direction, about two times more than in the radial
direction, see figure 3.
Table 1: Modulus of elasticity of lime wood (Jakiela, Bratasz and Kozlowski (2008))

126

RH [%]

Tangential direction [MPa]

Radial direction [MPa]

20

600

1120

35

490

900

50

450

820

65

420

770

dimensional change [%]

equilibrium moisture content [%]

Figure 3: Relation between the Equilibrium Moisture Content (EMC) and dimensional change for
shrinkage coefficients R = 0.13 (radial direction) and T = 0.28 (tangential direction).

Between two systems, whose concentration of moisture content differ, there is a natural
tendency for of moisture transfer. Due to moisture transfer, both systems seek for
equilibrium, minimizing the difference in moisture concentration, this process is called
diffusion. Table 2 represents the moisture diffusion coefficient of lime wood as a function
of equilibrium moisture content.

The Equilibrium Moisture Content (EMC) as a function of Relative Humidity (RH) and
temperature T is calculated according to equation (7) (Kollmann and Cote (1968)

EMC =

K K (RH) + 2K 1K 2 K 2 (RH)2
1800 K (RH)
+ 1

W 1 K (RH) 1 + K1K (RH) + K1K 2 K 2 (RH)2

(7)

with:
W = 345 + 1.29T + 0.0135 T 2
K = 0.805 + 0.000736 T 0.00000273 T 2

(8)
(9)

K 1 = 6.27 0.00938 T 0.000303 T 2

(10)

K 2 = 1.19 + 0.0407 T 0.000293 T 2

(11)

127

Table 2: Moisture diffusion coefficient as a function of equilibrium moisture content (EMC) (Jakiela,
Bratasz and Kozlowski (2008))
equilibrium

diffusion coefficient

moisture content

radial

tangential

longitudinal

[%]

[m/h]

[m/h]

[m/h]

0.0003888

0.0003888

0.0009000

0.0004751

0.0004751

0.0050400

5.5

0.0004841

0.0004841

0.0053500

0.0005137

0.0005137

0.0056700

8.5

0.0005461

0.0005461

0.0058500

0.0005572

0.0005572

0.0056700

13.5

0.0006690

0.0006690

0.0045400

18

0.0008026

0.0008026

0.0030700

23

0.0009690

0.0009690

0.0021000

28

0.0012029

0.0012029

0.0013500

equilibrium moisture content [%]

relative humidity [1/100 %]

Figure 4: Equilibrium moisture content EMC as a function of relative humidity RH at 10 C,

20 C, 30 C and 40 C
2.2

Results

At time t = 0, the wooden cylinder is in equilibrium with its surrounding. The relative
humidity at t = 0 equals 70%, resulting in a wood equilibrium moisture content of 14%, see
also figure 4. At t = 1 the surrounding conditions suddenly change. The relative humidity
128

at t = 1 drops from 70% to 30%, corresponding to 6% wood equilibrium moisture content.

After this event the environment is kept constant at RH1 = 30%. From t = 1 the wooden
cylinder is slowly drying as it releases moisture. From t = 1 the wooden cylinder seeks for
a new equilibrium with its surrounding condition. With help of ABAQUS finite element
model the time to reach complete equilibrium was calculated at 40 days.
2.2.1 Distribution of moisture content after 24 hours

After 24 hours the centre of the cylinder is still at the initial 14% moisture content. The
surface of the wooden cylinder shows a fast transition in moisture content, see figure 5,
which shows the changing moisture content at different distances from the surface as a
function of time. The selected distances from surface up to 10 mm inside the wooden
cylinder show a non-linear diffusion process. The first 1 mm to 5 mm from the surface
level instantaneously changes its moisture content. Figure 5 also shows that the core of the
cylinder lying deeper than 1 centimetre does not experience any change in the moisture
content in the first 3 hours.

equilibrium moister content [%]

time [s]
Figure 5: Distribution of moisture content at selected distances from the surface up to 10 mm into a
wooden cylinder with a step change of 14% MC to 6% MC which is equal to 70% RH to 30% RH
after 24 hours
129

2.2.2 Strain and stress development in tangential direction

Figure 6 shows the development of radial strain in compression perpendicular to the grain
at different depths from the surface as a function of time. Generally, the maximum elastic
strain at which wood starts to deform plastically perpendicular to the grain lies around the
0.004 (Ormarsson et al. (1997).

strain [-]

time [s]
Figure 6: Radial strain at selected distances from the surface up to 10 mm into a wooden cylinder
with a step change of 14% to 6% moisture content (MC) equal to 70% to 30% relative humidity
(RH) after 24 h
Figure 7 shows the development of tangential stress at different depths from the surface as
a function of time. The elastic range up to a tangential stress of 2.5 MPa and maximum
tangential strength of 5.5 MPa as denoted by Jakiela et al. (2008) is exceeded. With
continues drying of the interior layers, the stress slowly decreases. This slow decrease is
the result of the slow vanishing of the moisture gradient as the interior layers dry,
resulting in more evenly shrinkage. The content of this paper is limited to stress
development in tangential direction. Information about radial stress development,
tangential and radial strain development can be found in Jakiela et al. (2008) and Reijnen
(2012).

130

2.2.3 Verification

Comparison of results is part of finding the right environment for modelling. The results of
the ABAQUS model are compared to the results from Jakiela et al. (2008) and Schellen et al.
(2011). The results found by Schellen et al. (2011) were computed with help of a numerical
model within COMSOL as a comparative benchmark of the research done by Jakiela et al.
(2008). Comparing the results found by Jakiela et al. (2008) and Shellen et al. (2011), it can
be concluded that these results are very much the same as the results calculated by
ABAQUS, see Reijnen (2012).

Summary and conclusion

Within ABAQUS finite element software, it is shown that due to the analogy between
Fouriers law of heat conduction and Ficks law of mass diffusion, it is possible to solve a
moisture movement problem, which is a mass diffusion problem, with a heat transfer
analysis. Careful implementation, proper material data and using an appropriate driving
potential deserves full attention. Wood is most often only part of a bigger construction,
consisting of many other layers and different materials as in panel paintings. Different
materials come with different material properties. Because of this different behaviour, it is

stress [MPa]

time [s]
Figure 7: Tangential stress at selected distances from the surface up to 10 mm into a wooden
cylinder with a step change of 14% to 6% moisture content (MC) equal to 70% to 30% relative
humidity (RH) after 24 h
131

wise to use a driving potential, which is consistent for different materials. Water vapour
pressure P and water vapour content w are such driving potentials. Unfortunately,
moisture transport below fibre saturation point cannot be regarded as being a pure Fickian
process. In case of short time moisture movement simulations, one should use a nonFickian or multi-Fickian model.
A Lime wood cylinder exposed to a change in relative humidity from 70% to 30% shows
that moisture distribution is strongly non-linear. The first few millimetres from surface
level respond very quickly to changing relative humidity. As a result, cracking of the
surface level would be almost instantaneous. Practical research and experience has proven
that this is not the case. As a result, a non-Fickian or multi-Fickian model must be used for
describing the process of moisture distribution through wood more realistically, see figure
21. This non-Fickian or multi-Fickian model predicts a much less steeper gradient at
surface level, corresponding to practical experience.

Figure 21: Stress development due to changing wood moisture content simulated with so-called
Fickian and Non-Fickian model descriptions

132

Literature
Jakiela S., Bratasz L., Kozlowski R. (2008) Numerical modelling of moisture movement and
related stress field in lime wood subjected to changing climate conditions, Wood Science
Technology, 42,21-37
Kollmann F., Cote W. (1968) Principles of Wood science and Technology: 1 Solid Wood. Institute
fur Holzforschung und Holztechnik , Munchen, Germany / College of Forestry,
Syracuse, New York.
Mirianon F. Fortino S. and Toratti T. (2008) A method to model wood by using ABAQUS
finite element software, VTT (Technical Research Centre of Finland) Publications,
Laskut, Finland.
Ormarsson S., Dahlblom O., Petersson H. (1997) A numerical study of the shape stability of
sawn timber subjected to moisture variations, Lund Institute of Technology.
Reijnen S. (2012) Moisture transport and shape stability of wood exposed to humidity variations
A numerical study. Eindhoven University of Technology, Department of Building and
Architecture (master diploma thesis), Eindhoven, the Netherlands.
Schellen H., van Schijndel J. (2011) Numerical modelling of moisture related mechanical
stress in wooden cylindrical objects using COMSOL: a comparative benchmark,
Eindhoven University of Technology, Department of Building and Architecture,
Eindhoven, the Netherlands.

133

134