9. Summary of results#
For all the models, a difference between solution SCILAB and the Code_Aster solution is obtained that is less than \(\text{1.5 \%}\), which makes it possible to validate the implementation of the various drying laws in the code. Let’s just note that we observe a violation of the maximum principle at the beginning of the simulation with Aster for Mensi’s law. This can be explained (by analogy with thermal) by the significant « shock » due to the way in which boundary conditions are imposed (imposed water concentration). This problem should be able to be solved by using lumped elements in the same way as in thermal.
Scilab command file
Main.sci:
getf (“/home/xxxx/library.sci”);
//PARAMETRES FROM SIMULATION NUMERIQUE
//
//width discretization
x0 = 0.08;
x = [-0.080:0.001:+0.080]; [n1 n2] = size (x);
//initial water content
Cinit = 128.8;
Ci = Cinit*ones (1, n2);
//limit conditions at 50% RH
CL = [58.8 58.8];.
Ci (1) = CL (1); Ci ($) = CL (2);
Ci_bazant = Ci;
//no time
dt = 60;//[s]
//coefficients of Bazant’s law
D1 = 3.0E-10;//[m2/s]
a = 0.04;
n = 6;
TMAX = 5;//years
//______________________________________________________________________
//
//SIMULATION NUMERIQUE
//
j = 0;
u=file (“open”, “resultat_g”, “unknown”);
for year = 0: TMAX,
Year
for day = 0:364,
For time = 0:23,
minute = 0;
for minute = 0:59,
d_Bazant = diffusion_bazant (D1, a, n, n, Cinit,58.8, Ci_bazant,293,293*ones (Ci) ,4700);
CI_bazant = linear_drying (d_Bazant, CI_bazant, CL, dt, x, « polar »);
if ((year == 0 & day == 0 & hour == 0 & hour == 1 & minute == 0) |…
(year == 0 & day == 3 & hour == 0 & minute == 0) |…
(year == 0 & day == 28 & hour == 0 & minute == 0) |…
(year == 1 & day == 91 & hour == 0 & minute == 0) |…
(year == 3 & day == 0 & hour == 0 & minute == 0) |…
(year == 5 & day == 0 & hour == 0 & hour == 0)) then,
Year, day, hour
t= 81:1:161;
for tk=t,
fprintf (u, “%6.3f %6.3f”, x (tk), ci_bazant (tk));
End,
end,//if
End,//for minute
End,//for time
end,//for day
End,//for year
file (“close”, u);
Library.sci:
//
//COEFFICIENT FROM DIFFUSION NON LINEAIRE POUR LES SECHAGE FROM BETON
//LOI FROM MENSI D (C) = a.exp (b.C)
//ACTIVATION THERMIQUE D (C, T) = D (C, T0). (T/T0) .exp [-Q/R* (1/T-1/T0)]
//
//The coefficient of Mensi’s law
//b Mensi’s law coefficient
//C vector of water contents [-]
//T0 reference temperature [K]
//T temperature vector [K]
//Q_R Q/A (worth 4700 K)
function D = diffusion_mensi (a, b, C, C, T0, T, Q_R),
D = a*ones (C).*exp (B*c);
D = D.*(T./ (T0*ones (T)));
D = D.*exp (Q_R* ((ones (T). /T0) - (ones (T). /T)));
end function,
//
//________________________________________________________________
//
//COEFFICIENT FROM DIFFUSION NON LINEAIRE POUR LES SECHAGE FROM BETON
//LOI FROM BAZANT
//ACTIVATION THERMIQUE D (C, T) = D (C, T0). (T/T0) .exp [-Q/R* (1/T-1/T0)]
//
//D1 coefficient of Bazant’s law
//The coefficient of Bazant’s law (alpha)
//n coefficient of Bazant’s law
//C0 water content at 100% RH
//Cext water content of the surrounding environment
//C vector of water contents [-]
//T0 reference temperature [K]
//T temperature vector [K]
//Q_R Q/R (worth 4700 K-1)
function D = diffusion_bazant (D1, a, n, n, C0, Cext, C, T0, T, Q_R),
h = ones (C) -0.5*((C-C0*ones (C))/(Cext-C0)) **2;
D = (((1-a) *ones (C). /(ones (C) + (4**n) *(ones (C) -h)**n)) +a*ones (C))*D1;
D = D.*(T./ (T0*ones (T)));
D = D.*exp (Q_R* ((ones (T). /T0) - (ones (T). /T)));
end function,
//________________________________________________________________
//
//DIFFUSION
//Resolution using the finite difference method
//
//D vector of the diffusion coefficients
//Ci vector of the water contents at time j [-]
//CL boundary conditions in xmin and xmax of the Dirichlet type (C=C0)
//and no time [s]
//x x-axis vector [m]
//mode_ polar/Cartesian
function Cf = linear_drying (D, Ci, CL, dt, x, mode_),
[n1,n2] = size (Ci);
dx_ = zeros (1, n2-2); dx_ (1: $) = (x (3: $) -x (1: $-2))) *0.5;
//Cf_ = (d*dt* (ones (x_). /(dx_**2)).* (Ci (3: $) -2*Ci (2: $-1) +Ci (1: $-2))) +Ci (2: $-1);
x3 = ((…
(x (2: $-1) -x (1: $-2)) .*…
(x (3: $) -x (1: $-2))…
) .*…
(x (3: $) -x (2: $-1))…
);
D2c_dx2 = 2*(Ci (3: $) .* (x (2: $-1) -x (1: $-2)))…
-Ci (2: $-1) .* (x (3: $) -x (1: $-2))…
+Ci (1: $-2) .* (x (3: $) -x (2: $-1))));
d2c_dx2 = d2c_dx2. /dx3;
if (mode_ == « polar ») then,
dc_Dx = (Ci (3: $) .*(x (2: $-1) -x (1: $-2))**2…
-Ci (1: $-2) .*(x (3: $) -x (2: $-1))**2);
//-Ci (2: $-1). *((x (2: $-1) -x (1: $-2))**2 - (x (3: $) -x (2: $-1))) **2)…
dd_Dx = (D (3: $) .*(x (2: $-1) -x (1: $-2))**2…
D (1: $-2) .*(x (3: $) -x (2: $-1))**2);
//- D (2: $-1). *((x (2: $-1) -x (1: $-2))**2 - (x (3: $) -x (2: $-1))) **2)…
DC_Dx = DC_Dx. /dx3;
dd_Dx = dd_Dx. /dx3;
i = find (x==0); [k1 k2] = size (i);
if (~ (k1==0)) then, x (i) = x (i+1) /10, end,
//printf (« 1st order %s; 2nd order %s », string (min (dc_Dx)), string (min (D2c_dx2))));
d2c_dx2 = d2c_dx2 + dc_Dx. /x (2: $-1);
End,
Cf_ = Ci (2: $-1) +dt*(D (2: $-1).*d2c_dx2);
if (mode_ == « polar ») then,
End,
Cf = zeros (1, n2); Cf (2: $-1) = Cf_; Cf (1) = CL (1); Cf ($) = CL (2);
end function,
//
//________________________________________________________________
//
Aster/Scilab comparison
A2.1 SECH_MENSI/SECH_GRANGER/SECH_NAPPE
A2.2 SECH_BAZANT
