Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, February 2006 / Copyright © 2006 Japan Concrete Institute147

Scientific paper

Time-Dependent Structural Analysis Considering Mass Transfer to
Evaluate Deterioration Process of RC Structures
Hikaru Nakamura1, Worapong Srisoros 2, Ryosuke Yashiro3 and Minoru Kunieda4
Received 24 August 2005, accepted 1 November 2005

Abstract
A time-dependent structural analysis method under multi actions in consideration of drying shrinkage due to moisturetransfer and rebar corrosion due to chloride ions penetration as well as external load actions was developed. The
Rigid-Body-Spring Networks (RBSN) model and the truss networks model were used for structural analysis and masstransfer analysis, respectively. In addition, mass transfer through bulk concrete and mass transfer through cracks by setting truss networks on the boundaries of Voronoi particles, was also considered. The developed method was confirmed tosimulate well the deterioration process due to mass transfer for initial cracking behavior and ultimate behavior of concrete
structures.

1. Introduction
Cracks in concrete often occur due to the internal stressinduced by mass transfer such as moisture and heat or bycorrosion expansion of rebar under no external load.
Once cracks occur, mass transfer accelerates through thecracks, causing gradual degradation in structural performance and durability. It is, therefore, important toclarify the effect of mass transfer on structural performance and durability in order to assess the service lifeof concrete structures.
Although the effects of internal stress due to masstransfer have been estimated analytically, the main purpose was to simulate the occurrence of cracks throughtime. Further, the cracking behavior and other structural
performance characteristics such as ultimate behavior
have hardly been estimated. On the other hand, althoughnumerical methods for predicting the structural performance of concrete structures before and after the peakload have been advanced, most methods did not consider
the continuous effects of initial damage due to masstransfer.
This paper proposes a unified analytical method withload applying analysis for short term behavior and masstransfer analysis for long term behavior under multiactions in consideration of drying shrinkage due tomoisture transfer and rebar corrosion due to chloride ion

1Professor, Department of Civil Engineering, Nagoya
University, Nagoya, Japan.
E-mail: hikaru@civil.nagoya-u.ac.jp
2
Doctoral student, Graduate School, Department of Civil
Engineering, Nagoya University, Nagoya, Japan.
3
Master of Engineering, Graduate School, Department ofCivil Engineering, Nagoya University, Nagoya, Japan.
4
Associate Professor, Department of Civil Engineering,
Nagoya University, Nagoya, Japan.

penetration. The Rigid-Body-Spring Networks (RBSN),
developed by Kawai (Kawai 1978) is used for structural
analysis and the truss networks model is applied for masstransfer analysis. Several analyses in two-dimensions
were performed to confirm the effectiveness of the developed method and to simulate the deterioration process
induced by mass transfer.

2. Structural analysis
Nowadays, numerical methods for predicting the structural performance of concrete structures are becoming
more effective owing to advances in computing technology. The Finite Element Method (FEM) has generally
been applied for estimating mechanical behavior of
concrete structures. This method provides reasonableprediction for the totality of concrete structure behavior,such as the formation of cracks and load bearing capacity.
The deterioration process of concrete structures greatly
depends on cracking. Therefore, a model that can accurately simulate crack formation and propagation is required for predicting both the structural performance anddurability of concrete structures. Discrete method isbetter than the continuum models for modeling of material discontinuity. From this viewpoint, theRigid-Body-Spring Networks approach, which is one ofthe discrete approaches, is used for structural analysis,
since this method is an analytical technique based on
discrete mechanics that easily deals with crack propagation of concrete directly (Saito et al. 1999).
2.1. Rigid-body-spring networks
Rigid-Body-Spring Networks represent a continuum
material as an assembly of rigid particle elements interconnected by zero-size springs along their boundaries.
The response of the spring model provides comprehension of the interaction between particles instead of theinternal behavior of each particle based on continuum

148

kθ

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

ks

Voronoi meshknｈFig. 1 Rigid-Body-Spring Networks.

mechanics. Each rigid particle has two translations andone rotational degree of freedom defined at the nuclei
within. The interface between two particles consists ofthree individual springs in the normal, tangential androtations, as shown in Fig. 1. Since concrete cracks initiate and propagate along interparticle boundaries, thecrack pattern is strongly affected by the mesh design.
Therefore, random geometry using Voronoi diagrams(Bolander et al. 1998) is applied to partition the material
onto an assembly of rigid particle. The random geometryof the networks does not represent any structural feature
within the concrete material, but rather is used to reducemesh bias on potential crack directions.

σ
ft

E
ft1

