1 : C CGADXF SOURCE GOUNAND 05/12/21 21:17:04 5281 2 : SUBROUTINE CGADXF(JAC,IREL,IREF, 3 : $ SDETG,IG,IGH,IGHIG, 4 : $ NOMLOI, 5 : $ CONTRI, 6 : $ IMPR,IRET) 7 : IMPLICIT REAL*8 (A-H,O-Z) 8 : IMPLICIT INTEGER (I-N) 9 : C*********************************************************************** 10 : C NOM : CGADXF 11 : C DESCRIPTION : Calcul de la loi de comportement aux points de Gauss : 12 : C Partie en Fortran 77 13 : C 14 : C LANGAGE : ESOPE 15 : C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF) 16 : C mél : gounand@semt2.smts.cea.fr 17 : C*********************************************************************** 18 : C APPELES : 19 : C APPELE PAR : 20 : C*********************************************************************** 21 : C ENTREES : 22 : C ENTREES/SORTIES : 23 : C SORTIES : - 24 : C TRAVAIL : 25 : C*********************************************************************** 26 : C VERSION : v1, 04/08/04, version initiale 27 : C HISTORIQUE : v1, 04/08/04, création 28 : C HISTORIQUE : 29 : C HISTORIQUE : 30 : C*********************************************************************** 31 : C Prière de PRENDRE LE TEMPS de compléter les commentaires 32 : C en cas de modification de ce sous-programme afin de faciliter 33 : C la maintenance ! 34 : C*********************************************************************** 35 : -INC CCOPTIO 36 : REAL*8 JAC(IREL,IREF) 37 : REAL*8 SDETG 38 : REAL*8 IG(IREF,IREF) 39 : REAL*8 IGH(IREF,IREF) 40 : REAL*8 IGHIG(IREF,IREF) 41 : CHARACTER*8 NOMLOI 42 : REAL*8 CONTRI 43 : * 44 : INTEGER IMPR,IRET 45 : * 46 : * Executable statements 47 : * 48 : IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans cgadxf' 49 : * 50 : 50 : CALL MARE(IGH,IREF,IREF,'TRACE ', 51 : $ TRIGH,IMPR,IRET) 52 : IF (IRET.NE.0) GOTO 9999 53 : * 54 : * Calcul du déterminant 55 : * 56 : IF (NOMLOI.EQ.'ADCD ') THEN 57 : CONTRI=SDETG 58 : * 59 : * Calcul de la fonctionnelle 60 : * 61 : ELSEIF (NOMLOI.EQ.'ADCF ') THEN 62 : CONTRI=0.5D0*TRIGH*SDETG 63 : * 64 : * Calcul du résidu 65 : * 66 : ELSEIF (NOMLOI(1:4).EQ.'ADCR') THEN 66 : CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET) 67 : IF (IRET.NE.0) GOTO 9999 67 : CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET) 68 : IF (IRET.NE.0) GOTO 9999 69 : TR1=0.D0 70 : DO I=1,IREF 71 : DO J=1,IREF 72 : TR1=TR1+ 73 : $ (IG(I,J)*JAC(IDIM1,J)*IGH(IDIM2,I)) 74 : ENDDO 75 : ENDDO 76 : TR2=0.D0 77 : DO I=1,IREF 78 : DO L=1,IREF 79 : TR2=TR2+ 80 : $ (IG(I,IDIM2)*JAC(IDIM1,L)*IGH(L,I)) 81 : ENDDO 82 : ENDDO 83 : CONTR1=-0.5D0*(TR1+TR2)*SDETG 84 : TR3=0.D0 85 : DO J=1,IREF 86 : TR3=TR3+ 87 : $ (IG(IDIM2,J)*JAC(IDIM1,J)) 88 : ENDDO 89 : TR4=0.D0 90 : DO I=1,IREF 91 : TR4=TR4+ 92 : $ (IG(I,IDIM2)*JAC(IDIM1,I)) 93 : ENDDO 94 : CONTR2=0.25D0*(TR3+TR4)*TRIGH*SDETG 95 : CONTRI=CONTR1+CONTR2 96 : * CONTRI=CONTR1 97 : C WRITE(IOIMP,*) 'ADCR' 98 : C WRITE(IOIMP,*) 'TR1=',TR1 99 : C WRITE(IOIMP,*) 'TR2=',TR2 100 : C WRITE(IOIMP,*) 'TR3=',TR3 101 : C WRITE(IOIMP,*) 'TR4=',TR4 102 : * 103 : * Calcul du jacobien 104 : * 105 : ELSEIF (NOMLOI(1:4).EQ.'ADCJ') THEN 105 : CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET) 106 : IF (IRET.NE.0) GOTO 9999 106 : CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET) 107 : IF (IRET.NE.0) GOTO 9999 107 : CALL CH2INT(NOMLOI(7:7),IDIM3,IMPR,IRET) 108 : IF (IRET.NE.0) GOTO 9999 108 : CALL CH2INT(NOMLOI(8:8),IDIM4,IMPR,IRET) 109 : IF (IRET.NE.0) GOTO 9999 110 : TR0=0.D0 111 : DO K=1,IREF 112 : TR0=TR0+ 113 : $ (JAC(IDIM1,K)*IG(K,IDIM2)) 114 : ENDDO 115 : * 116 : TR1=0.D0 117 : DO I=1,IREF 118 : DO J=1,IREF 119 : TR1=TR1+ 120 : $ (IG(I,J)*JAC(IDIM3,J)*IGH(IDIM4,I)) 121 : ENDDO 122 : ENDDO 123 : TR2=0.D0 124 : DO I=1,IREF 125 : DO M=1,IREF 126 : TR2=TR2+ 127 : $ (IG(I,IDIM4)*JAC(IDIM3,M)*IGH(M,I)) 128 : ENDDO 129 : ENDDO 130 : * 131 : * 132 : TR3=0.D0 133 : DO K=1,IREF 134 : DO L=1,IREF 135 : TR3=TR3+ 136 : $ (JAC(IDIM1,K)*IG(K,L)*JAC(IDIM3,L) 137 : $ *IGHIG(IDIM4,IDIM2)) 138 : ENDDO 139 : ENDDO 140 : TR3=TR3*(-1.D0) 141 : * 142 : TR4=0.D0 143 : DO K=1,IREF 144 : DO IO=1,IREF 145 : TR4=TR4+ 146 : $ (JAC(IDIM1,K)*IGHIG(K,IDIM4) 147 : $ *JAC(IDIM3,IO)*IG(IO,IDIM2)) 148 : ENDDO 149 : ENDDO 150 : TR4=TR4*(-1.D0) 151 : * 152 : TR5=0.D0 153 : IF (IDIM1.EQ.IDIM3) THEN 154 : TR5=TR5+ 155 : $ (IGHIG(IDIM4,IDIM2)) 156 : ENDIF 157 : * 158 : TR6=0.D0 159 : DO K=1,IREF 160 : DO M=1,IREF 161 : TR6=TR6+ 162 : $ (JAC(IDIM1,K)*IGHIG(K,M) 163 : $ *JAC(IDIM3,M)*IG(IDIM4,IDIM2)) 164 : ENDDO 165 : ENDDO 166 : TR6=TR6*(-1.D0) 167 : * 168 : TR7=0.D0 169 : DO K=1,IREF 170 : DO N=1,IREF 171 : TR7=TR7+ 172 : $ (JAC(IDIM1,K)*IG(K,IDIM4)*JAC(IDIM3,N) 173 : $ *IGHIG(N,IDIM2)) 174 : ENDDO 175 : ENDDO 176 : TR7=TR7*(-1.D0) 177 : CONTR1=(TR1+TR2)*TR0 178 : CONTR2=(TR3+TR4+TR5+TR6+TR7)*2.D0 179 : CONT11=(CONTR1+CONTR2)*SDETG 180 : **** 181 : **** 182 : TR10=0.D0 183 : DO K=1,IREF 184 : TR10=TR10+ 185 : $ (JAC(IDIM1,K)*IGHIG(K,IDIM2)) 186 : ENDDO 187 : TR10=TR10*(-2.D0) 188 : * 189 : TR0=0.D0 190 : DO K=1,IREF 191 : TR0=TR0+ 192 : $ (JAC(IDIM1,K)*IG(K,IDIM2)) 193 : ENDDO 194 : * 195 : TR11=0.D0 196 : DO J=1,IREF 197 : TR11=TR11+ 198 : $ (IG(IDIM4,J)*JAC(IDIM3,J)) 199 : ENDDO 200 : TR12=0.D0 201 : DO I=1,IREF 202 : TR12=TR12+ 203 : $ (IG(I,IDIM4)*JAC(IDIM3,I)) 204 : ENDDO 205 : 206 : * 207 : * 208 : TR13=0.D0 209 : DO K=1,IREF 210 : DO L=1,IREF 211 : TR13=TR13+ 212 : $ (JAC(IDIM1,K)*IG(K,L)*JAC(IDIM3,L) 213 : $ *IG(IDIM4,IDIM2)) 214 : ENDDO 215 : ENDDO 216 : TR13=TR13*(-1.D0) 217 : TR14=0.D0 218 : DO K=1,IREF 219 : DO N=1,IREF 220 : TR14=TR14+ 221 : $ (JAC(IDIM1,K)*IG(K,IDIM4) 222 : $ *JAC(IDIM3,N)*IG(N,IDIM2)) 223 : ENDDO 224 : ENDDO 225 : TR14=TR14*(-1.D0) 226 : * 227 : TR15=0.D0 228 : IF (IDIM1.EQ.IDIM3) THEN 229 : TR15=TR15+IG(IDIM4,IDIM2) 230 : ENDIF 231 : * 232 : CONTR3=(TR11+TR12)*(TR10+(TR0*TRIGH))*SDETG 233 : CONTR4=(TR13+TR14+TR15)*2.D0*TRIGH*SDETG 234 : CONT12=(CONTR3+CONTR4) 235 : ***** 236 : ***** Finalement ! 237 : ***** 238 : * WRITE(IOIMP,*) 'CONT11=',CONT11 239 : * WRITE(IOIMP,*) 'CONT12=',CONT12 240 : CONTRI=(CONT11*(-0.5D0)) 241 : $ +(CONT12*0.25D0) 242 : * CONTRI=(CONT11*(-0.5D0)) 243 : C* 244 : C* Calcul du résidu (test) 245 : C* 246 : C ELSEIF (NOMLOI(1:4).EQ.'ADCS') THEN 246 : C CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET) 247 : C IF (IRET.NE.0) GOTO 9999 247 : C CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET) 248 : C IF (IRET.NE.0) GOTO 9999 249 : CC TR1=JAC(IDIM1,IDIM2) 250 : C TR1=0.D0 251 : C DO I=1,IREF 252 : C DO J=1,IREF 253 : C TR1=TR1+ 254 : C $ (IG(I,J)*JAC(IDIM1,J)*IGH(IDIM2,I)) 255 : C* $ (IG(I,J)*JAC(IDIM1,J)*IG(IDIM2,I)) 256 : C ENDDO 257 : C ENDDO 258 : C TR2=0.D0 259 : C DO I=1,IREF 260 : C DO M=1,IREF 261 : C TR2=TR2+ 262 : C $ (IG(I,IDIM2)*JAC(IDIM1,M)*IGH(M,I)) 263 : C ENDDO 264 : C ENDDO 265 : C CONTRI=(TR1+TR2)*SDETG 266 : C* CONTRI=(TR1+TR2) 267 : C* CONTRI=TR1 268 : C* 269 : C* Calcul du jacobien (test) 270 : C* 271 : C ELSEIF (NOMLOI(1:4).EQ.'ADCK') THEN 271 : C CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET) 272 : C IF (IRET.NE.0) GOTO 9999 272 : C CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET) 273 : C IF (IRET.NE.0) GOTO 9999 273 : C CALL CH2INT(NOMLOI(7:7),IDIM3,IMPR,IRET) 274 : C IF (IRET.NE.0) GOTO 9999 274 : C CALL CH2INT(NOMLOI(8:8),IDIM4,IMPR,IRET) 275 : C IF (IRET.NE.0) GOTO 9999 276 : CC IF ((IDIM1.EQ.IDIM3).AND.(IDIM2.EQ.IDIM4)) THEN 277 : CC CONTRI=1.D0 278 : CC ELSE 279 : CC CONTRI=0.D0 280 : CC ENDIF 281 : C TR0=0.D0 282 : C DO K=1,IREF 283 : C TR0=TR0+ 284 : C $ (JAC(IDIM1,K)*IG(K,IDIM2)) 285 : C ENDDO 286 : C* 287 : C TR1=0.D0 288 : C DO I=1,IREF 289 : C DO J=1,IREF 290 : C TR1=TR1+ 291 : C $ (IG(I,J)*JAC(IDIM3,J)*IGH(IDIM4,I)) 292 : C ENDDO 293 : C ENDDO 294 : C TR2=0.D0 295 : C DO I=1,IREF 296 : C DO M=1,IREF 297 : C TR2=TR2+ 298 : C $ (IG(I,IDIM4)*JAC(IDIM3,M)*IGH(M,I)) 299 : C ENDDO 300 : C ENDDO 301 : C* 302 : C* 303 : C TR3=0.D0 304 : C DO K=1,IREF 305 : C DO L=1,IREF 306 : C DO IO=1,IREF 307 : C TR3=TR3+ 308 : C $ (JAC(IDIM1,K)*IG(K,L)*JAC(IDIM3,L) 309 : C $ *IG(IDIM4,IO)*HIG(IO,IDIM2)) 310 : C* $ (JAC(IDIM1,K)*IG(K,L)*JAC(IDIM3,L) 311 : C* $ *IG(IDIM4,IO)*IG(IO,IDIM2)) 312 : C ENDDO 313 : C ENDDO 314 : C ENDDO 315 : C TR3=TR3*(-1.D0) 316 : C* 317 : C TR4=0.D0 318 : C DO K=1,IREF 319 : C DO L=1,IREF 320 : C DO IO=1,IREF 321 : C TR4=TR4+ 322 : C $ (JAC(IDIM1,K)*IGH(K,L)*IG(L,IDIM4) 323 : C $ *JAC(IDIM3,IO)*IG(IO,IDIM2)) 324 : C* $ (JAC(IDIM1,K)*IG(K,L)*IG(L,IDIM4) 325 : C* $ *JAC(IDIM3,IO)*IG(IO,IDIM2)) 326 : C ENDDO 327 : C ENDDO 328 : C ENDDO 329 : C TR4=TR4*(-1.D0) 330 : C* 331 : C TR5=0.D0 332 : C IF (IDIM1.EQ.IDIM3) THEN 333 : C DO L=1,IREF 334 : C TR5=TR5+ 335 : C $ (IGH(IDIM4,L)*IG(L,IDIM2)) 336 : C* $ (IG(IDIM4,L)*IG(L,IDIM2)) 337 : C ENDDO 338 : C ENDIF 339 : C* 340 : C TR6=0.D0 341 : C DO K=1,IREF 342 : C DO L=1,IREF 343 : C DO M=1,IREF 344 : C TR6=TR6+ 345 : C $ (JAC(IDIM1,K)*IGH(K,L)*IG(L,M) 346 : C $ *JAC(IDIM3,M)*IG(IDIM4,IDIM2)) 347 : C ENDDO 348 : C ENDDO 349 : C ENDDO 350 : C TR6=TR6*(-1.D0) 351 : C* 352 : C TR7=0.D0 353 : C DO K=1,IREF 354 : C DO N=1,IREF 355 : C DO IO=1,IREF 356 : C TR7=TR7+ 357 : C $ (JAC(IDIM1,K)*IG(K,IDIM4)*JAC(IDIM3,N) 358 : C $ *IG(N,IO)*HIG(IO,IDIM2)) 359 : C ENDDO 360 : C ENDDO 361 : C ENDDO 362 : C TR7=TR7*(-1.D0) 363 : C* 364 : C CONTR1=(TR1+TR2)*TR0 365 : C CONTR2=(TR3+TR4+TR5+TR6+TR7)*2.D0 366 : C CONTRI=(CONTR1+CONTR2)*SDETG 367 : C* CONTRI=(TR3+TR4+TR5+TR6+TR7)+(TR3B+TR4B+TR5B+TR6B+TR7B) 368 : C* CONTRI=(TR3+TR4+TR5+TR6+TR7) 369 : C* 370 : C* WRITE(IOIMP,*) 'CONTRI=',CONTRI 371 : C* 372 : C* Calcul du résidu (test) 373 : C* 374 : C ELSEIF (NOMLOI(1:4).EQ.'ADCS') THEN 374 : C CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET) 375 : C IF (IRET.NE.0) GOTO 9999 375 : C CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET) 376 : C IF (IRET.NE.0) GOTO 9999 377 : C TR11=0.D0 378 : C DO J=1,IREF 379 : C TR11=TR11+ 380 : C $ (IG(IDIM2,J)*JAC(IDIM1,J)) 381 : C ENDDO 382 : C TR12=0.D0 383 : C DO I=1,IREF 384 : C TR12=TR12+ 385 : C $ (IG(I,IDIM2)*JAC(IDIM1,I)) 386 : C ENDDO 387 : C CONTRI=(TR11+TR12)*TRIGH*SDETG 388 : C* 389 : C* Calcul du jacobien (test) 390 : C* 391 : C ELSEIF (NOMLOI(1:4).EQ.'ADCK') THEN 391 : C CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET) 392 : C IF (IRET.NE.0) GOTO 9999 392 : C CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET) 393 : C IF (IRET.NE.0) GOTO 9999 393 : C CALL CH2INT(NOMLOI(7:7),IDIM3,IMPR,IRET) 394 : C IF (IRET.NE.0) GOTO 9999 394 : C CALL CH2INT(NOMLOI(8:8),IDIM4,IMPR,IRET) 395 : C IF (IRET.NE.0) GOTO 9999 396 : C* 397 : C TR10=0.D0 398 : C DO K=1,IREF 399 : C DO L=1,IREF 400 : C TR10=TR10+ 401 : C $ (JAC(IDIM1,K)*IG(K,L)*HIG(L,IDIM2)) 402 : C ENDDO 403 : C ENDDO 404 : C TR10=TR10*(-2.D0) 405 : C* 406 : C TR0=0.D0 407 : C DO K=1,IREF 408 : C TR0=TR0+ 409 : C $ (JAC(IDIM1,K)*IG(K,IDIM2)) 410 : C ENDDO 411 : C* 412 : C TR11=0.D0 413 : C DO J=1,IREF 414 : C TR11=TR11+ 415 : C $ (IG(IDIM4,J)*JAC(IDIM3,J)) 416 : C ENDDO 417 : C TR12=0.D0 418 : C DO I=1,IREF 419 : C TR12=TR12+ 420 : C $ (IG(I,IDIM4)*JAC(IDIM3,I)) 421 : C ENDDO 422 : C 423 : C* 424 : C* 425 : C TR13=0.D0 426 : C DO K=1,IREF 427 : C DO L=1,IREF 428 : C TR13=TR13+ 429 : C $ (JAC(IDIM1,K)*IG(K,L)*JAC(IDIM3,L) 430 : C $ *IG(IDIM4,IDIM2)) 431 : C ENDDO 432 : C ENDDO 433 : C TR13=TR13*(-1.D0) 434 : C TR14=0.D0 435 : C DO K=1,IREF 436 : C DO N=1,IREF 437 : C TR14=TR14+ 438 : C $ (JAC(IDIM1,K)*IG(K,IDIM4) 439 : C $ *JAC(IDIM3,N)*IG(N,IDIM2)) 440 : C ENDDO 441 : C ENDDO 442 : C TR14=TR14*(-1.D0) 443 : C* 444 : C TR15=0.D0 445 : C IF (IDIM1.EQ.IDIM3) THEN 446 : C TR15=TR15+IG(IDIM4,IDIM2) 447 : C ENDIF 448 : C* 449 : CC CONTRI=(TR13+TR14+TR15)*2.D0 450 : C CONTR3=(TR11+TR12)*(TR10+(TR0*TRIGH))*SDETG 451 : C CONTR4=(TR13+TR14+TR15)*2.D0*TRIGH*SDETG 452 : CC CONTR1=(TR1+TR2)*TR0 453 : CC CONTR2=(TR3+TR4+TR5+TR6+TR7)*2.D0 454 : CC CONTRI=(CONTR1+CONTR2)*SDETG 455 : C CONTRI=(CONTR3+CONTR4) 456 : C* 457 : C WRITE(IOIMP,*) 'CONTRI=',CONTRI 458 : ELSE 459 : WRITE(IOIMP,*) 'Erreur grave' 460 : GOTO 9999 461 : ENDIF 462 : * 463 : * Normal termination 464 : * 465 : IRET=0 466 : RETURN 467 : * 468 : * Format handling 469 : * 470 : * 471 : * Error handling 472 : * 473 : 9999 CONTINUE 474 : IRET=1 475 : WRITE(IOIMP,*) 'An error was detected in subroutine cgadxf' 476 : RETURN 477 : * 478 : * End of subroutine CGADXF 479 : * 480 : END 481 : 482 :
© Cast3M 2003 - All rights reserved.
Disclaimer