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