2.2. Concrete material models
The fracture criterion in RBSN is not based on a tensorial
measure of stress, but it uses the average stresses acting
normally and tangentially to the particle interface. Normal springs are set to represent the tensile and compressive behavior of concrete.
In this study, the tensile behavior of concrete up to thetensile strength is modeled by using a linear elastic, whilea bilinear softening branch of 1/4 model is assumed after
cracking, as shown in Fig. 2 (a) which is represented by
the tensile strength, ft, the tensile fracture energy, Gft, andthe distance between nuclei, h. In this paper, Gft is set to
0.1 N/mm. The behavior of concrete under compressive
strength f’c is modeled using a parabolic curve up to thecompressive strength and a linear softening branch isassumed thereafter, as shown in Fig. 2 (b). The slope ofthe linear softening branch is defined by considering thecompressive fracture energy, Gfc, of Equation (1) (Nakamura et al. 2001) to avoid mesh size dependence aswell as tension behavior.
(1)

G fc = 8.8 f c'Tangential springs represent the shear transfer mechanism of uncracked and cracked concrete. The shear
strength is assumed to the Mohr-Coulomb type criterionwith tension and compression caps, as shown in Fig. 3(a). The cohesion, c, and internal friction angle, φ, are set
to 0.138f'c and 37°, respectively (Ueda et al. 1984). It is
assumed that the shear strength is constant when normal
stress is greater than ψf'c, with ψ set to 0.5. After shear
stress reaches the yield strength, the stress moves on theft1 = 0.25ft
εt1 = 0.75Gf/fth
εt1 = 0.75Gf/fth

τ
φ

Gft/h

εt0

c

ψf’c

εt1

εtu

ε

σ

f’c

ft

a) Tension (1/4 model)a) Mohr-Coulomb type criteria

σ

Gcr

f’c

G

Ec Ec
ε0/2 ε0

Gfc/h
εcu

ε

wmax wb) Compression (parabolic)b) Reduction of shear stiffness

Fig. 2 Stress-strain relationship of concrete.

Fig. 3 Shear model in concrete.

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

yield surface until the shear strain reaches the ultimate
strain, γu. The force in the shear spring is released and the
local stiffness, kt, is set to zero when the shear strain
exceeds the ultimate strain. The ultimate strain is set to0.004 µ in this study. The shear transfer capacity at crack
interfaces depends on crack opening. Thus, the shear
stress, τ, is calculated by Equation (2) (Saito et al. 1999)
with the function of the strain of normal spring, ε,Young’s modulus, E, and shear strain, γ.

τ = (1 − λ )βGγ

149

τ
τmax

τmax/100.2

s (mm)

0.4(2)
Fig. 5 Bond stress-slip relationship.

where
0


λ=
 γ
 k

1 − 0 exp  − (γ − γ 0 ) 
γ
γ
 u




1

β =
ε t 0 0.3

(ε − ε t 0 ) exp −
ε
 ε
 tu


for γ ≤ γ 0

for γ > γ 0for ε ≤ ε t 0
for ε > ε t 0

where G = shear stiffness and γ0 = τf/G. The λ and βparameters represent the degree of damage and reduction
factor, respectively, as shown in Fig. 3 (b).
2.3. Reinforcement model
Each reinforcement element is represented by a series ofregular beam elements. The beam nodes are attached tothe concrete particles without regard to the concrete
mesh design through a zero-size link element, as shownin Fig. 4 (Saito et al. 1997). The stress-strain relation ofreinforcement is used as a bilinear model. The bondinteraction between the concrete and reinforcement
strongly affects crack development. The bond stress-slip
relation is introduced into the spring parallel to the rebarof the linked elements. Figure 5 shows the relationshipmodeled by Equation (3) (Shima et al. 1989) up to thebond strength and the same curve as that proposed by

CEB-FIB is assumed thereafter (CEB-FIB 1990). Shimaet al., however, proposed Equation (3) for deformed bars
embedded in massive concrete. For application to reinforced concrete beams, which have smaller cover concrete compared to the cover of massive concrete, thereduction factor due to the effect of the concrete covering
should also be considered (Kitoku et al. 2005). In thispaper, this reduction factor, α was set to 0.6.

τ = α 0.9 f c'

2/3

{1 − exp(− 40(s / D ) )}0 .5

where s = slip and D = diameter of reinforcement.
Concrete corrosion decreases the effective area of thereinforcement. Furthermore, the bond property between
the reinforcement bars and the surrounding concrete is
influenced by the corrosion of rebars. However, the reductions in effective area and the change of bond property have not been considered in this paper.

3. Mass transfer analysis
3.1. Truss networks model
Mass transfer is a continuous flow and it is usually analyzed using the same continuum model as that use for the
Finite Element approach, which is represented by a partial differential equation, whereas RBSN is used forstructural analysis, which does not require continuality.
In this section, the truss networks model is used for masstransfer analysis to obtain the initial strain for structural
analysis by using RBSN with Voronoi diagrams. Each of
the Voronoi elements is linked by truss elements with the

Zero size link
Beam element

Reinforcement
direction

Fig. 4 Reinforcement arrangement.

(3)

