1 : C YLAP2C SOURCE MAGN 09/07/30 21:22:31 6458
2 : SUBROUTINE YLAP2C(ICOGRV,MPROC,MPVITC,
3 : $ MPVOLU,MPNORM,MPSURF,MELEFL,
4 : $ KRFACE,KRCENT,LCLIMV,KRVIMP,LCLITO,KRTOIM,
5 : $ NOMINC,
6 : $ MU,
7 : $ IJACO,
8 : $ IMPR,IRET)
9 : IMPLICIT INTEGER(I-N)
10 : IMPLICIT REAL*8 (A-H,O-Z)
11 : C***********************************************************************
12 : C NOM : YLAP2C
13 : C DESCRIPTION : Calcul de la matrice jacobienne du résidu du laplacien
14 : C VF 3D.
15 : C Ici, on ne calcule que les contributions 'simples'
16 : C à la matrice jacobienne faisant intervenir les
17 : C coefficients pour le calcul des gradients de vitesse
18 : C (ICOGRV).
19 : C (contributions à (d Res_{\rho u} / d var) et (d
20 : C Res_{\rho v} / d var)
21 : C var prenant successivement les valeurs :
22 : C \rho, \rho u, \rho v, \rho e_t )
23 : C
24 : C
25 : C NOTE PERSO : On pourrait programmer ça plus lisiblement en stockant
26 : C les contributions dans un tableau de pointeurs (2
27 : C indices, c'est possible ?) et en effectuant des produits
28 : C matriciels (coeff. x matrices de dérivées).
29 : C On n'y coupera pas si on veut regrouper 2D et 3D...
30 : C On ne va pas le faire.
31 : C
32 : C
33 : C LANGAGE : ESOPE
34 : C AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
35 : C mél : gounand@semt2.smts.cea.fr
36 : C***********************************************************************
37 : C APPELES (UTIL) : AJMTK : ajoute un objet de type MATSIM (non
38 : C standard) à un objet de type MATRIK.
39 : C APPELE PAR : YLAP2A : Calcul de la matrice jacobienne du
40 : C résidu du laplacien VF 3D.
41 : C***********************************************************************
42 : C ENTREES : ICOGRV (type MCHELM) : coefficients pour le
43 : C calcul du gradient de la vitesse aux
44 : C interfaces.
45 : C MPROC (type MPOVAL) : masse volumique par
46 : C élément.
47 : C MPVITC (type MPOVAL) : vitesse par élément.
48 : C MPVOLU (type MPOVAL) : volume des éléments.
49 : C MPNORM (type MPOVAL) : normale aux faces.
50 : C MPSURF (type MPOVAL) : surface des faces.
51 : C MELEFL (type MELEME) : connectivités face-(centre
52 : C gauche, centre droit).
53 : C KRFACE (type MLENTI) : tableau de repérage dans
54 : C le maillage des faces des éléments.
55 : C KRCENT (type MLENTI) : tableau de repérage dans
56 : C le maillage des centres des éléments.
57 : C LCLIMV (type logique) : .TRUE. => CL de Dirichlet
58 : C sur la vitesse.
59 : C KRVIMP (type MLENTI) : tableau de repérage dans
60 : C maillage des CL de Dirichlet sur la
61 : C vitesse.
62 : C LCLITO (type logique) : .TRUE. => CL de Dirichlet
63 : C sur le tenseur des contraintes.
64 : C KRTOIM (type MLENTI) : tableau de repérage dans
65 : C maillage des CL de Dirichlet sur le tenseur des
66 : C contraintes
67 : C NOMINC (type MLMOTS) : noms des inconnues.
68 : C MU (type réel) : viscosité dynamique (SI).
69 : C ENTREES/SORTIES : IJACO (type MATRIK) : matrice jacobienne du
70 : C résidu du laplacien VF 3D.
71 : C SORTIES : -
72 : C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
73 : C***********************************************************************
74 : C VERSION : v1, 28/08/2001, version initiale
75 : C HISTORIQUE : v1, 28/08/2001, création
76 : C HISTORIQUE :
77 : C HISTORIQUE :
78 : C***********************************************************************
79 : C Prière de PRENDRE LE TEMPS de compléter les commentaires
80 : C en cas de modification de ce sous-programme afin de faciliter
81 : C la maintenance !
82 : C***********************************************************************
83 : -INC CCOPTIO
84 : -INC SMCHPOI
85 : POINTEUR MPROC.MPOVAL ,MPVITC.MPOVAL
86 : POINTEUR MPSURF.MPOVAL,MPNORM.MPOVAL,MPVOLU.MPOVAL
87 : -INC SMCHAML
88 : POINTEUR ICOGRV.MCHELM,JCOGRV.MCHAML
89 : POINTEUR KDUNDX.MELVAL,KDUNDY.MELVAL,KDUNDZ.MELVAL
90 : -INC SMELEME
91 : POINTEUR MELEFL.MELEME
92 : POINTEUR MCOGRV.MELEME
93 : -INC SMLENTI
94 : POINTEUR KRVIMP.MLENTI,KRTOIM.MLENTI
95 : POINTEUR KRCENT.MLENTI,KRFACE.MLENTI
96 : -INC SMLMOTS
97 : POINTEUR NOMINC.MLMOTS
98 : -INC SMMATRIK
99 : POINTEUR IJACO.MATRIK
100 : *
101 : * Objet matrice élémentaire simplifié
102 : *
103 : SEGMENT GMATSI
104 : INTEGER POIPR1(NPP1,NEL1)
105 : INTEGER POIDU1(1,NEL1)
106 : INTEGER POIPR2(NPP2,NEL2)
107 : INTEGER POIDU2(2,NEL2)
108 : POINTEUR LMATSI(0).MATSIM
109 : ENDSEGMENT
110 : * Contributions simples de la part du gradient de
111 : * vitesse (CTSGRV)
112 : POINTEUR CTSGRV.GMATSI
113 : SEGMENT MATSIM
114 : CHARACTER*8 NOMPRI,NOMDUA
115 : REAL*8 VALMA1(1,NPP1,NEL1)
116 : REAL*8 VALMA2(2,NPP2,NEL2)
117 : ENDSEGMENT
118 : POINTEUR ROURHO.MATSIM
119 : POINTEUR ROVRHO.MATSIM
120 : POINTEUR ROWRHO.MATSIM
121 : POINTEUR ROUROU.MATSIM
122 : POINTEUR ROVROU.MATSIM
123 : POINTEUR ROWROU.MATSIM
124 : POINTEUR ROUROV.MATSIM
125 : POINTEUR ROVROV.MATSIM
126 : POINTEUR ROWROV.MATSIM
127 : POINTEUR ROUROW.MATSIM
128 : POINTEUR ROVROW.MATSIM
129 : POINTEUR ROWROW.MATSIM
130 : *
131 : REAL*8 MU
132 : *
133 : INTEGER IMPR,IRET
134 : *
135 : LOGICAL LCLIMV,LCLITO
136 : LOGICAL LMUR
137 : LOGICAL LCTRB1
138 : *
139 : INTEGER IELEM,IPD,IPP,ISOUCH,IEL1,IEL2
140 : INTEGER NELEM,NPD,NPP,NSOUCH,NEL1,NEL2,NPP1,NPP2
141 : INTEGER NGCDRO,NGCGAU,NGFACE,NPPRIM,NPDUAL
142 : INTEGER NLCENP,NLCEND,NLFACE,NLCLTO,NLCLV
143 : INTEGER NPTEL
144 : *
145 : REAL*8 ALPHAX,ALPHAY,ALPHAZ,CNX,CNY,CNZ
146 : REAL*8 SIGNOR,SURFFA,VOLUEL
147 : REAL*8 RHOP,UP,VP,WP
148 : REAL*8 FACTOR
149 : REAL*8 DUDRHO,DUDROU
150 : REAL*8 DVDRHO,DVDROV
151 : REAL*8 DWDRHO,DWDROW
152 : REAL*8 BROUDU,BROUDV,BROUDW
153 : REAL*8 BROVDU,BROVDV,BROVDW
154 : REAL*8 BROWDU,BROWDV,BROWDW
155 : C
156 : *
157 : * Executable statements
158 : *
159 : IF (IMPR.GT.2) WRITE(IOIMP,*) 'Entrée dans ylap2c.eso'
160 : * On calcule les contributions à (d Res_{\rho u} / d var), (d
161 : * Res_{\rho v} / d var) et (d Res_{\rho w} / d var) ; var prenant
162 : * successivement les valeurs :
163 : * \rho, \rho u, \rho v, \rho w, \rho e_t
164 : * Note :
165 : * - (d Res_{\rho u} / d \rho e_t) = 0
166 : * - (d Res_{\rho v} / d \rho e_t) = 0
167 : * - (d Res_{\rho w} / d \rho e_t) = 0
168 : * On dérive les termes : \tens{\tau} \prod \vect{n}
169 : * Les noms de matrices élémentaires (type MATSIM) associées sont :
170 : * ROURHO, ROUROU, ROUROV, ROUROW,
171 : * ROVRHO, ROVROU, ROVROV, ROVROW,
172 : * ROWRHO, ROWROU, ROWROV, ROWROW.
173 : IF (LCLIMV) THEN
174 : SEGACT KRVIMP
175 : ENDIF
176 : IF (LCLITO) THEN
177 : SEGACT KRTOIM
178 : ENDIF
179 : SEGACT NOMINC
180 : SEGACT KRCENT
181 : SEGACT KRFACE
182 : SEGACT MELEFL
183 : SEGACT MPSURF
184 : SEGACT MPNORM
185 : SEGACT MPVOLU
186 : SEGACT MPROC
187 : SEGACT MPVITC
188 : SEGACT ICOGRV
189 : NSOUCH=ICOGRV.IMACHE(/1)
190 : DO 1 ISOUCH=1,NSOUCH
191 : MCOGRV=ICOGRV.IMACHE(ISOUCH)
192 : JCOGRV=ICOGRV.ICHAML(ISOUCH)
193 : SEGACT JCOGRV
194 : KDUNDX=JCOGRV.IELVAL(1)
195 : KDUNDY=JCOGRV.IELVAL(2)
196 : KDUNDZ=JCOGRV.IELVAL(3)
197 : SEGDES JCOGRV
198 : SEGACT KDUNDX
199 : SEGACT KDUNDY
200 : SEGACT KDUNDZ
201 : SEGACT MCOGRV
202 : NELEM=MCOGRV.NUM(/2)
203 : NPTEL=MCOGRV.NUM(/1)
204 : NPP1=NPTEL-1
205 : NPP2=NPTEL-1
206 : NEL1=NELEM
207 : NEL2=NELEM
208 : IEL1=1
209 : IEL2=1
210 : SEGINI ROURHO
211 : SEGINI ROVRHO
212 : SEGINI ROWRHO
213 : SEGINI ROUROU
214 : SEGINI ROVROU
215 : SEGINI ROWROU
216 : SEGINI ROUROV
217 : SEGINI ROVROV
218 : SEGINI ROWROV
219 : SEGINI ROUROW
220 : SEGINI ROVROW
221 : SEGINI ROWROW
222 : SEGINI CTSGRV
223 : ROURHO.NOMPRI(1:4)=NOMINC.MOTS(1)
224 : ROURHO.NOMPRI(5:8)=' '
225 : ROURHO.NOMDUA(1:4)=NOMINC.MOTS(2)
226 : ROURHO.NOMDUA(5:8)=' '
227 : ROVRHO.NOMPRI(1:4)=NOMINC.MOTS(1)
228 : ROVRHO.NOMPRI(5:8)=' '
229 : ROVRHO.NOMDUA(1:4)=NOMINC.MOTS(3)
230 : ROVRHO.NOMDUA(5:8)=' '
231 : ROWRHO.NOMPRI(1:4)=NOMINC.MOTS(1)
232 : ROWRHO.NOMPRI(5:8)=' '
233 : ROWRHO.NOMDUA(1:4)=NOMINC.MOTS(4)
234 : ROWRHO.NOMDUA(5:8)=' '
235 : ROUROU.NOMPRI(1:4)=NOMINC.MOTS(2)
236 : ROUROU.NOMPRI(5:8)=' '
237 : ROUROU.NOMDUA(1:4)=NOMINC.MOTS(2)
238 : ROUROU.NOMDUA(5:8)=' '
239 : ROVROU.NOMPRI(1:4)=NOMINC.MOTS(2)
240 : ROVROU.NOMPRI(5:8)=' '
241 : ROVROU.NOMDUA(1:4)=NOMINC.MOTS(3)
242 : ROVROU.NOMDUA(5:8)=' '
243 : ROWROU.NOMPRI(1:4)=NOMINC.MOTS(2)
244 : ROWROU.NOMPRI(5:8)=' '
245 : ROWROU.NOMDUA(1:4)=NOMINC.MOTS(4)
246 : ROWROU.NOMDUA(5:8)=' '
247 : ROUROV.NOMPRI(1:4)=NOMINC.MOTS(3)
248 : ROUROV.NOMPRI(5:8)=' '
249 : ROUROV.NOMDUA(1:4)=NOMINC.MOTS(2)
250 : ROUROV.NOMDUA(5:8)=' '
251 : ROVROV.NOMPRI(1:4)=NOMINC.MOTS(3)
252 : ROVROV.NOMPRI(5:8)=' '
253 : ROVROV.NOMDUA(1:4)=NOMINC.MOTS(3)
254 : ROVROV.NOMDUA(5:8)=' '
255 : ROWROV.NOMPRI(1:4)=NOMINC.MOTS(3)
256 : ROWROV.NOMPRI(5:8)=' '
257 : ROWROV.NOMDUA(1:4)=NOMINC.MOTS(4)
258 : ROWROV.NOMDUA(5:8)=' '
259 : ROUROW.NOMPRI(1:4)=NOMINC.MOTS(4)
260 : ROUROW.NOMPRI(5:8)=' '
261 : ROUROW.NOMDUA(1:4)=NOMINC.MOTS(2)
262 : ROUROW.NOMDUA(5:8)=' '
263 : ROVROW.NOMPRI(1:4)=NOMINC.MOTS(4)
264 : ROVROW.NOMPRI(5:8)=' '
265 : ROVROW.NOMDUA(1:4)=NOMINC.MOTS(3)
266 : ROVROW.NOMDUA(5:8)=' '
267 : ROWROW.NOMPRI(1:4)=NOMINC.MOTS(4)
268 : ROWROW.NOMPRI(5:8)=' '
269 : ROWROW.NOMDUA(1:4)=NOMINC.MOTS(4)
270 : ROWROW.NOMDUA(5:8)=' '
271 : DO 12 IELEM=1,NELEM
272 : * Le premier point du support de ICOGRV est un point FACE
273 : NGFACE=MCOGRV.NUM(1,IELEM)
274 : NLFACE=KRFACE.LECT(NGFACE)
275 : IF (NLFACE.EQ.0) THEN
276 : WRITE(IOIMP,*) 'Erreur de programmation n°1'
277 : GOTO 9999
278 : ENDIF
279 : * On calcule la contribution à la matrice jacobienne IJACO de la face
280 : * NGFAC (points duaux : centres à gauche et à droite de la face)
281 : * (points primaux : une partie (bicoz conditions aux limites)
282 : * de ceux du stencil pour le calcul du gradient
283 : * à la face, ils doivent être des points centres)
284 : * Si le tenseur des contraintes sur la face est imposé par les
285 : * conditions aux limites, la contribution de la face à IJACO est
286 : * nulle.
287 : LCTRB1=.TRUE.
288 : IF (LCLITO) THEN
289 : NLCLTO=KRTOIM.LECT(NGFACE)
290 : IF (NLCLTO.NE.0) THEN
291 : LCTRB1=.FALSE.
292 : ENDIF
293 : ENDIF
294 : IF (LCTRB1) THEN
295 : NGCGAU=MELEFL.NUM(1,NLFACE)
296 : NGCDRO=MELEFL.NUM(3,NLFACE)
297 : LMUR=(NGCGAU.EQ.NGCDRO)
298 : * On distingue le cas où la face est un bord du maillage (mur)
299 : * du cas où la face est interne au maillage
300 : IF (.NOT.LMUR) THEN
301 : NPD=2
302 : ELSE
303 : NPD=1
304 : ENDIF
305 : NPP=NPTEL-1
306 : * IPD=1 : point à gauche du point NGFACE
307 : * IPD=2 : point à droite du point NGFACE
308 : DO 122 IPD=1,NPD
309 : NPDUAL=MELEFL.NUM((2*IPD)-1,NLFACE)
310 : IF (.NOT.LMUR) THEN
311 : CTSGRV.POIDU2(IPD,IEL2)=NPDUAL
312 : ELSE
313 : CTSGRV.POIDU1(IPD,IEL1)=NPDUAL
314 : ENDIF
315 : NLCEND=KRCENT.LECT(NPDUAL)
316 : IF (NLCEND.EQ.0) THEN
317 : WRITE(IOIMP,*) 'Erreur grave n°1'
318 : GOTO 9999
319 : ENDIF
320 : DO 124 IPP=1,NPP
321 : NPPRIM=MCOGRV.NUM(IPP+1,IELEM)
322 : NLCENP=KRCENT.LECT(NPPRIM)
323 : C
324 : C AB We do not check BC anymore
325 : C
326 : IF(NLCENP .EQ. 0)THEN
327 : * Lorsque une contribution est nulle, on fixe artificiellement le
328 : * point primal égal au point dual.
329 : IF (.NOT.LMUR) THEN
330 : CTSGRV.POIPR2(IPP,IEL2)=NPDUAL
331 : ROURHO.VALMA2(IPD,IPP,IEL2)=0.D0
332 : ROVRHO.VALMA2(IPD,IPP,IEL2)=0.D0
333 : ROWRHO.VALMA2(IPD,IPP,IEL2)=0.D0
334 : ROUROU.VALMA2(IPD,IPP,IEL2)=0.D0
335 : ROVROU.VALMA2(IPD,IPP,IEL2)=0.D0
336 : ROWROU.VALMA2(IPD,IPP,IEL2)=0.D0
337 : ROUROV.VALMA2(IPD,IPP,IEL2)=0.D0
338 : ROVROV.VALMA2(IPD,IPP,IEL2)=0.D0
339 : ROWROV.VALMA2(IPD,IPP,IEL2)=0.D0
340 : ROUROW.VALMA2(IPD,IPP,IEL2)=0.D0
341 : ROVROW.VALMA2(IPD,IPP,IEL2)=0.D0
342 : ROWROW.VALMA2(IPD,IPP,IEL2)=0.D0
343 : ELSE
344 : CTSGRV.POIPR1(IPP,IEL1)=NPDUAL
345 : ROURHO.VALMA1(IPD,IPP,IEL1)=0.D0
346 : ROVRHO.VALMA1(IPD,IPP,IEL1)=0.D0
347 : ROWRHO.VALMA1(IPD,IPP,IEL1)=0.D0
348 : ROUROU.VALMA1(IPD,IPP,IEL1)=0.D0
349 : ROVROU.VALMA1(IPD,IPP,IEL1)=0.D0
350 : ROWROU.VALMA1(IPD,IPP,IEL1)=0.D0
351 : ROUROV.VALMA1(IPD,IPP,IEL1)=0.D0
352 : ROVROV.VALMA1(IPD,IPP,IEL1)=0.D0
353 : ROWROV.VALMA1(IPD,IPP,IEL1)=0.D0
354 : ROUROW.VALMA1(IPD,IPP,IEL1)=0.D0
355 : ROVROW.VALMA1(IPD,IPP,IEL1)=0.D0
356 : ROWROW.VALMA1(IPD,IPP,IEL1)=0.D0
357 : ENDIF
358 : ELSE
359 : * Les contributions valent :
360 : * 1. (d Res_{\rho u})_d / (d var)_p =
361 : * +/-1 (normale sortante, rentrante) (1/V_d) * (S_f) * \mu
362 : * * [ [ ( 4/3 (n_x * \alpha_x) + (n_y * \alpha_y)
363 : * + (n_z * \alpha_z)) *
364 : * ((du)_p / (d var)_p)]
365 : * + [ (-2/3 (n_x * \alpha_y) + (n_y * \alpha_x)) *
366 : * ((dv)_p / (d var)_p)]
367 : * + [ (-2/3 (n_x * \alpha_z) + (n_z * \alpha_x)) *
368 : * ((dw)_p / (d var)_p)]
369 : * ]
370 : * 2. (d Res_{\rho v})_d / (d var)_p =
371 : * +/-1 (normale sortante, rentrante) (1/V_d) * (S_f) * \mu
372 : * * [ [ (-2/3 (n_y * \alpha_x) + (n_x * \alpha_y)) *
373 : * ((du)_p / (d var)_p)]
374 : * + [ ( 4/3 (n_y * \alpha_y) + (n_x * \alpha_x)
375 : * + (n_z * \alpha_z)) *
376 : * ((dv)_p / (d var)_p)]
377 : * + [ (-2/3 (n_y * \alpha_z) + (n_z * \alpha_y)) *
378 : * ((dw)_p / (d var)_p)]
379 : * ]
380 : * 3. (d Res_{\rho w})_d / (d var)_p =
381 : * +/-1 (normale sortante, rentrante) (1/V_d) * (S_f) * \mu
382 : * * [ [ (-2/3 (n_z * \alpha_x) + (n_x * \alpha_z)) *
383 : * ((du)_p / (d var)_p)]
384 : * + [ (-2/3 (n_z * \alpha_y) + (n_y * \alpha_z)) *
385 : * ((dv)_p / (d var)_p)]
386 : * + [ ( 4/3 (n_z * \alpha_z) + (n_x * \alpha_x)
387 : * + (n_y * \alpha_y)) *
388 : * ((dw)_p / (d var)_p)]
389 : * ]
390 : *
391 : * avec :
392 : * (du)_p / (d \rho)_p = - u / \rho_p
393 : * (du)_p / (d \rho u)_p = 1 / \rho_p
394 : * (du)_p / (d \rho v)_p = 0
395 : * (du)_p / (d \rho w)_p = 0
396 : * (dv)_p / (d \rho)_p = - v / \rho_p
397 : * (dv)_p / (d \rho u)_p = 0
398 : * (dv)_p / (d \rho v)_p = 1 / \rho_p
399 : * (dv)_p / (d \rho w)_p = 0
400 : * (dv)_p / (d \rho)_p = - w / \rho_p
401 : * (dv)_p / (d \rho u)_p = 0
402 : * (dv)_p / (d \rho v)_p = 0
403 : * (dv)_p / (d \rho w)_p = 1 / \rho_p
404 : *
405 : * normale sortante pour IPD=1, rentrante pour IPD=2
406 : SIGNOR=(-1.D0)**(IPD+1)
407 : VOLUEL=MPVOLU.VPOCHA(NLCEND,1)
408 : SURFFA=MPSURF.VPOCHA(NLFACE,1)
409 : CNX =MPNORM.VPOCHA(NLFACE,1)
410 : CNY =MPNORM.VPOCHA(NLFACE,2)
411 : CNZ =MPNORM.VPOCHA(NLFACE,3)
412 : ALPHAX=KDUNDX.VELCHE(IPP+1,IELEM)
413 : ALPHAY=KDUNDY.VELCHE(IPP+1,IELEM)
414 : ALPHAZ=KDUNDZ.VELCHE(IPP+1,IELEM)
415 : RHOP =MPROC.VPOCHA(NLCENP,1)
416 : UP =MPVITC.VPOCHA(NLCENP,1)
417 : VP =MPVITC.VPOCHA(NLCENP,2)
418 : WP =MPVITC.VPOCHA(NLCENP,3)
419 : FACTOR=SIGNOR*(1.D0/VOLUEL)*SURFFA*MU
420 : BROUDU=(( 4.D0/3.D0)*(CNX*ALPHAX))
421 : $ + (CNY*ALPHAY) + (CNZ*ALPHAZ)
422 : BROUDV=((-2.D0/3.D0)*(CNX*ALPHAY))
423 : $ + (CNY*ALPHAX)
424 : BROUDW=((-2.D0/3.D0)*(CNX*ALPHAZ))
425 : $ + (CNZ*ALPHAX)
426 : BROVDU=((-2.D0/3.D0)*(CNY*ALPHAX))
427 : $ + (CNX*ALPHAY)
428 : BROVDV=(( 4.D0/3.D0)*(CNY*ALPHAY))
429 : $ + (CNX*ALPHAX) + (CNZ*ALPHAZ)
430 : BROVDW=((-2.D0/3.D0)*(CNY*ALPHAZ))
431 : $ + (CNZ*ALPHAY)
432 : BROWDU=((-2.D0/3.D0)*(CNZ*ALPHAX))
433 : $ + (CNX*ALPHAZ)
434 : BROWDV=((-2.D0/3.D0)*(CNZ*ALPHAY))
435 : $ + (CNY*ALPHAZ)
436 : BROWDW=(( 4.D0/3.D0)*(CNZ*ALPHAZ))
437 : $ + (CNX*ALPHAX) + (CNY*ALPHAY)
438 : DUDRHO=-UP /RHOP
439 : DUDROU=1.D0/RHOP
440 : DVDRHO=-VP /RHOP
441 : DVDROV=1.D0/RHOP
442 : DWDRHO=-WP /RHOP
443 : DWDROW=1.D0/RHOP
444 : IF (.NOT.LMUR) THEN
445 : CTSGRV.POIPR2(IPP,IEL2)=NPPRIM
446 : ROURHO.VALMA2(IPD,IPP,IEL2)=
447 : $ FACTOR*((BROUDU*DUDRHO)
448 : $ +(BROUDV*DVDRHO)+(BROUDW*DWDRHO))
449 : ROVRHO.VALMA2(IPD,IPP,IEL2)=
450 : $ FACTOR*((BROVDU*DUDRHO)
451 : $ +(BROVDV*DVDRHO)+(BROVDW*DWDRHO))
452 : ROWRHO.VALMA2(IPD,IPP,IEL2)=
453 : $ FACTOR*((BROWDU*DUDRHO)
454 : $ +(BROWDV*DVDRHO)+(BROWDW*DWDRHO))
455 : ROUROU.VALMA2(IPD,IPP,IEL2)=
456 : $ FACTOR* (BROUDU*DUDROU)
457 : ROVROU.VALMA2(IPD,IPP,IEL2)=
458 : $ FACTOR* (BROVDU*DUDROU)
459 : ROWROU.VALMA2(IPD,IPP,IEL2)=
460 : $ FACTOR* (BROWDU*DUDROU)
461 : ROUROV.VALMA2(IPD,IPP,IEL2)=
462 : $ FACTOR* (BROUDV*DVDROV)
463 : ROVROV.VALMA2(IPD,IPP,IEL2)=
464 : $ FACTOR* (BROVDV*DVDROV)
465 : ROWROV.VALMA2(IPD,IPP,IEL2)=
466 : $ FACTOR* (BROWDV*DVDROV)
467 : ROUROW.VALMA2(IPD,IPP,IEL2)=
468 : $ FACTOR* (BROUDW*DWDROW)
469 : ROVROW.VALMA2(IPD,IPP,IEL2)=
470 : $ FACTOR* (BROVDW*DWDROW)
471 : ROWROW.VALMA2(IPD,IPP,IEL2)=
472 : $ FACTOR* (BROWDW*DWDROW)
473 : ELSE
474 : CTSGRV.POIPR1(IPP,IEL1)=NPPRIM
475 : ROURHO.VALMA1(IPD,IPP,IEL1)=
476 : $ FACTOR*((BROUDU*DUDRHO)
477 : $ +(BROUDV*DVDRHO)+(BROUDW*DWDRHO))
478 : ROVRHO.VALMA1(IPD,IPP,IEL1)=
479 : $ FACTOR*((BROVDU*DUDRHO)
480 : $ +(BROVDV*DVDRHO)+(BROVDW*DWDRHO))
481 : ROWRHO.VALMA1(IPD,IPP,IEL1)=
482 : $ FACTOR*((BROWDU*DUDRHO)
483 : $ +(BROWDV*DVDRHO)+(BROWDW*DWDRHO))
484 : ROUROU.VALMA1(IPD,IPP,IEL1)=
485 : $ FACTOR* (BROUDU*DUDROU)
486 : ROVROU.VALMA1(IPD,IPP,IEL1)=
487 : $ FACTOR* (BROVDU*DUDROU)
488 : ROWROU.VALMA1(IPD,IPP,IEL1)=
489 : $ FACTOR* (BROWDU*DUDROU)
490 : ROUROV.VALMA1(IPD,IPP,IEL1)=
491 : $ FACTOR* (BROUDV*DVDROV)
492 : ROVROV.VALMA1(IPD,IPP,IEL1)=
493 : $ FACTOR* (BROVDV*DVDROV)
494 : ROWROV.VALMA1(IPD,IPP,IEL1)=
495 : $ FACTOR* (BROWDV*DVDROV)
496 : ROUROW.VALMA1(IPD,IPP,IEL1)=
497 : $ FACTOR* (BROUDW*DWDROW)
498 : ROVROW.VALMA1(IPD,IPP,IEL1)=
499 : $ FACTOR* (BROVDW*DWDROW)
500 : ROWROW.VALMA1(IPD,IPP,IEL1)=
501 : $ FACTOR* (BROWDW*DWDROW)
502 : ENDIF
503 : ENDIF
504 : 124 CONTINUE
505 : 122 CONTINUE
506 : IF (.NOT.LMUR) THEN
507 : IEL2=IEL2+1
508 : ELSE
509 : IEL1=IEL1+1
510 : ENDIF
511 : ENDIF
512 : 12 CONTINUE
513 : NPP1=NPTEL-1
514 : NPP2=NPTEL-1
515 : NEL1=IEL1-1
516 : NEL2=IEL2-1
517 : SEGADJ ROURHO
518 : SEGADJ ROVRHO
519 : SEGADJ ROWRHO
520 : SEGADJ ROUROU
521 : SEGADJ ROVROU
522 : SEGADJ ROWROU
523 : SEGADJ ROUROV
524 : SEGADJ ROVROV
525 : SEGADJ ROWROV
526 : SEGADJ ROUROW
527 : SEGADJ ROVROW
528 : SEGADJ ROWROW
529 : SEGADJ CTSGRV
530 : CTSGRV.LMATSI(**)=ROURHO
531 : CTSGRV.LMATSI(**)=ROVRHO
532 : CTSGRV.LMATSI(**)=ROWRHO
533 : CTSGRV.LMATSI(**)=ROUROU
534 : CTSGRV.LMATSI(**)=ROVROU
535 : CTSGRV.LMATSI(**)=ROWROU
536 : CTSGRV.LMATSI(**)=ROUROV
537 : CTSGRV.LMATSI(**)=ROVROV
538 : CTSGRV.LMATSI(**)=ROWROV
539 : CTSGRV.LMATSI(**)=ROUROW
540 : CTSGRV.LMATSI(**)=ROVROW
541 : CTSGRV.LMATSI(**)=ROWROW
542 : * On accumule les matrices résultantes dans IJACO
542 : CALL AJMTK(CTSGRV,IJACO,IMPR,IRET)
543 : IF (IRET.NE.0) GOTO 9999
544 : SEGSUP ROURHO
545 : SEGSUP ROVRHO
546 : SEGSUP ROWRHO
547 : SEGSUP ROUROU
548 : SEGSUP ROVROU
549 : SEGSUP ROWROU
550 : SEGSUP ROUROV
551 : SEGSUP ROVROV
552 : SEGSUP ROWROV
553 : SEGSUP ROUROW
554 : SEGSUP ROVROW
555 : SEGSUP ROWROW
556 : SEGSUP CTSGRV
557 : *
558 : SEGDES MCOGRV
559 : SEGDES KDUNDZ
560 : SEGDES KDUNDY
561 : SEGDES KDUNDX
562 : 1 CONTINUE
563 : SEGDES ICOGRV
564 : SEGDES MPVITC
565 : SEGDES MPROC
566 : SEGDES MPVOLU
567 : SEGDES MPNORM
568 : SEGDES MPSURF
569 : SEGDES MELEFL
570 : SEGDES KRFACE
571 : SEGDES KRCENT
572 : SEGDES NOMINC
573 : IF (LCLITO) THEN
574 : SEGDES KRTOIM
575 : ENDIF
576 : IF (LCLIMV) THEN
577 : SEGDES KRVIMP
578 : ENDIF
579 : *
580 : * Normal termination
581 : *
582 : IRET=0
583 : RETURN
584 : *
585 : * Format handling
586 : *
587 : *
588 : * Error handling
589 : *
590 : 9999 CONTINUE
591 : IRET=1
592 : WRITE(IOIMP,*) 'An error was detected in subroutine ylap2c'
593 : RETURN
594 : *
595 : * End of subroutine YLAP2C
596 : *
597 : END
598 :
599 :
600 :
601 :
602 :
603 :
604 :
605 :
© Cast3M 2003 - All rights reserved.
Disclaimer