1 : C CCGTDI SOURCE GOUNAND 07/07/05 21:15:15 5784 2 : SUBROUTINE CCGTDI(LCOF, 3 : $ FC, 4 : $ IMPR,IRET) 5 : IMPLICIT REAL*8 (A-H,O-Z) 6 : IMPLICIT INTEGER (I-N) 7 : C*********************************************************************** 8 : C NOM : CCGTDI 9 : C DESCRIPTION : Calcul de la loi de comportement aux points de Gauss : 10 : C Diamètre interne d'un élément suivant une direction donnée 11 : C 12 : C LANGAGE : ESOPE 13 : C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF) 14 : C mél : gounand@semt2.smts.cea.fr 15 : C*********************************************************************** 16 : C APPELES : 17 : C APPELE PAR : 18 : C*********************************************************************** 19 : C ENTREES : 20 : C ENTREES/SORTIES : 21 : C SORTIES : - 22 : C TRAVAIL : 23 : C*********************************************************************** 24 : C VERSION : v1, 13/09/06, version initiale 25 : C HISTORIQUE : v1, 13/09/06, création 26 : C HISTORIQUE : 27 : C HISTORIQUE : 28 : C*********************************************************************** 29 : C Prière de PRENDRE LE TEMPS de compléter les commentaires 30 : C en cas de modification de ce sous-programme afin de faciliter 31 : C la maintenance ! 32 : C*********************************************************************** 33 : -INC CCOPTIO 34 : -INC CCREEL 35 : CBEGININCLUDE SMCHAEL 36 : SEGMENT MCHAEL 37 : POINTEUR IMACHE(N1).MELEME 38 : POINTEUR ICHEVA(N1).MCHEVA 39 : ENDSEGMENT 40 : SEGMENT MCHEVA 41 : REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM) 42 : ENDSEGMENT 43 : SEGMENT LCHEVA 44 : POINTEUR LISCHE(NBCHE).MCHEVA 45 : ENDSEGMENT 46 : CENDINCLUDE SMCHAEL 47 : INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM,N1 48 : POINTEUR FC.MCHEVA 49 : POINTEUR LCOF.LCHEVA 50 : POINTEUR JMAJAC.MCHEVA 51 : POINTEUR JMIJAC.MCHEVA 52 : POINTEUR JDTJAC.MCHEVA 53 : POINTEUR JMAREG.MCHEVA 54 : POINTEUR JDIAMA.MCHEVA 55 : POINTEUR JPC.MCHEVA 56 : POINTEUR JVIT.MCHEVA 57 : CHARACTER*8 NOMLOI 58 : INTEGER ICOF 59 : * 60 : -INC TMXMAT 61 : POINTEUR A.MXMAT 62 : POINTEUR AP.MXMAT 63 : POINTEUR JMA.MXMAT 64 : POINTEUR KJMA.MXMAT 65 : POINTEUR JM1.MXMAT,J.MXMAT 66 : POINTEUR K.MXMAT 67 : * 68 : SEGMENT MVIT 69 : POINTEUR MVCOMP(IDIM).MCHEVA 70 : ENDSEGMENT 71 : * 72 : LOGICAL LREGP,LRELP 73 : INTEGER IMPR,IRET 74 : * 75 : * Executable statements 76 : * 77 : IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans ccgtdi' 78 : * XPETI=SQRT(XPETIT) 79 : NLFC=FC.VELCHE(/6) 80 : NPFC=FC.VELCHE(/5) 81 : ICOF=0 82 : * 83 : * Récupération des coefficients de la metrique 84 : * 85 : SEGINI MVIT 86 : DO IIDIM=1,IDIM 87 : ICOF=ICOF+1 88 : JVIT=LCOF.LISCHE(ICOF) 89 : IF (ICOF.EQ.1) THEN 90 : NLJV=JVIT.VELCHE(/6) 91 : NPJV=JVIT.VELCHE(/5) 92 : ELSE 93 : NLJV2=JVIT.VELCHE(/6) 94 : NPJV2=JVIT.VELCHE(/5) 95 : IF (NLJV2.NE.NLJV.OR.NPJV2.NE.NPJV) THEN 96 : WRITE(IOIMP,*) 'Erreur grave dims JVIT' 97 : GOTO 9999 98 : ENDIF 99 : ENDIF 100 : MVIT.MVCOMP(IIDIM)=JVIT 101 : ENDDO 102 : * 103 : ICOF=ICOF+1 104 : JMAJAC=LCOF.LISCHE(ICOF) 105 : C NLJA=JMAJAC.VELCHE(/6) 106 : C NPJA=JMAJAC.VELCHE(/5) 107 : IREF=JMAJAC.VELCHE(/4) 108 : IREL=JMAJAC.VELCHE(/3) 109 : * 110 : IF (IREL.NE.IDIM) THEN 111 : WRITE(IOIMP,*) 'Erreur dims JMAJAC' 112 : GOTO 9999 113 : ENDIF 114 : * 115 : ICOF=ICOF+1 116 : JMIJAC=LCOF.LISCHE(ICOF) 117 : IF (JMIJAC.EQ.0) THEN 118 : WRITE(IOIMP,*) 'Erreur JMIJAC=0' 119 : GOTO 9999 120 : ENDIF 121 : NLJI=JMIJAC.VELCHE(/6) 122 : NPJI=JMIJAC.VELCHE(/5) 123 : IREL2=JMIJAC.VELCHE(/4) 124 : IREF2=JMIJAC.VELCHE(/3) 125 : * 126 : IF (IREL2.NE.IREL.OR.IREF2.NE.IREF) THEN 127 : WRITE(IOIMP,*) 'Erreur dims JMIJAC' 128 : GOTO 9999 129 : ENDIF 130 : * 131 : ICOF=ICOF+1 132 : ICOF=ICOF+1 133 : JMAREG=LCOF.LISCHE(ICOF) 134 : NLJR=JMAREG.VELCHE(/6) 135 : NPJR=JMAREG.VELCHE(/5) 136 : I1 =JMAREG.VELCHE(/4) 137 : I2 =JMAREG.VELCHE(/3) 138 : IF ((NLJR.NE.1).OR.(NPJR.NE.1).OR.(I1.NE.IREF).OR.(I2.NE.IREF)) 139 : $ THEN 140 : WRITE(IOIMP,*) 'Erreur dims JMAREG' 141 : GOTO 9999 142 : ENDIF 143 : ICOF=ICOF+1 144 : JDIAMA=LCOF.LISCHE(ICOF) 145 : NLJD=JDIAMA.VELCHE(/6) 146 : NPJD=JDIAMA.VELCHE(/5) 147 : I1 =JDIAMA.VELCHE(/4) 148 : I2 =JDIAMA.VELCHE(/3) 149 : IF ((NLJD.NE.1).OR.(NPJD.NE.1).OR.(I1.NE.1).OR.(I2.NE.1)) 150 : $ THEN 151 : WRITE(IOIMP,*) 'Erreur dims JDIAMA' 152 : GOTO 9999 153 : ENDIF 154 : XDIAMA=JDIAMA.VELCHE(1,1,1,1,1,1) 155 : * 156 : * Objets temporaires 157 : * 158 : LDIM1=IREL 159 : LDIM2=1 160 : SEGINI,A 161 : SEGINI,AP 162 : LDIM1=IREF 163 : LDIM2=1 164 : SEGINI,JMA 165 : LDIM1=IREF 166 : LDIM2=1 167 : SEGINI,KJMA 168 : LDIM1=IREL 169 : LDIM2=IREF 170 : SEGINI,J 171 : LDIM1=IREF 172 : LDIM2=IREL 173 : SEGINI,JM1 174 : LDIM1=IREF 175 : LDIM2=IREF 176 : SEGINI,K 177 : * Copie de la matrice jacobienne ref -> regulier 177 : CALL MAMA(JMAREG.VELCHE,IREF,IREF, 178 : $ 'COPIE ',K.XMAT,IREF,IREF, 179 : $ IMPR,IRET) 180 : IF (IRET.NE.0) GOTO 9999 181 : * SEGPRT,H 182 : DO ILFC=1,NLFC 183 : IF (NLJV.EQ.1) THEN 184 : ILJV=1 185 : ELSE 186 : ILJV=ILFC 187 : ENDIF 188 : C IF (NLJA.EQ.1) THEN 189 : C ILJA=1 190 : C ELSE 191 : C ILJA=ILFC 192 : C ENDIF 193 : IF (NLJI.EQ.1) THEN 194 : ILJI=1 195 : ELSE 196 : ILJI=ILFC 197 : ENDIF 198 : * 199 : DO IPFC=1,NPFC 200 : IF (NPJV.EQ.1) THEN 201 : IPJV=1 202 : ELSE 203 : IPJV=IPFC 204 : ENDIF 205 : C IF (NPJA.EQ.1) THEN 206 : C IPJA=1 207 : C ELSE 208 : C IPJA=IPFC 209 : C ENDIF 210 : IF (NPJI.EQ.1) THEN 211 : IPJI=1 212 : ELSE 213 : IPJI=IPFC 214 : ENDIF 215 : * 216 : * Copie de la vitesse 217 : * 218 : DO IIDIM=1,IDIM 219 : JVIT=MVIT.MVCOMP(IIDIM) 220 : A.XMAT(IIDIM,1)=JVIT.VELCHE(1,1,1,1,IPJV,ILJV) 221 : ENDDO 222 : * SEGPRT,A 223 : * 224 : * Copie de l'inverse (ou pseudo-inverse) du jacobien ref->reel 225 : * 225 : CALL MAMA(JMIJAC.VELCHE(1,1,1,1,IPJI,ILJI),IREF,IREL, 226 : $ 'COPIE ', 227 : $ JM1.XMAT,IREF,IREL, 228 : $ IMPR,IRET) 229 : IF (IRET.NE.0) GOTO 9999 230 : C SEGPRT,JM1 231 : * 232 : * Copie du jacobien ref->reel 233 : * 233 : CALL MAMA(JMAJAC.VELCHE(1,1,1,1,IPJI,ILJI),IREL,IREF, 234 : $ 'COPIE ', 235 : $ J.XMAT,IREL,IREF, 236 : $ IMPR,IRET) 237 : IF (IRET.NE.0) GOTO 9999 238 : C SEGPRT,JM1 239 : * 240 : * Vecteur A dans le repère de référence 241 : * 241 : CALL MAMAMA(JM1.XMAT,IREF,IREL,A.XMAT,IREL,1, 242 : $ 'FOIS ',JMA.XMAT,IREF,1, 243 : $ IMPR,IRET) 244 : IF (IRET.NE.0) GOTO 9999 245 : * 246 : * Vecteur A' dans le repère réel (projection sur le domaine considéré) 247 : * 247 : CALL MAMAMA(J.XMAT,IREL,IREF,JMA.XMAT,IREF,1, 248 : $ 'FOIS ',AP.XMAT,IREL,1, 249 : $ IMPR,IRET) 250 : IF (IRET.NE.0) GOTO 9999 251 : * 252 : * Vecteur A dans le repère de l'élément régulier 253 : * 253 : CALL MAMAMA(K.XMAT,IREF,IREF,JMA.XMAT,IREF,1, 254 : $ 'FOIS ',KJMA.XMAT,IREF,1, 255 : $ IMPR,IRET) 256 : IF (IRET.NE.0) GOTO 9999 257 : * 258 : * Normes de A 259 : * 260 : XAREG2=XZERO 261 : DO IIREF=1,IREF 262 : XAREG2=XAREG2+(KJMA.XMAT(IIREF,1)**2) 263 : ENDDO 264 : XA2=XZERO 265 : DO IIDIM=1,IDIM 266 : XA2=XA2+(A.XMAT(IIDIM,1)**2) 267 : ENDDO 268 : XAP2=XZERO 269 : DO IIDIM=1,IDIM 270 : XAP2=XAP2+(AP.XMAT(IIDIM,1)**2) 271 : ENDDO 272 : XAPA=XZERO 273 : DO IIDIM=1,IDIM 274 : XAPA=XAPA+(AP.XMAT(IIDIM,1)*A.XMAT(IIDIM,1)) 275 : ENDDO 276 : C WRITE(IOIMP,*) 'XA2=',XA2 277 : C WRITE(IOIMP,*) 'XAP2=',XAP2 278 : C WRITE(IOIMP,*) 'XAREG2=',XAREG2 279 : C WRITE(IOIMP,*) 'XDIAMA=',XDIAMA 280 : * 281 : IF (XA2.LT.XPETIT) THEN 282 : WRITE(IOIMP,*) 'The given direction is 0' 283 : GOTO 9999 284 : ENDIF 285 : IF (XAREG2.LT.XPETIT) THEN 286 : CONTRI=0.D0 287 : ELSE 288 : CONTRI=SQRT((XAP2*XAPA)/(XAREG2*XA2))*XDIAMA 289 : * CONTRI=SQRT(XA2/XAREG2)*XDIAMA 290 : ENDIF 291 : FC.VELCHE(1,1,1,1,IPFC,ILFC)= 292 : $ FC.VELCHE(1,1,1,1,IPFC,ILFC)+ 293 : $ CONTRI 294 : ENDDO 295 : ENDDO 296 : SEGSUP,K 297 : SEGSUP,J 298 : SEGSUP,JM1 299 : SEGSUP,KJMA 300 : SEGSUP,JMA 301 : SEGSUP,AP 302 : SEGSUP,A 303 : SEGSUP,MVIT 304 : * 305 : * Normal termination 306 : * 307 : IRET=0 308 : RETURN 309 : * 310 : * Format handling 311 : * 312 : * 313 : * Error handling 314 : * 315 : 9999 CONTINUE 316 : IRET=1 317 : WRITE(IOIMP,*) 'An error was detected in subroutine ccgtdi' 318 : RETURN 319 : * 320 : * End of subroutine CCGTDI 321 : * 322 : END 323 : 324 : 325 : 326 : 327 :
© Cast3M 2003 - All rights reserved.
Disclaimer