Node
Truss

Fig. 6 Truss networks model.

150

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

nodes at the Voronoi nuclei and the intermediate points ofparticle boundaries, as shown in Fig. 6. Then, a simplified one-dimensional diffusion equation using truss
elements is employed to carry out mass transfer. A masstransfer analysis based on RBSN and the truss networks
model, in which truss elements are generated between
the Voronoi nuclei, was proposed by Bolander (Bolander
et al. 2003). This is one of the differences with the proposed model. A feature of this study is that it assumes the
position of the node of truss elements on intermediate
points of the particle boundaries, as mentioned in section
3.5.
3.2. Basic concept of diffusion equation
In the truss networks model, nonstationary potential flow
problems are governed by Equation (4) below, aone-dimensional diffusion equation.

∂φ ∂ ⎛ ∂φ ⎞= ⎜D ⎟∂t ∂x ⎝ ∂x ⎠(4)

where φ = some potential quantity such as the relative
humidity and D = diffusivity.
D is variously referred to in the case of Fourier’s,
Fick’s and Darcy’s law, depending on the types of diffusivity. This can be expressed as a matrix,⎧ ∂φ ⎫⎧ φ1 − φ1 ⎫ ⎧0⎫
AD ⎡ 1 − 1⎤ ⎧φ1 ⎫ 1 AL ⎡2 1 ⎤ ⎪ 1 ∂t ⎪+⎨ ⎬⎨∂φ ⎬ + A⎨
⎬=⎨ ⎬
⎥
⎢
⎥
⎢
L ⎣− 1 1 ⎦ ⎩φ2 ⎭ ω 6 ⎣1 2⎦ ⎪ 2 ⎪φ2 − φ2 ⎭ ⎩0⎭
⎩
∂t ⎭
⎩

(5)
where φ1, φ2 = potential quantity at each node, A = area of
each truss element which is equal to the area of the corresponding Voronoi facet, and L = length of each truss
element. The parameter ω in the second term is added toconvert the volume of Voronoi to the volume of truss
elements due to the overlap volume of adjacent truss
elements, as shown in Fig. 7 and Equation (6). Thisparameter is set as 2.0 for two-dimensional arrangement.
Truss element
Voronoielement

Overlap area

Effective area oftruss element
Fig. 7 Definition of overlap area and parameter omega
(ω).

n

ω=

∑A ⋅L
i =1

i

i

Vreal

(6)

where n = total number of elements, Ai = cross-sectional
area of each truss element, Li = length of each truss
element, i = element number and Vreal = volume of thestructure.
When several mass transfers are considered, they haverelations to each other. Thus, the truss networks for eachmass transfer are superposed and the various properties
relate to each other. For example, the effective diffusion
coefficient of chloride ions in the truss networks is assumed to depend on the relative humidity obtained fromthe truss networks of moisture transfer.
3.3. Moisture transfer
External drying shrinkage generally occurs in concrete
members as soon as the formwork is removed. The difference between water loss by external drying andchemical hydration is, therefore, extreme when water
diffusion of early age concrete is the highest. Since thearea where the moisture content decreases shrinks, internal stress is induced, which in turn causes tensile stressat the surface of the specimen depending on the internaland external restraints. A diffusion equation using
Darcy’s law can describe the moisture transfer as follows.
∂R ∂ ⎛∂R ⎞= ⎜ Dw⎟∂t ∂x ⎝∂x ⎠(7)

where R = relative humidity, t = time and Dw = moisture
conductivity. The equation is solved by the initial condition and the following boundary condition.

∂R+ αm ( Hs + H0 ) = 0∂n(8)where n = normal vector of drying surface, αm = moisturetransmission coefficient, and HS, H0 = relative humidityof the drying surface and atmosphere, respectively. In
this study, Dw is assumed as a constant value for simplicity.
3.4. Chloride ions penetration
Chloride ion penetration is one of the major problems
that affect the durability of reinforced concrete structures.
Although chloride ions in concrete do not directly cause
severe damage to the concrete, they contribute to the
corrosion of embedded rebars in the structure. The general diffusion equation can be written from Equation (4)using Fick’s law, as shown in Equation (9). Thus, it isassumed that the diffusion coefficient changes according
to the relative humidity (Saeki et al. 1996).

151

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

∂C ∂  ' ∂C 
=  Dcl∂t ∂x 
∂x 

(9)

Dcl'
= 0.0032 × 100.025 R
Dcl 0

(10)

