1
%-----------------------------------------------------------------------
3
% This file is part of the Code_Saturne Kernel, element of the
4
% Code_Saturne CFD tool.
6
% Copyright (C) 1998-2008 EDF S.A., France
8
% contact: saturne-support@edf.fr
10
% The Code_Saturne Kernel is free software; you can redistribute it
11
% and/or modify it under the terms of the GNU General Public License
12
% as published by the Free Software Foundation; either version 2 of
13
% the License, or (at your option) any later version.
15
% The Code_Saturne Kernel is distributed in the hope that it will be
16
% useful, but WITHOUT ANY WARRANTY; without even the implied warranty
17
% of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18
% GNU General Public License for more details.
20
% You should have received a copy of the GNU General Public License
21
% along with the Code_Saturne Kernel; if not, write to the
22
% Free Software Foundation, Inc.,
23
% 51 Franklin St, Fifth Floor,
24
% Boston, MA 02110-1301 USA
26
%-----------------------------------------------------------------------
32
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
33
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
35
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
36
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
37
Ce sous-programme, appel� entre autre par \fort{preduv}, \fort{turbke}, \fort{covofi},
38
\fort{resrij}, \fort{reseps}, ..., r�sout les �quations de convection-diffusion
39
d'un scalaire $a$ avec termes sources du type :
40
\begin{equation}\label{Base_Codits_eq_ref}
42
\displaystyle f_s^{\,imp} (a^{n+1} - a^{n}) +
43
\theta \ \underbrace{\dive((\rho \underline{u})\,a^{n+1})}_{\text{convection implicite}}
44
-\theta \ \underbrace{\dive(\mu_{\,tot}\,\grad a^{n+1})}_{\text{diffusion implicite}}
46
= f_s^{\,exp}-(1-\theta) \ \underbrace{\dive((\rho \underline{u})\,a^{n})}_{\text{convection explicite}}
47
+ (1-\theta) \ \underbrace{\dive(\mu_{\,tot}\,\grad a^{n})}_{\text{diffusion explicite}}
50
o� $\rho \underline{u}$, $f_s^{exp}$ et $f_s^{imp}$ d�signent respectivement le flux de masse, les termes sources explicites et les termes lin�aris�s en $a^{n+1}$.
51
$a$ est un scalaire d�fini sur toutes les cellules\footnote{$a$, sous forme discr\`ete en espace, correspond \`a un vecteur dimensionn\'e \`a \var{NCELET} de composante $a_I$, I d\'ecrivant l'ensemble des cellules.}.
52
Par souci de clart� on suppose, en l'absence d'indication, les propri�tes
53
physiques $\Phi$ (viscosit� totale $\mu_{tot}$,...) et le flux de masse $(\rho
54
\underline{u})$ pris respectivement aux instants $n+\theta_\Phi$ et
55
$n+\theta_F$, o� $\theta_\Phi$ et $\theta_F$ d�pendent des sch�mas en temps
56
sp�cifiquement utilis�s pour ces grandeurs\footnote{cf. \fort{introd}}.
58
L'�criture des termes de convection et diffusion en maillage non orthogonal
59
engendre des difficult�s (termes de reconstruction et test de pente) qui sont
60
contourn�es en utilisant une m\'ethode it\'erative dont la limite, si elle
61
existe, est la solution de l'�quation pr�c�dente.
63
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
64
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
65
\section{Discr\'etisation}
66
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
67
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
68
Afin d'expliquer la proc�dure utilis�e pour traiter les difficult�s dues aux
69
termes de reconstruction et de test de pente dans les termes de
70
convection-diffusion, on note, de fa\c con analogue \`a ce qui est d\'efini dans
71
\fort{navsto} mais sans discr\'etisation spatiale associ\'ee, $\mathcal{E}_{n}$ l'op\'erateur :
72
\begin{equation}\label{Base_Codits_Eq_ref_small}
74
\mathcal{E}_{n}(a) = f_s^{\,imp}\,a + \theta\,\, \dive((\rho
75
\underline{u})\,a) - \theta\,\, \dive(\mu_{\,tot}\,\grad a)\\
76
- f_s^{\,exp} - f_s^{\,imp}\,a^{n} +(1-\theta)\,\,\dive((\rho
77
\underline{u})\, a^n) - (1-\theta)\,\, \dive(\mu_{\,tot}\,\grad a^n)
80
L'\'equation (\ref{Base_Codits_eq_ref}) s'\'ecrit donc :
82
\mathcal{E}_{n}(a^{n+1}) = 0
84
La quantit\'e $\mathcal{E}_{n}(a^{n+1})$ comprend donc :\\
85
\hspace*{1.cm} $\rightsquigarrow$ $f_s^{\,imp}\,a^{n+1}$, contribution des
86
termes diff\'erentiels d'ordre $0$ lin\'eaire en $a^{n+1}$,\\
87
\hspace*{1.cm} $\rightsquigarrow$ $\theta\,\,
88
\dive((\rho\underline{u})\,a^{n+1})
89
- \theta\,\, \dive(\mu_{\,tot}\,\grad a^{n+1})$, termes de convection-diffusion
90
implicites complets (termes non reconstruits + termes de reconstruction)
91
lin\'eaires\footnote{Lors de la discr�tisation en espace, le caract�re lin�aire
92
de ces termes pourra cependant �tre perdu, notamment � cause du test de pente.}
94
\hspace*{1.cm} $\rightsquigarrow$ $f_s^{\,exp}- f_s^{\,imp}\,a^n$ et
95
$(1-\theta)\,\,\dive((\rho
96
\underline{u})\,a^n) - (1-\theta)\,\, \dive(\mu_{\,tot}\,\grad a^n)$ l'ensemble
97
des termes explicites (y compris la partie explicite provenant du sch\'ema en
98
temps appliqu\'e \`a la convection diffusion).\\\\
100
De m\^eme, on introduit un op\'erateur $\mathcal{EM}_{n}$ approch\'e de
101
$\mathcal{E}_{n}$, lin\'eaire et simplement inversible, tel que son
102
expression contient :\\
103
\hspace*{1.cm}$\rightsquigarrow$ la prise en compte des termes lin\'eaires en $a$,\\
104
\hspace*{1.cm}$\rightsquigarrow$ la convection int�gr�e par un sch�ma d�centr� amont
105
(upwind) du premier ordre en espace,\\
106
\hspace*{1.cm}$\rightsquigarrow$ les flux diffusifs non reconstruits.\\
108
\mathcal{EM}_{n}(a) = f_s^{\,imp}\,a + \theta\,\,[\dive((\rho
109
\underline{u})\,a)]^{\textit{amont}} - \theta\,\, [\dive(\mu_{\,tot}\,\grad a)]^{\textit{N Rec}}
111
Cet op\'erateur permet donc de contourner la difficult\'e induite par la pr\'esence d'\'eventuelles non lin\'earit\'es introduites par l'activation du test de pente lors du sch\'ema convectif, et par le remplissage important de la structure de la matrice d\'ecoulant de la pr\'esence des gradients propres \`a la reconstruction.\\
112
On a la relation\footnote{On pourra se reporter au sous-programme
113
\fort{matrix} pour plus de d\'etails relativement \`a
114
$\mathcal{EM_{\it{disc}}}$, op\'erateur discret agissant sur un scalaire $a$.}, pour toute cellule $\Omega_I$ de centre $I$ :
115
\begin{equation}\notag
116
\mathcal{EM_{\it{disc}}}(a,I) = \int_{\Omega_i}\mathcal{EM}_{n}(a) \, d\Omega
118
On cherche � r�soudre :
120
0 =\mathcal{E}_{n}(a^{n+1}) = \mathcal{EM}_{n}(a^{n+1}) + \mathcal{E}_{n}(a^{n+1}) - \mathcal{EM}_{n}(a^{n+1})
124
\mathcal{EM}_{n}(a^{n+1}) = \mathcal{EM}_{n}(a^{n+1}) - \mathcal{E}_{n}(a^{n+1})
126
On va pour cela utiliser un algorithme de type point fixe en d�finissant la
127
suite $(a^{n+1,\,k})_{k\in \grandN}$\footnote{Dans le cas ou le point fixe en
128
vitesse-pression est utilis� (\var{NTERUP}$>$ 1) $a^{n+1,0}$ est initialis� par
129
la derni�re valeur obtenue de $a^{n+1}$.}:
130
\begin{equation}\notag
131
\left\{\begin{array}{l}
132
a^{n+1,\,0} = a^{n}\\
133
a^{n+1,\,k+1} = a^{n+1,\,k} + \delta a^{n+1,\,k+1}
136
o� $\delta a^{n+1,\,k+1}$ est solution de :
138
\mathcal{EM}_{n}(a^{n+1,\,k} + \delta a^{n+1,\,k+1}) = \mathcal{EM}_{n}(a^{n+1,\,k}) - \mathcal{E}_{n}(a^{n+1,\,k})
140
Soit encore, par lin�arit� de $\mathcal{EM}_{n}$ :
142
\mathcal{EM}_{n}(\delta a^{n+1,\,k+1}) = - \mathcal{E}_{n}(a^{n+1,\,k})
143
\label{Base_Codits_Eq_Codits}
146
Cette suite, coupl\'ee avec le choix de l'op\'erateur $\mathcal{E}_{n}$, permet donc de lever la difficult� induite par la
147
pr\'esence de la convection (discr\'etis\'ee \`a l'aide de sch\'emas num\'eriques
148
qui peuvent introduire des non lin\'earit\'es) et les termes de
149
reconstruction. Le sch\'ema r\'eellement choisi par l'utilisateur pour la
150
convection (donc \'eventuellement non lin\'eaire si le test de pente est activ\'e) ainsi que les termes de
151
reconstruction vont \^etre pris \`a l'it�ration $k$ et trait\'es au second membre {\it via} le sous-programme \fort{bilsc2}, alors que les termes
152
non reconstruits sont pris � l'it�ration $k+1$ et repr\'esentent donc les
153
inconnues du syst\`eme lin\'eaire r\'esolu par \fort{codits}\footnote{cf. le sous-programme
156
On suppose de plus que cette suite $(a^{n+1,\,k})_k$ converge vers la solution
157
$a^{n+1}$ de l'\'equation (\ref{Base_Codits_Eq_ref_small}), {\it i.e.}
158
$\lim\limits_{k\rightarrow\infty} \delta a^{n+1,\,k}\,=\,0$, ceci pour tout $n$ donn\'e.\\
159
(\ref{Base_Codits_Eq_Codits}) correspond \`a l'\'equation r\'esolue par \fort{codits}. La
160
matrice $\tens{EM}_{\,n}$, matrice associ\'ee \`a $\mathcal{EM}_{n}$ est
161
\`a inverser, les termes non lin\'eaires sont mis au second membre mais sous forme
162
explicite (indice $k$ de $a^{n+1,\,k}$) et ne posent donc plus de probl\`eme.
164
\minititre{Remarque 1}
165
La viscosit� $\mu_{\,tot}$ prise dans $\mathcal{EM}_{n}$ et dans
166
$\mathcal{E}_{n}$ d�pend du mod�le de turbulence utilis�. Ainsi on a
167
$\mu_{\,tot}=\mu_{\,laminaire} + \mu_{\,turbulent}$
168
dans $\mathcal{EM}_{n}$ et dans $\mathcal{E}_{n}$ sauf lorsque l'on
169
utilise un mod�le $R_{ij}-\varepsilon$, auquel cas on a
170
$\mu_{\,tot}=\mu_{\,laminaire}$.\\
171
Le choix de $\mathcal{EM}_{n}$ �tant {\it a
172
priori} arbitraire ($\mathcal{EM}_{n}$ doit \^etre lin\'eaire et la suite
173
$(a^{n+1,\,k})_{k\in\grandN}$ doit converger pour tout $n$ donn\'e), une option des mod�les
174
$R_{ij}-\varepsilon$ ($\var{IRIJNU}=1$) consiste \`a for\c cer $\mu_{\,tot}^n$
175
dans l'expression de $\mathcal{EM}_{n}$ \`a la
176
valeur $\mu_{\,laminaire}^n + \mu_{\,turbulent}^n$ lors de l'appel \`a
177
\fort{codits} dans le sous-programme \fort{navsto}, pour l'\'etape de
178
pr\'ediction de la vitesse. Ceci n'a pas de sens
179
physique (seul $\mu_{\,laminaire}^n$ �tant cens� intervenir), mais cela peut
180
dans certains cas avoir un effet stabilisateur, sans que cela modifie pour
181
autant les valeurs de la limite de la suite $(a^{n+1,\,k})_k$.\\
183
\minititre{Remarque 2}
184
Quand \fort{codits} est utilis\'e pour le couplage instationnaire renforc�
185
vitesse-pression (\var{IPUCOU}=1), on fait une seule it�ration $k$ en initialisant la suite $(a^{n+1,\,k})_{k\in\grandN}$ � z�ro. Les conditions de type Dirichlet sont
186
annul�es (on a $\var{INC}\,=\,0$) et le second membre est �gal � $\rho |\Omega_i|$.
187
Ce qui permet d'obtenir une approximation de type diagonal de
189
n\'ecessaire lors de l'�tape de correction de la vitesse\footnote{cf. le sous-programme \fort{resolp}.}.
191
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
192
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
193
\section{Mise en \oe uvre}
194
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
195
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
196
L'algorithme de ce sous-programme est le suivant :\\
197
- d�termination des propri�t�s de la matrice $\tens{EM}_{n}$ (sym\'etrique
198
si pas de convection, non sym\'etrique sinon)\\
199
- choix automatique de la m�thode de r�solution pour l'inverser si l'utilisateur
200
ne l'a pas sp\'ecifi\'e pour la variable trait\'ee. La m\'ethode de Jacobi est utilis\'ee par d\'efaut pour toute variable scalaire $a$ convect\'ee. Les m\'ethodes
201
disponibles sont la m\'ethode du gradient conjugu\'e, celle de Jacobi, et le
202
bi-gradient conjugu\'e stabilis\'e ($BICGStab$) pour les matrices non
203
sym\'etriques. Un pr\'econditionnement diagonal est possible et utilis\'e par
204
d\'efaut pour tous ces solveurs except� Jacobi.\\
205
- prise en compte de la p\'eriodicit\'e (translation ou rotation d'un scalaire, vecteur ou tenseur),\\
206
- construction de la matrice $\tens{EM}_{n}$ correspondant � l'op\'erateur
207
lin\'eaire $\mathcal{EM}_{n}$ par appel au sous-programme
208
\fort{matrix}\footnote{ On rappelle que dans \fort{matrix}, la convection est
209
trait�e, quelque soit le choix de l'utilisateur, avec un sch�ma d�centr� amont d'ordre 1 en espace et qu'il n'y a pas de reconstruction
210
pour la diffusion. Le choix de l'utilisateur quant au
211
sch�ma num�rique pour la convection intervient uniquement lors de l'int�gration
212
des termes de convection de $\mathcal{E}_{n}$, au second membre de
213
(\ref{Base_Codits_Eq_Codits}) dans le sous-programme \fort{bilsc2}.}. Les termes implicites correspondant �
214
la partie diagonale de la matrice et donc aux contributions diff\'erentielles
215
d'ordre $0$ lin\'eaires en $a^{n+1}$,({\it i.e} $f_s^{imp}$), sont stock�s dans le tableau \var{ROVSDT} (r\'ealis\'e en amont du sous-programme appelant \fort{codits}).\\
216
- cr�ation de la hi�rarchie de maillage si on utilise le multigrille
217
($ \var{IMGRP}\,>0 $).\\
218
- appel � \fort{bilsc2} pour une �ventuelle prise en compte de la
219
convection-diffusion explicite lorsque $\theta \ne 0$.\\
220
- boucle sur le nombre d'it�rations de 1 \`a $\var{NSWRSM}$ (appel� $\var{NSWRSP}$ dans \fort{codits}).
221
Les it�rations sont repr�sent�es par $k$ appel�
222
\var{ISWEEP} dans le code et d\'efinissent les indices de la suite $(a^{n+1,\,k})_k$
223
et de $(\delta a^{n+1,\,k})_k$.\\
224
Le second membre est scind\'e en deux parties :\\
225
\hspace*{1cm}{\tiny$\blacksquare$}\ un terme, affine
226
en $a^{n+1,\,k-1}$, facile \`a mettre \`a jour dans le cadre de la r\'esolution
227
par incr\'ement, et qui s'\'ecrit :
228
\begin{equation}\notag
229
-f_s^{\,imp} \left(\,a^{n+1,\,k-1} - a^{n+1,0}\right) + f_s^{\,exp}- (1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^{n+1,0}) - \dive(\mu_{\,tot}\,\grad a^{n+1,0})\,\right]\\
232
\hspace*{1cm}{\tiny$\blacksquare$}\ les termes issus de la
233
convection/diffusion (avec reconstruction) calcul�e par \fort{bilsc2}.\\
234
\begin{equation}\notag
235
- \theta\,\left[\,\dive\left((\rho \underline{u})\,a^{n+1,\,k-1}\right)- \dive\left(\mu_{\,tot}\,\grad a^{n+1,\,k-1}\right)\right]
238
La boucle en $k$ est alors la suivante :
240
\item Calcul du second membre, hors contribution des termes de
241
convection-diffusion explicite $\var{SMBINI}$; le second membre complet correspondant
242
\`a $\mathcal{E}_{n}(a^{n+1,\,k-1})$ est, quant \`a lui, stock\'e dans le
243
tableau $\var{SMBRP}$, initialis\'e par $\var{SMBINI}$ et compl\'et\'e par les
244
termes reconstruits de convection-diffusion par appel au sous-programme
246
\`A l'it\'eration $k$, $\var{SMBINI}$ not\'e $\var{SMBINI}^{\,k}$ vaut :\\
247
\begin{equation}\notag
248
\var{SMBINI}^{\,k}\ = f_s^{\,exp}-(1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^n) - \dive(\mu_{\,tot}\,\grad a^n)\,\right]-\,f_s^{\,imp}(\,a^{n+1,\,k-1} - a^n\,) \\
251
$\bullet$ Avant de boucler sur $k$,un premier appel au sous-programme \fort{bilsc2} avec $\var{THETAP}=1-\theta$ permet de prendre
252
en compte la partie explicite des termes de convection-diffusion provenant du sch\'ema en temps.
253
\begin{equation}\notag
255
\var{SMBRP}^{\,0} = f_s^{\,exp} -(1-\theta)\,[\,\dive((\rho \underline{u})\,a^n) - \dive(\mu_{\,tot}\,\grad a^n)\,]\\
257
Avant de boucler sur $k$, le second membre $\var{SMBRP}^{\,0}$ est stock\'e dans le tableau $\var{SMBINI}^{\,0}$ et sert pour l'initialisation du reste du calcul.
258
\begin{equation}\notag
259
\var{SMBINI}^{\,0} =\var{SMBRP}^{\,0}
262
$\bullet$ pour $k = 1$,
263
\begin{equation}\notag
265
\var{SMBINI}^{\,1}\ &=f_s^{\,exp}-(1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^n) - \dive(\mu_{\,tot}\,\grad a^n)\,\right]-\,f_s^{\,imp}\,(\,a^{n+1,\,0} - a^n\,)\\
266
&=f_s^{\,exp}- (1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^{n+1,\,0}) - \dive(\mu_{\,tot}\,\grad a^{n+1,\,0})\,\right]-f_s^{\,imp}\,\delta a^{n+1,\,0} \\
270
\begin{equation}\notag
271
\var{SMBINI}^{\,1}\ =\ \var{SMBINI}^{\,0} - \var{ROVSDT}\,*(\,\var{PVAR}-\,\var{PVARA})
273
et $\var{SMBRP}^{\,1}$ est compl\'et\'e par un second appel au sous-programme \fort{bilsc2} avec $\var{THETAP}=\theta$, de mani\`ere \`a ajouter dans le second membre la partie de la convection-diffusion implicite.
274
\begin{equation}\notag
276
\var{SMBRP}^{\,1} & = \var{SMBINI}^{\,1} -\theta\,\left[\,\dive((\rho \underline{u})\,a^{n+1,\,0}) - \dive(\mu_{\,tot}\,\grad a^{n+1,\,0})\,\right]\\
277
& = f_s^{\,exp}\ - (1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^{n}) - \dive(\mu_{\,tot}\,\grad a^{n})\,\right]- f_s^{\,imp}\,(a^{n+1,\,0} -a^{n}) \\
278
& -\theta\,\left[\,\dive((\rho \underline{u})\,a^{n+1,\,0}) - \dive(\mu_{\,tot}\,\grad a^{n+1,\,0})\,\right]\\
281
$\bullet$ pour $k = 2$,\\
282
de fa\c con analogue, on obtient :
283
\begin{equation}\notag
285
\var{SMBINI}^{\,2}\ &=f_s^{\,exp}-(1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^n) - \dive(\mu_{\,tot}\,\grad a^n)\,\right]-\,f_s^{\,imp}\,(\,a^{n+1,\,1} - a^n\,)\\
289
\begin{equation}\notag
290
\var{SMBINI}^{\,2}\ =\ \var{SMBINI}^{\,1} - \var{ROVSDT}\,*\,\var{DPVAR}^{\,1}
292
l'appel au sous-programme \fort{bilsc2}, \'etant syst\'ematiquement fait par la suite avec $\var{THETAP}=\theta$, on obtient de m\^eme :
293
\begin{equation}\notag
295
\var{SMBRP}^{\,2}\ &=\ \var{SMBINI}^{\,2}-\theta\left[\dive\left((\rho \underline{u})\,a^{n+1,\,1}\right)- \dive\left(\mu_{\,tot}\,\grad \,a^{n+1,\,1}\right)\right]\\
299
\begin{equation}\notag
300
a^{n+1,\,1}=\var{PVAR}^{\,1}=\var{PVAR}^{\,0}+\var{DPVAR}^{\,1}=a^{n+1,\,0}+\delta a^{n+1,\,1}
302
$\bullet$ pour l'it\'eration $k+1$,\\
303
Le tableau $\var{SMBINI}^{\,k+1}$ initialise le second membre complet
304
$\var{SMBRP}^{\,k+1}$ auquel vont \^etre rajout\'ees les contributions
305
convectives et diffusives {\it via} le sous-programme \fort{bilsc2}.\\
307
\begin{equation}\notag
309
\var{SMBINI}^{\,k+1}\ &= \var{SMBINI}^{\,k} - \var{ROVSDT}\,*\,\var{DPVAR}^{\,k}\\
312
Puis suit le calcul et l'ajout des termes de convection-diffusion reconstruits de $-\ \mathcal{E}_{n}(a^{n+1,\,k})$, par appel au sous-programme
313
\fort{bilsc2}. On rappelle que la convection est prise en compte � cette �tape
314
par le sch�ma num�rique choisi par l'utilisateur (sch�ma d�centr� amont du
315
premier ordre en espace, sch�ma centr� du second ordre en espace, sch�ma
316
d�centr� amont du second ordre S.O.L.U. ou une
317
pond�ration (blending) des sch�mas dits du second ordre (centr� ou S.O.L.U.) avec le sch�ma
318
amont du premier ordre, utilisation \'eventuelle d'un test de pente).\\
319
Cette contribution (convection-diffusion) est alors
320
ajout�e dans le second membre $\var{SMBRP}^{\,k+1}$ (initialis\'e par $\var{SMBINI}^{\,k+1}$).
321
\begin{equation}\notag
323
\var{SMBRP}^{\,k+1}\ &= \var{SMBINI}^{\,k+1} - \theta\,\left[\,\dive\left((\rho \underline{u})\,a^{n+1,\,k}\right)- \dive\left(\mu_{\,tot}\,\grad a^{n+1,\,k}\right)\right]\\
324
& = f_s^{\,exp}-(1-\theta)\,\left[\,\dive((\rho \underline{u})\,a^n) - \dive(\mu_{\,tot}\,\grad a^n)\,\right]- f_s^{\,imp}\,(a^{n+1,\,k} -a^{n}) \\
325
&-\theta\,\left[\,\dive((\rho \underline{u})\,a^{n+1,k}) - \dive(\mu_{\,tot}\,\grad a^{n+1,k})\,\right]\\
329
\item R\'esolution du syst�me lin�aire en $\delta a^{n+1,\,k+1}$ correspondant
330
\`a l'\'equation (\ref{Base_Codits_Eq_Codits}) par inversion de la matrice
331
$\tens{EM}_{n}$, en appelant le sous programme \fort{invers}.
332
On calcule $a^{n+1,\,k+1}$ gr\^ace \`a la formule :
333
\begin{equation}\notag
334
a^{n+1,\,k+1} = a^{n+1,\,k} + \delta a^{n+1,\,k+1}
337
\begin{equation}\notag
338
\var{PVAR}^{\,k+1} = \var{PVAR}^{\,k} + \var{DPVAR}^{\,k+1}
341
\item Traitement de la p\'eriodicit\'e et du parall�lisme.
342
\item Test de convergence :\\
343
Il porte sur la quantit\'e $||\var{SMBRP}^{\,k+1}|| < \varepsilon
344
||\tens{EM}_{n}(a^{n}) + \var{SMBRP}^{\,1}|| $, o\`u $||\,.\,||$ repr\'esente la
346
Si le test est v\'erifi\'e, la convergence est atteinte et on sort de la
347
boucle sur les it�rations. La solution recherch�e est $a^{\,n+1} = a^{n+1,\,k+1}$.\\
348
Sinon, on continue d'it\'erer dans la limite des it�rations impos�es par $\var{NSWRSM}$ dans \fort{usini1}.\\
349
En ``continu'' ce test de convergence s'�crit aussi :
350
\begin{equation}\notag
352
||\var{SMBRP}^{\,k+1}||& < \varepsilon ||f_s^{\,exp}\ - \dive((\rho \underline{u})\,a^{n}) + \dive(\mu_{\,tot}\,\grad a^{n}) \\
353
& +[\dive((\rho \underline{u})\,a^{n})]^{\textit{amont}} + [\dive(\mu_{\,tot}\,\grad a^{n})]^{\textit{N Rec}}||\\
356
Si bien que sur maillage orthogonal avec sch�ma de convection upwind et en l'absence de terme source, la suite converge en th�orie en une unique it�ration puisque par construction~:
357
\begin{equation}\notag
359
||\var{SMBRP}^{\,2}||=\,0\,& < \varepsilon ||f_s^{\,exp}||
366
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
367
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
368
\section{Points \`a traiter}\label{Base_Codits_section4}
369
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
370
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
371
\etape{Approximation $\mathcal{EM}_{n}$ de l'op\'erateur
373
D'autres approches visant soit \`a modifier la d\'efinition de l'approxim\'ee, prise en compte du sch\'ema centr\'e sans reconstruction par exemple,
374
soit \`a abandonner cette voie seraient \`a \'etudier.\\
376
\etape{Test de convergence}
377
La quantit\'e d\'efinissant le test de convergence est �galement � revoir, �ventuellement � simplifier.
379
\etape{Prise en compte de $T_s^{imp}$}
380
Lors de la r�solution de l'�quation par \fort{codits}, le tableau \var{ROVSDT} a
381
deux fonctions : il sert � calculer la diagonale de la matrice (par appel de
382
\fort{matrix}) et il sert � mettre � jour le second membre � chaque
383
sous-it�ration de la r�solution en incr�ments. Or, dans le cas o� $T_s^{imp}$
384
est positif, on ne l'int�gre pas dans \var{ROVSDT}, afin de ne pas affaiblir la
385
diagonale de la matrice. De ce fait, on ne l'utilise pas pour mettre � jour le
386
second membre, alors que ce serait tout � fait possible. Au final, on obtient
387
donc un terme source trait� totalement en explicite ($T_s^{exp}+T_s^{imp}a^n$),
388
alors que la r�solution en incr�ments nous permettrait justement de l'impliciter
389
quasiment totalement ($T_s^{exp}+T_s^{imp}a^{n+1,k_{fin}-1}$, o� $k_{fin}$ est
390
la derni�re sous-it�ration effectu�e).\\
391
Pour ce faire, il faudrait d�finir deux tableaux \var{ROVSDT} dans \fort{codits}.
added
removed
doc/theory/codits.tex
