1 : C ZCVFPU SOURCE CHAT 05/01/13 04:21:45 5004 2 : SUBROUTINE ZCVFPU(NEL,K0,NPF,IES,NPTD,IAXI,IPADL, 3 : & LEF,XCOOR, 4 : & VOLF, 5 : & UN,TK,TE, 6 : & F, 7 : & DK,DE, 8 : & ANU,IKC,UET,YP) 9 : IMPLICIT INTEGER(I-N) 10 : IMPLICIT REAL*8 (A-H,O-Z) 11 : 12 : C************************************************************************ 13 : C 14 : C SYNTAXE : 15 : C 16 : C FPU (NU,UET,YP <,VPAROI> ) 17 : C 18 : C 1------2 19 : C (R1,AL1) LEF FLUIDE NOEUDS 1 2 20 : C 21 : C 22 : C ANU VISCOSITE CINEMATIQUE 23 : C UET U* 24 : C YP DISTANCE A LA PAROI 25 : C VPAROI VITESSE DE LA PAROI (PAR DEFAUT 0.D0) 26 : C 27 : C CAS TRIDIMENSIONNEL 28 : C 4 ________ 3 29 : C / FLUIDE / 30 : C 1 /________/2 31 : C 32 : C 33 : C************************************************************************ 34 : -INC CCOPTIO 35 : -INC CCREEL 36 : DIMENSION UN(NPTD,IES),TK(*),TE(*) 37 : 38 : REAL*8 F(NPTD,IES) 39 : REAL*8 DK(*),DE(*) 40 : 41 : DIMENSION UET(*) 42 : DIMENSION AF(2,2),CF(2,2) 43 : DIMENSION ANU(*),VOLF(*) 44 : DIMENSION XCOOR(*),XYZ(3,4) 45 : DIMENSION LEF(NPF,NEL),IPADL(*) 46 : DIMENSION LE(4) 47 : DIMENSION UIX(4),UIY(4),UIZ(4),UFI(4) 48 : 49 : DIMENSION UPIU(101) 50 : 51 : C -------------------CAROLI----I------------------------------------ 52 : DATA UPIU( 1) /0.000000D+00/ 53 : DATA UPIU( 2) /0.999954D+00/ 54 : DATA UPIU( 3) /0.199858D+01/ 55 : DATA UPIU( 4) /0.298974D+01/ 56 : DATA UPIU( 5) /0.395957D+01/ 57 : DATA UPIU( 6) /0.488802D+01/ 58 : DATA UPIU( 7) /0.575502D+01/ 59 : DATA UPIU( 8) /0.654720D+01/ 60 : DATA UPIU( 9) /0.726016D+01/ 61 : DATA UPIU( 10) /0.789663D+01/ 62 : DATA UPIU( 11) /0.846323D+01/ 63 : DATA UPIU( 12) /0.896802D+01/ 64 : DATA UPIU( 13) /0.941906D+01/ 65 : DATA UPIU( 14) /0.982372D+01/ 66 : DATA UPIU( 15) /0.101884D+02/ 67 : DATA UPIU( 16) /0.105188D+02/ 68 : DATA UPIU( 17) /0.108193D+02/ 69 : DATA UPIU( 18) /0.110941D+02/ 70 : DATA UPIU( 19) /0.113464D+02/ 71 : DATA UPIU( 20) /0.115790D+02/ 72 : DATA UPIU( 21) /0.117944D+02/ 73 : DATA UPIU( 22) /0.119944D+02/ 74 : DATA UPIU( 23) /0.121809D+02/ 75 : DATA UPIU( 24) /0.123553D+02/ 76 : DATA UPIU( 25) /0.125189D+02/ 77 : DATA UPIU( 26) /0.126727D+02/ 78 : DATA UPIU( 27) /0.128178D+02/ 79 : DATA UPIU( 28) /0.129549D+02/ 80 : DATA UPIU( 29) /0.130849D+02/ 81 : DATA UPIU( 30) /0.132082D+02/ 82 : DATA UPIU( 31) /0.133256D+02/ 83 : DATA UPIU( 32) /0.134375D+02/ 84 : DATA UPIU( 33) /0.135443D+02/ 85 : DATA UPIU( 34) /0.136465D+02/ 86 : DATA UPIU( 35) /0.137444D+02/ 87 : DATA UPIU( 36) /0.138383D+02/ 88 : DATA UPIU( 37) /0.139285D+02/ 89 : DATA UPIU( 38) /0.140153D+02/ 90 : DATA UPIU( 39) /0.140988D+02/ 91 : DATA UPIU( 40) /0.141794D+02/ 92 : DATA UPIU( 41) /0.142573D+02/ 93 : DATA UPIU( 42) /0.143325D+02/ 94 : DATA UPIU( 43) /0.144052D+02/ 95 : DATA UPIU( 44) /0.144757D+02/ 96 : DATA UPIU( 45) /0.145440D+02/ 97 : DATA UPIU( 46) /0.146102D+02/ 98 : DATA UPIU( 47) /0.146746D+02/ 99 : DATA UPIU( 48) /0.147371D+02/ 100 : DATA UPIU( 49) /0.147979D+02/ 101 : DATA UPIU( 50) /0.148570D+02/ 102 : DATA UPIU( 51) /0.149147D+02/ 103 : DATA UPIU( 52) /0.149708D+02/ 104 : DATA UPIU( 53) /0.150256D+02/ 105 : DATA UPIU( 54) /0.150790D+02/ 106 : DATA UPIU( 55) /0.151311D+02/ 107 : DATA UPIU( 56) /0.151821D+02/ 108 : DATA UPIU( 57) /0.152319D+02/ 109 : DATA UPIU( 58) /0.152806D+02/ 110 : DATA UPIU( 59) /0.153282D+02/ 111 : DATA UPIU( 60) /0.153749D+02/ 112 : DATA UPIU( 61) /0.154206D+02/ 113 : DATA UPIU( 62) /0.154653D+02/ 114 : DATA UPIU( 63) /0.155092D+02/ 115 : DATA UPIU( 64) /0.155522D+02/ 116 : DATA UPIU( 65) /0.155944D+02/ 117 : DATA UPIU( 66) /0.156358D+02/ 118 : DATA UPIU( 67) /0.156765D+02/ 119 : DATA UPIU( 68) /0.157164D+02/ 120 : DATA UPIU( 69) /0.157556D+02/ 121 : DATA UPIU( 70) /0.157942D+02/ 122 : DATA UPIU( 71) /0.158321D+02/ 123 : DATA UPIU( 72) /0.158694D+02/ 124 : DATA UPIU( 73) /0.159060D+02/ 125 : DATA UPIU( 74) /0.159421D+02/ 126 : DATA UPIU( 75) /0.159776D+02/ 127 : DATA UPIU( 76) /0.160126D+02/ 128 : DATA UPIU( 77) /0.160470D+02/ 129 : DATA UPIU( 78) /0.160809D+02/ 130 : DATA UPIU( 79) /0.161143D+02/ 131 : DATA UPIU( 80) /0.161472D+02/ 132 : DATA UPIU( 81) /0.161797D+02/ 133 : DATA UPIU( 82) /0.162117D+02/ 134 : DATA UPIU( 83) /0.162433D+02/ 135 : DATA UPIU( 84) /0.162744D+02/ 136 : DATA UPIU( 85) /0.163051D+02/ 137 : DATA UPIU( 86) /0.163354D+02/ 138 : DATA UPIU( 87) /0.163653D+02/ 139 : DATA UPIU( 88) /0.163949D+02/ 140 : DATA UPIU( 89) /0.164240D+02/ 141 : DATA UPIU( 90) /0.164528D+02/ 142 : DATA UPIU( 91) /0.164813D+02/ 143 : DATA UPIU( 92) /0.165094D+02/ 144 : DATA UPIU( 93) /0.165371D+02/ 145 : DATA UPIU( 94) /0.165646D+02/ 146 : DATA UPIU( 95) /0.165917D+02/ 147 : DATA UPIU( 96) /0.166185D+02/ 148 : DATA UPIU( 97) /0.166450D+02/ 149 : DATA UPIU( 98) /0.166712D+02/ 150 : DATA UPIU( 99) /0.166971D+02/ 151 : DATA UPIU(100) /0.167228D+02/ 152 : DATA UPIU(101) /0.167481D+02/ 153 : C -----------------------CAROLI---------F----------------------- 154 : 155 : C CONSTANTES PHYSIQUES 156 : C 157 : DATA CNU,AKER /0.09D0,0.41D0/ 158 : 159 : C 160 : IF(IKC.NE.1)GO TO 9990 161 : R1=1.D0 162 : C 163 : C 164 : SQCNU=SQRT(CNU) 165 : KC=1 166 : C --------------------CAROLI-----------I----------------------- 167 : ANUL=ANU(KC) 168 : C --------------------CAROLI-----------F----------------------- 169 : 170 : IF(IES.EQ.3)GO TO 300 171 : 172 : C ********* 173 : C * 2D * 174 : C ********* 175 : 176 : CF(1,1)=R1/3.D0 177 : CF(2,1)=CF(1,1)/2.D0 178 : CF(1,2)=CF(2,1) 179 : CF(2,2)=CF(1,1) 180 : 181 : C IF(IPT.EQ.100)WRITE(IOIMP,*)' YCVFPU : ANUL,R1 ',ANUL,R1 182 : 183 : DO 50 K=1,NEL 184 : NK=K0+K 185 : C WRITE(IOIMP,*)' XCVFPS NK=',NK 186 : KC=1+(1-IKC)*(NK-1) 186 : CALL RSETD(AF,CF,4) 187 : C 188 : C 189 : IP=LEF(1,K) 190 : XYZ(1,1)=XCOOR((IP-1)*(IES+1) +1) 191 : XYZ(2,1)=XCOOR((IP-1)*(IES+1) +2) 192 : IP=LEF(2,K) 193 : XYZ(1,2)=XCOOR((IP-1)*(IES+1) +1) 194 : XYZ(2,2)=XCOOR((IP-1)*(IES+1) +2) 195 : 196 : LE(1)=IPADL(LEF(1,K)) 197 : LE(2)=IPADL(LEF(2,K)) 198 : 199 : C WRITE(IOIMP,*)' LE.... ',(LE(MM),MM=1,4) 200 : 201 : DO 5 I=1,2 202 : 203 : NF=LE(I) 204 : 205 : UIX(I)= UN(NF,1) 206 : UIY(I)= UN(NF,2) 207 : 208 : 5 CONTINUE 209 : 210 : C 211 : C CALCUL DES LONGUEURS DE L'ELEMENT : AL1 212 : C 213 : IF (IAXI.EQ.0) THEN 214 : R1=1.D0 215 : ELSEIF (IAXI.EQ.1) THEN 216 : R1=(XYZ(2,1)+XYZ(2,2))*XPI 217 : ELSEIF (IAXI.EQ.2) THEN 218 : R1=(XYZ(1,1)+XYZ(1,2))*XPI 219 : ENDIF 220 : 221 : AL1X=XYZ(1,2)-XYZ(1,1) 222 : AL1Y=XYZ(2,2)-XYZ(2,1) 223 : AL1=SQRT(AL1X*AL1X+AL1Y*AL1Y) 224 : 225 : C 226 : C CALCUL DES COMPOSANTES TANGENTIELLES A LA PAROI DES VITESSES 227 : C ET DE LA MOYENNE DE LEURS DIFFERENCES FLUIDE-PAROI : UTM 228 : C 229 : 230 : UT1=UIX(1)*UIX(1) + UIY(1)*UIY(1) 231 : UT2=UIX(2)*UIX(2) + UIY(2)*UIY(2) 232 : UT1=SQRT(UT1) 233 : UT2=SQRT(UT2) 234 : UTM=(UT1+UT2)/2.D0 235 : 236 : C 237 : C CALCUL DE U ETOILE A PARTIR UTM 238 : C 239 : C --------------------CAROLI---I-------------------------- 240 : DO 57 I=1,2 241 : YPLUS=YP*UET(K)/ANUL 242 : YPLUS=ABS(YPLUS)+1.D-5 243 : IF(YPLUS.GT.100.0D0) THEN 244 : UPLUS=5.5D0+(LOG(YPLUS))/0.41D0 245 : ELSE 246 : IY=INT(YPLUS) 247 : RY=YPLUS-IY 248 : UPLUS=UPIU(IY+1)+RY*(UPIU(IY+2)-UPIU(IY+1)) 249 : ENDIF 250 : C ---- RELAX SU UET ------------------------------------ 251 : UET(NK)=0.5D0*UET(NK)+0.5D0*(ABS(UTM/(UPLUS+1.D-5))) 252 : 57 CONTINUE 253 : C --------------------CAROLI----F------------------ 254 : 255 : C************************************************************************ 256 : C CALCUL Q D M 257 : C************************************************************************ 258 : 259 : C CALCUL DU COEFFICIENT AK A PARTIR DE UTM 260 : 261 : AUTM=ABS(UTM) 262 : AK=0.D0 263 : IF(AUTM.GT.1.D-10)AK=UET(NK)*UET(NK)/AUTM 264 : 265 : C CALCUL DE LA MATRICE VITESSE 266 : DO 82 I=1,2 267 : DO 72 J=1,2 268 : AF(J,I)=AK*AL1*AF(J,I)*R1 269 : 72 CONTINUE 270 : 82 CONTINUE 271 : C 272 : C SATURATION VITESSE 273 : C 274 : C IF(IPAT.LE.2)WRITE(IOIMP,*)' UTM=',UTM,AL1 275 : DO 54 J=1,2 276 : AF(1,J)=AF(1,J)*UT1 277 : AF(2,J)=AF(2,J)*UT2 278 : 54 CONTINUE 279 : 280 : DO 65 I=1,2 281 : NF=LE(I) 282 : FF1=(AF(1,I)+AF(2,I))*AL1X/AL1 283 : FF2=(AF(1,I)+AF(2,I))*AL1Y/AL1 284 : C if(k.ge.10.AND.k.le.15)then 285 : C WRITE(IOIMP,*)' FF1,FF2=',ff1,ff2 286 : C endif 287 : F(NF,1)=F(NF,1)+SIGN(FF1,UN(NF,1)) 288 : F(NF,2)=F(NF,2)+SIGN(FF2,UN(NF,2)) 289 : 290 : C------------------------------------------------------------------------ 291 : C CALCUL Q D M FIN 292 : C------------------------------------------------------------------------ 293 : 294 : 65 CONTINUE 295 : 296 : C 297 : C CALCUL DE KP ET DE EPSILONP 298 : C 299 : YPLUS2=-(YPLUS+0.01D0)*(YPLUS+0.01D0)*0.0017D0 300 : DO 66 I=1,2 301 : NF=LE(I) 302 : DK(NF)=XGRAND 303 : DE(NF)=XGRAND 304 : TK(NF)=UET(NK)*UET(NK)/SQCNU 305 : TE(NF)=UET(NK)*UET(NK)*UET(NK)/ 306 : & (AKER*YP*(1.D0-EXP(YPLUS2))) 307 : 66 CONTINUE 308 : 309 : 50 CONTINUE 310 : RETURN 311 : 312 : C ********* 313 : C * 3D * 314 : C ********* 315 : 300 CONTINUE 316 : 317 : DO 350 K=1,NEL 318 : NK=K0+K 319 : KC=1+(1-IKC)*(NK-1) 320 : C 321 : UTM=0.D0 322 : 323 : DO 35 I=1,NPF 324 : NFA=IPADL(LEF(I,K)) 325 : UIX(I)=UN(NFA,1)+XPETIT 326 : UIY(I)=UN(NFA,2)+XPETIT 327 : UIZ(I)=UN(NFA,3)+XPETIT 328 : UFI(I)=UIX(I)*UIX(I)+UIY(I)*UIY(I)+UIZ(I)*UIZ(I) 329 : UFI(I)=SQRT(UFI(I)) 330 : UTM=UTM+UFI(I) 331 : 35 CONTINUE 332 : 333 : UTM=UTM/NPF 334 : 335 : AIRF=VOLF(NK) 336 : C 337 : C CALCUL DE U ETOILE A PARTIR UTM 338 : C 339 : C --------------------CAROLI---I-------------------------- 340 : DO 357 I=1,2 341 : YPLUS=YP*UET(NK)/ANUL 342 : YPLUS=ABS(YPLUS)+1.D-5 343 : IF(YPLUS.GT.100.0D0) THEN 344 : UPLUS=5.5D0+(LOG(YPLUS))/0.41D0 345 : ELSE 346 : IY=INT(YPLUS) 347 : RY=YPLUS-IY 348 : UPLUS=UPIU(IY+1)+RY*(UPIU(IY+2)-UPIU(IY+1)) 349 : ENDIF 350 : C ---- RELAX SU UET ------------------------------------ 351 : UET(NK)=0.5D0*UET(NK)+0.5D0*(ABS(UTM/(UPLUS+1.D-5))) 352 : 357 CONTINUE 353 : C --------------------CAROLI----F------------------ 354 : 355 : C************************************************************************ 356 : C CALCUL Q D M 357 : C************************************************************************ 358 : 359 : C CALCUL DU COEFFICIENT AK A PARTIR DE UTM 360 : 361 : AUTM=ABS(UTM) 362 : AK=0.D0 363 : IF(AUTM.GT.1.D-10)AK=UET(NK)*UET(NK)/AUTM 364 : C 365 : DO 365 I=1,NPF 366 : NFA=IPADL(LEF(I,K)) 367 : F(NFA,1)=F(NFA,1)+AK/NPF*UIX(I)/UFI(I)*AIRF 368 : F(NFA,2)=F(NFA,2)+AK/NPF*UIY(I)/UFI(I)*AIRF 369 : F(NFA,3)=F(NFA,3)+AK/NPF*UIZ(I)/UFI(I)*AIRF 370 : 365 CONTINUE 371 : C------------------------------------------------------------------------ 372 : C CALCUL Q D M FIN 373 : C------------------------------------------------------------------------ 374 : C 375 : C CALCUL DE KP ET DE EPSILONP 376 : C 377 : YPLUS2=-(YPLUS+0.01D0)*(YPLUS+0.01D0)*0.0017D0 378 : DO 366 I=1,NPF 379 : NFA=IPADL(LEF(I,K)) 380 : DK(NFA)=XGRAND 381 : DE(NFA)=XGRAND 382 : TK(NFA)=UET(NK)*UET(NK)/SQCNU 383 : TE(NFA)=UET(NK)*UET(NK)*UET(NK)/ 384 : & (AKER*YP*(1.D0-EXP(YPLUS2))) 385 : 366 CONTINUE 386 : 387 : 350 CONTINUE 388 : 389 : RETURN 390 : 391 : 9990 CONTINUE 392 : WRITE(IOIMP,*)' IKC ',IKC 393 : WRITE(IOIMP,*) 394 : $ ' COEF D''UN MAUVAIS TYPE : ON VEUT (SCAL OU VECT DOMA)' 394 : CALL ARRET(0) 395 : RETURN 396 : 1002 FORMAT(10(1X,1PE11.4)) 397 : 1001 FORMAT(20(1X,I5)) 398 : END 399 : 400 : 401 : 402 : 403 : 404 :
© Cast3M 2003 - All rights reserved.
Disclaimer