where C = concentration of chloride ions, t = time, D'cl =effective chloride ion diffusion coefficient, Dcl0 = diffusion coefficient in water logged condition, and R = relative humidity.
3.5. Mass transfer through crack
The truss networks connected between the Voronoi nuclei and the intermediate points of the particle boundaries
are expressed for the mass transfer within bulk concrete.
Steel corrosion due to mass transfer through cracks,
especially chloride ion diffusion, is, however, more important for structural performance. When surface cracks
occur in reinforced concrete structures due to drying
shrinkage or load application, chloride ions penetrate
into the structure through these cracks. As described in
the previous section, cracks will occur and propagate
along the particle boundaries. To consider mass transfer
through cracks in addition to through the bulk concrete,
the truss networks are generated on particle boundaries inwhich a diffusion coefficient from that of the concretematerial is assumed and the cross sectional area is variedby the crack width, as shown in Fig. 8. The diffusion
coefficient of mass transfer through cracks may be quite
large compared with that of bulk concrete.
Before cracking, the cross sectional area of truss elements on the particle boundaries is 0 and no mass transfer
occurs in these truss elements. After cracking, the crack
width can be calculated directly from the relative displacement between rigid particles. Therefore, the crosssectional area of the truss elements on the particle
boundaries can be evaluated by multiplying the relative
displacement with the thickness of the Voronoi element

and mass can transfer along these truss elements in theconcrete.
3.6. Checking efficiency of truss networks
modelThis section describes the application of the JCMAX
thermal stress analysis program developed by the JCI
Committee on Computer Code Development for Crack
Control in Massive Concrete Using Finite Element
Method (JCI committee, 2005) to analyze heat transfer inconcrete specimens. The results of temperature distribution obtained from JCMAX were compared with theresults of the truss networks model. The thermal diffusion can be described from Equation (11) using Fourier’s
law as follows:∂T ∂  ∂T 
= K 
∂t ∂x  ∂x 

where T = temperature, t = time, and K = temperature
diffusion coefficient.
Figure 9 displays the geometry of a specimen with
1500 Voronoi elements for the truss networks model and
100 rectangular elements for FEM analysis. Table 1 lists
the thermal properties of the concrete (thermal conductivity, heat transfer coefficient, specific heat and massdensity). The initial and surrounding temperatures were
set to 20°C and 0°C, respectively. For the boundary
condition of heat transfer, higher temperatures inside thespecimen can be transferred to the outside only via theheat transfer boundary using the truss networks generTable 1 Thermal properties of concrete.

Thermal conductivity (W/m⋅K)

2.7

2Heat transfer coefficient (W/m K)

14

Specific heat (KJ/kg⋅K)

1.1

3

Mass density (kg/m )Zero area truss onparticle boundary

Adiabatic boundary
Heat transfer boundaryA

before cracking

NodeTruss onboundary
Truss insideof concrete

Area depending on
crack widthFig. 8 Truss networks at crack.

10000 mm

B

after cracking

(11)

2650

Initial: 20℃
Surrounding: 0℃

2000 mm2000 mm

C

1500 elements

10×10 elements

10000 mm
Truss network model

FEM (JCMAC)

Fig. 9 Analytical model.

152

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006
Point C

Temperature (℃）

20

Point B
Finite Element Method
100 elements (JCMAC)
Truss network modelvoronoi 1500 elements, ω = 1
voronoi 1500 elements, ω = 2

10Structural model

Truss model

Point A0

1020
Age (days)

30(a) Different ϖ values

20
Temperature (℃）

Point C

Combined model

Point B
Finite Element Method
100 elements (JCMAC)
Truss network modelvoronoi 1500 elements, ω = 2voronoi 1000 elements, ω = 2voronoi 500 elements, ω = 2

10Point A0

1020
Age (days)

30(b) Different element numbers
Fig. 10 Comparison results.

ated between the Voronoi nuclei.
Figure 10 (a) compares the temperature distribution
obtained from JCMAX and the truss networks model ofpoints A, B and C. Two values of the conversion coefficient, ω (ω = 1 and 2), were used for analysis and comparison with the JCMAX program. The temperature
distributions in which ω is set to 2.0 obtained from thetruss networks model were found to closely agree withthe results from the JCMAX program. Therefore, thetruss networks model with Voronoi random mesh can beapplied to analyze mass transfer for which ω is set to 2.0
for two-dimensional arrangements.
To investigate the effect of the number of elements,
specimens divided into 500, 1000 and 1500 Voronoi
elements were analyzed by using the truss networks
model and setting ω to 2.0. Figure 10 (b) shows thecomparison results between 500, 1000 and 1500 Voronoi
elements for the truss networks model and 100 rectangular elements for the JCMAC program. Based on the
results, a small mesh for calculating mass transfer byusing the truss networks model is necessary to obtain thesame results with FEM due to the difference in shapefunction.

4. Unification with structural analysis andmass transfer analysis
4.1. Internal stress due to drying shrinkage
Internal strain occurs in concrete due to the distribution

Fig. 11 Concept of unified model.

of the mass as a result of mass transfer. In the case ofdrying shrinkage, a linear relationship between thechange in relative humidity, ∆R, and the corresponding
change in shrinkage strain, ∆εn, is assumed:∆εn = αsh∆R(12)where αsh = shrinkage coefficient.
Internal stress also occurs due to internal restraint fromthe internal strain distribution and external restraint fromthe boundary condition. The internal strain of each truss
element obtained from mass transfer analysis is introduced to the corresponding normal springs of RBSN as
structural analysis. Figure 11 shows the unified model ofthe structural analysis and the mass transfer analysis usedin the analysis. The matrix form of structural analysis is
as follows. It solves the initial stress program taking into
consideration time-dependent strain.

∫∫∫[B] [D ][B]{d }dV =T

ep∫∫ [N ] { f }dS + ∫∫∫[B] [D ]{ε }dV
T

T

ep

t

(13)

where [Dep] = elastic-plastic constitutive matrix and εt =time-dependent strain obtained from the mass transfer
analysis.
4.2. Internal stress due to corrosion expansionof rebar
Cracks due to corrosion expansion of rebars are one of
the major causes of deterioration. It is usually assumed
that corrosion starts when the concentration of chloride
ions at the surface of rebars exceeds a threshold value. In
this study, the threshold value is set to 1.2 × 10-6 g/mm3(JSCE 2002) by unit volume of concrete and the corrosion velocity is set to 3.0 × 10-8 g/mm2 day for simplicity.

153

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

Table 2 Concrete properties.

Compelling
displacement

A1

Corrosion
expansion
Fig. 12 Corrosion expansion.

Then, the compelling displacement induced by corrosion
expansion is calculated by the following equation (Matsuo et al. 1997),U=

Wr (dV − 1)Tensile strength, ft (MPa)Young modulus, E (MPa)Moisture diffusivity withinbulk
concrete,Dw
2(mm /day)
Water transfer coefficient,αm (mm/day)
Shrinkage coefficient, αsh
Moisture diffusivity between cracks, Dcr (mm2/day)A2

3.0
2.00 × 104
10
30

30.0
0.001
1.0 × 105 1.0 × 105

0

(14)

ρs

where U = compelling displacement, Wr = corrosion
volume of rebar, ρS = density of steel, and dV = volume
expansion ratio of rust, set to 2.5.
The compelling displacement is translated to the displacement in the direction of normal springs at adjacentpoints, as shown in Fig. 12, and the internal strain is
introduced in the springs which are calculated by dividing the displacement depending on the corrosion expansion by the normal spring length.

5. Crack propagation analysis due to drying
shrinkage
5.1. Analysis outline
Drying shrinkage analysis is performed by considering
the moisture transfer. Figure 13 displays dog bone shapespecimens (160 × 100 × 1000 mm) for free drying
shrinkage analysis (free analysis) and 1-axial restrained
analysis (restrained analysis) to investigate the effect ofinternal and external restraint, respectively. In theanalysis, both specimens were divided by 1000 Voronoi

elements.
The analysis conditions are as follows. For boundary
condition of free analysis and restrained analysis, verticaldisplacement at the bottom of the specimen was fixed.
Humidity conditions of the atmosphere as boundary
conditions for mass transfer analysis were set uniformly
along both sides for free analysis and only along thebottleneck area for restrained analysis. The relative humidity of specimens was 100% at the initial stage. The
relative humidity of the atmosphere was constant at 80%
until 80 days and gradually decreased to 40% in the final
stage (total testing period of 200 days). Table 2 lists thematerial properties of A1 and A2. In the free analysis
model, A1 and A2 were used to investigate the effect ofinternal restraint by setting different values of moisture
conductivity, namely 10 and 30 mm2/day, respectively.
For the restrained analysis, the condition of A1 was applied and the restraint of the steel tendon was set to 100%
and 200% by the rigidity ratio shown in the following
equation. It should be noted that the moisture conductivity through cracks was assumed as the high value.

rigidity ratio =

ES AS
EC AC

(15)

Moisture boundaries

1000

where Es and As = Young’s modulus and cross-sectional
area of steel, respectively, Ec and Ac = Young’s modulus
and cross-sectional area of concrete at the middle sectionof the specimen, respectively.
Restrained
tendon

160

160

Free
analysis

Restrained
analysis

Unit: mmFig.13 Analysis model.

5.2. Free drying shrinkage analysis (free analysis)
Figure 14 shows the change in free drying shrinkage
strain for specimens A1 and A2 used for the free analysis.
The average strains rapidly occurred at the beginning due
to the different relative humidity values inside and outside the specimens. The rate of average strain thengradually decreased over time until 80 days. To accelerate the drying shrinkage strain, the relative humidity of
the atmosphere was set to decline after 80 days. There-

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

Shrinkage strain (µ)

40Relative humidityof atmosphere

80
60

A2

20

A1

40
20

0

Relative humidity (%)

154

A1
1.0

0.7

0.4

0

100

200

Age (days)

A2
1.0

Fig. 14 Free drying shrinkage strain.

Normal stress (MPa)20.7

40 days
80 days
160 days

1
0

0.4

40 days

-1

0

20
40
60
80
Length of specimen (mm)160 days

Relative moisture content

Specimen A1-2

80 days100A1a) Specimen A1

3.5 MPa

Normal stress (MPa)240 days
80 days
160 days10.0

0-0.35 MPa

-1Specimen A2
-2

0

A2
3.5 MPa

20
40
60
80
100
Length of specimen (mm)b) Specimen A2
0.0

Fig. 16 Profile of normal stress.

fore, normal shrinkage strains still occurred continuously.
Specimen A2, which had a higher moisture conductivity
value within bulk concrete obtained a higher averagestrain than specimen A1.
The relative humidity and normal stress distributions
that occurred inside the specimens are indicated in Fig.
15. Note that due to the difference in relative humidity
inside the specimens, tensile stresses occurred on thesurface and many micro cracks occurred at 160 days dueto the effect of internal restraint. For the normal stressprofile of both specimens at 40, 80, and 160 days, as
shown in Fig. 16, normal stress of specimen A1 washigher than that of A2 because specimen A1 has smaller
moisture conductivity and shows a larger gradient ofrelative humidity content. Therefore, low moisture conductivity generates higher internal restraint.

-0.35 MPa

40 days

80 daysNormal stress

160 days

160 days

Crack

Fig. 15 Conditions inside specimens.

5.3. Uniaxial restrained analysis (restrained
analysis)
The averaged strain that occurred in restrained tendons isshown in Fig. 17. In case of the high-restrained rate of
200%, through crack occurred at earlier age and theaveraged strain has a small value. When the throughcrack occurs, the strain of the tendon becomes almost
zero due to the release of restrained stress. Figure 18
shows the relative humidity distributions and crack pat-

155

80

100

Relative humidity
60 of atmosphere

80

40

60

a

100%

40

20
200%
0

50

100

20

b

c

150

200

0.0

Specimen A1

Relat ive humidity (%)

Shrinkage strain (µ)

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

Rebar
Specimen A2

Unit: mmFig. 19 Analytical model.
Table 3 Material properties.

1.0

Concrete
(MPa)

0.7

Rebar
(MPa)
0.4
c

（Voronoi1500） 200
1000

Fig.17 Average strain.

b（Voronoi1000）
Reinforcement ratio 2.0％

Rebar

Age (days)aReinforcement ratio 1.0％

a

b

c

Fig.18 Relative humidity and crack pattern.

terns before and after through crack occurrence for the
200% restrained specimen, with a, b, c corresponding tomarks a, b and c in Fig. 17.
Many cracks occurred on the surface, then cracks
propagated into the specimen and the crack widths
gradually increased. The specimen dehydrated according
to the distance between both sides and the crack plane.
Further, irregular relative humidity distribution on thesurface occurred due to the occurrence of many surface
cracks through which moisture was transferred.

6. Crack propagation Analysis due todrying shrinkage and corrosion
6.1. Analysis outline
Reinforced concrete bars as shown in Fig. 19 were
simulated by time-dependent structural analysis considering moisture transfer and chloride ion penetration withcorrosion expansion of rebars. The specimen size is 200
× 200 × 1000 mm. A single rebar was placed at the center,with a reinforcement ratio of 1.0% for specimen A1, andtwo rebars were placed near the specimen surfaces with areinforcement ratio of 2.0% for specimen A2. Specimen
A1 and A2 were analyzed in two-dimensions with 1000
and 1500 Voronoi elements, respectively. The material
properties are listed in Table 3.
For steel corrosion analysis (corrosion analysis) inwhich only chloride ion penetration is assumed, specimen A1 was analyzed. For moisture transfer and chloride
ion penetration (coupled analysis), specimens A1 and A2

Moisturetransfer(mm2/day)
Chlorideion
penetration(mm2/day)Tensile strength, ft
Compressive strength, fc
Yield strength, fy
Young's modulus, Es
Moisture diffusivity withinbulk concrete, Dw
Moisture diffusivity between cracks, Dcrw
Diffusion coefficient at water logged condition, Dcl0
Diffusion coefficient
between cracks, Dcrc

2.0
30.0
300.0
2.0 × 105
10.0
1.0 × 105
8.00
8.0 × 105

were analyzed. The threshold value was set to 1.2 × 10-6
g/mm3 at the initial corrosion. The diffusion coefficientsof chloride ions and moisture conductivity throughcracks were assumed to be high values. The diffusion
coefficient of chloride ions within bulk concrete was set
to 10 times the typical value in order to accelerate corrosion. For the boundary condition of structural analysis,
the displacement in the horizontal direction for one end
was fixed. On the other hand, only displacement in thevertical direction along specimen sides was assumed as
the boundary condition for mass transfer analysis. Relative humidity and chloride ion concentration of specimens at the initial stage were 1.0 and 0.0 g/mm3, respectively. The relative humidity of the atmosphere was
set to 0.8 during the period of 0 to 200 days and 0.6
afterward until 600 days. The chloride ion concentration
at the surface was set to 2.5 × 10-6 g/mm3 for all periods.
6.2. Steel corrosion analysis (corrosion analysis)
In this case, only chloride ion penetration was considered.
The distribution of stress in the lateral direction andcrack patterns of specimen A1 are shown in Fig. 20.
Cracks occurred along the rebar on both sides of the rebar
because chloride ions transfer into the specimen uniformly due to the absence of surface cracks. Tensilestress and micro-cracks occurred around the steel bar at
360 days and major cracks occurred at 400 days.

156

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

-3.0

0.0

3.0 Unit: MPa

360 daysStress distribution

360 days

Crack pattern

400 days

Crack pattern

50 days

250 days

600 daysa) Crack patterns
0.5

1.0

Fig. 20 Steel corrosion analysis.

6.3. Humidity and steel corrosion analysis(coupled analysis)
Specimen A1: Figure 21 shows cracking behaviors, thedistributions of relative humidity and the concentrationof chloride ions at 50, 250 and 600 days obtained fromthe developed analysis method. At 50 days, many micro-cracks occurred at the specimen surface due to theeffect of drying shrinkage. The soluble chloride ions
gradually penetrated into the specimen via the boundary
surface and the surface cracks. At 250 days, moreshrinkage cracks occurred and major cracks reached therebar, which led to diffusion of moisture and chloride
ions through the cracks, causing corrosion of the rebar.
At 600 days, mass transfer through cracks accelerated
because of the increase in crack widths.
Moreover the cracks propagated along the rebar, since
the specific concentration of chloride ions at the surfaceof rebar exceeded the threshold value. The analysis
simulated well structural damages due to mass transfer.
The crack patterns were different from the corrosion
analysis due to the effect of drying shrinkage, because
drying shrinkage can induce cracks on the surface andsoluble chloride ions can penetrate through existing
cracks. Thus, corrosion started locally near the cracks.
Specimen A2: In this case, proper covering and theeffect of chloride ions on the concrete covering wereconsidered, as shown in Fig. 22. Cracks occurred not
only at the specimen surface due to drying shrinkage but
also along the rebar due to the corrosion expansion of therebar. When the concentration of chloride ions along
both surfaces became higher than the threshold value,
spalling of the concrete covering occurred. As shown in
Fig. 22, spalling of concrete occurred at 262, 323 and374 days for 2, 2.5 and 3 cm covering, respectively.
Therefore, the thicker the covering, the longer the required slip off period and the longer the length.

7. Load application analysis
7.1 Analysis outline
The geometry of the beam specimen, which was modeledb) Relative humidity
Unit: 10-6g/mm3
0.0

2.5

c) Density of chloride ionsFig. 21 Combined analysis of drying shrinkage and corrosion.
262 days covering 20 mm

323 days

covering 25 mm

374 days

covering 30 mm

Fig. 22 Crack patterns of Specimen A2.

with 1500 Voronoi elements, is shown in Fig. 23. The
corrosion rate was set to the constant of 3.0 × 10-8g/mm2/day for all elements.

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

157

Loading plate C
L

2D10

250

2D29

500
Loading plate200

40 220 40

2D6

Unit : mm(a) Initial state at 400 days(b) No damage casea) RC beam(c) Drying caseb) Analytical model(d) Combination caseFig. 23 Geometry and dimensions of RC beam.

Beam specimens were analyzed considering the combined effect of shrinkage and steel corrosion before load
application (coupled-load analysis). The whole surfaceof the specimen except the surface of the loading plate
was assumed as the boundary condition of mass transfer
analysis. The relative humidity of the atmosphere was set
to 0.4 during the period of 0 to 200 days and 0.6 afterward until 400 days. The other conditions were the sameas described in section 6.1. It should be noted that thedeterioration of material properties was not consideredeven in the damaged specimens. The load application
analysis was performed for three cases: 1) Load application at 0 days for undamaged beam (undamaged), 2)
Load application at 400 days for damaged beams, with
only moisture transfer (drying) taken into consideration,and 3) Consideration of both moisture transfer and chloride ion penetration (combination).
7.2 Coupled-load analysis
Figure 24 shows the load displacement relation for threecases and Fig. 25 shows deformations and crack patterns
at the maximum load and at the initial state before
loading at 400 days. For damaged beams, the load

Fig. 25 Crack patterns and deformations.

gradually decreased during mass transfer, since the displacement of the loading plate was numerically fixed by
the displacement control. Different behaviors can beobserved in each of the three cases. The initial slope ofthe load displacement relationship for the damaged beam
becomes smaller than that of the undamaged beam. Thisis the effect of the initial cracks due to the mass transfer.
Figure 25 (a) shows the initial state at 400 days for the
combination case. Many cracks occurred on the surfaceof the specimen due to drying shrinkage and corrosion
expansion of the rebar. For the drying case, the maximum
load shows a larger value than that of the undamaged
case. This is the reason why the initial cracks obstructed
the propagation of diagonal cracks. This situation is often
observed in the test in which initial damage is introduced.
For the combination case, the load dropped suddenly
when the diagonal crack propagated. And, spalling of thecover concrete occurred at the same time. The deterioration of structural performance was simulated in the developed method.

8. Conclusions

Load (kN)

400
300
200
100
0
-50
0

Undamaged
DryingCombination
6
2
4
Displacement (mm)Fig.24 Load-displacement relationship.

This paper reports the development of a newly developed
unified analytical method with load application analysis
for short-term behavior and mass transfer analysis forlong-term behavior. The Rigid-Body-Spring Networks
model was used for structural analysis and the truss
networks model was utilized for mass transfer analysis.
As the mass transfer through bulk concrete and cracks,
moisture transfer and soluble chloride ion penetration
were adopted. Then, deterioration processes such as thecracking behavior under multi actions in consideration ofdrying shrinkage and rebar corrosion expansion were158

H. Nakamura, W. Srisoros R. Yashiro and M. Kunieda / Journal of Advanced Concrete Technology Vol. 4, No. 1, 147-158, 2006

simulated, taking into consideration the effect of masstransfer through cracks. The analysis was found tosimulate well structural damages as well as deterioration
process due to mass transfer. Moreover, a damaged
structure due to mass transfer was analyzed under external load, and the deterioration of structural performance was simulated with the developed method.
As the proposed model is under the first stage of development, the prime objectives of this paper are to propose and to verify the possibility of the new numericalmodel. Realizing the importance of creep and autogenous shrinkage as time-dependent behavior and degradation of material properties such as bond strength, theauthors will conduct further studies to consider theseissues in greater detail.
References
Bolander, J. E. and Berton, S. (2003). “Simulation ofshrinkage induced cracking in cement composite
overlays.” Cement and Concrete Composites, 26 (7),861-871.
Bolander, J. E. and Saito, S. (1998). “Fracture analysis
using spring networks with random geometry.”
Engineering Fracture Mechanic, 61. 569-591.
Comite Euro-International du Beton, CEB-FIP Model
Code 1990.
JCI. (2005). “Japan Concrete Institute Computer
Program of Thermal Stress Analysis for Massive
Concrete Structure (JCMAC).” Tokyo: Japan Concrete
Institute.
JSCE. (2002). “Standard specifications for concrete
structure-Structural
Performance
Verification.”
Tokyo: Japan Society of Civil Engineers.
Kawai, T. (1978). “New discrete models and theirapplication to seismic response analysis of structures.”

Nuclear Engineering and Design, 48. 207-229.
Kitoku, T., Higai, T and Saito, S. (2005). “The effect ofcover concrete and bar diameter on the bond stress-slip
relation of deformed bars.” Proceedings of the Japan
Concrete Institute, 27 (2). 757-762. (in Japanese)
Matsuo, T. et al. (1997). “Crack propagation analysis ofconcrete dependent on expansion of reinforcement.”
Proceedings of the Japan Concrete Institute, 19 (2).
99-104. (in Japanese)
Nakamura, H. and Higai, T. (2001). “Compressive
fracture energy and fracture zone length of concrete.”
In: Modeling of inelastic behavior of RC structures
under seismic loads. ASCE: 471-487.
Saeki, T. and Niki, H. (1996). “Migration of chloride ionsin non-saturated mortar.” Proceedings of the Japan
Concrete Institute, 18 (1). 969-974. (in Japanese)
Saito, S. and Bolander, J. E. (1997). “Numerical analysisof shrinkage cracking controlled by carbon fiber nets.”
Proceedings of 3rd International Symposium on
Non-Metallic (FRP) Reinforcement for Concrete
Structures, Sapporo, Japan, 2. 275-282.
Saito, S. and Hikosaka, H. (1999). “Numerical analysesof reinforced concrete structures using spring
networks models.” Journal of Materials, Concrete
Structures and Pavements, JSCE, 672 (44). 289-303.
Shima, H. and Okamura, H. (1987). “Bond-slip-strain
relationship of deformed bars embedded in massive
concrete.” Journal of Materials, Concrete Structures
and Pavements, JSCE, 378 (6). 165-174. (in Japanese)
Ueda, M., Takeuchi, H., Higuchi, H. and Kawai, T.
(1984). “A discrete limit analysis of reinforced
concrete structures.” In Computer-Aided Analysis andDesign of Concrete Structure, eds. Damjanic, F.,
Hinton, E., Owen, D. R. J., Bicanic, N. and Simovic,
V., Pineridge Press, Swansea, UK, 1369-1384.