1 : C TESJA1 SOURCE GOUNAND 05/12/21 21:37:04 5281 2 : SUBROUTINE TESJA1(MYLRFS,MYPGS, 3 : $ IMPR,IRET) 4 : IMPLICIT REAL*8 (A-H,O-Z) 5 : IMPLICIT INTEGER (I-N) 6 : C*********************************************************************** 7 : C NOM : TESJA1 8 : C PROJET : Noyau linéaire NLIN 9 : C DESCRIPTION : On vérifie le calcul de la matrice jacobienne et de son 10 : C déterminant sur un carré... 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 APPELES (E/S) : 18 : C APPELES (BLAS) : 19 : C APPELES (CALCUL) : 20 : C APPELE PAR : TESTJA 21 : C*********************************************************************** 22 : C SYNTAXE GIBIANE : 23 : C ENTREES : 24 : C ENTREES/SORTIES : 25 : C SORTIES : 26 : C*********************************************************************** 27 : C VERSION : v1, 16/08/99, version initiale 28 : C HISTORIQUE : v1, 16/08/99, création 29 : C HISTORIQUE : 30 : C HISTORIQUE : 31 : C*********************************************************************** 32 : C Prière de PRENDRE LE TEMPS de compléter les commentaires 33 : C en cas de modification de ce sous-programme afin de faciliter 34 : C la maintenance ! 35 : C*********************************************************************** 36 : -INC CCOPTIO 37 : -INC CCREEL 38 : CBEGININCLUDE SELREF 39 : SEGMENT ELREF 40 : CHARACTER*(LNNOM) NOMLRF 41 : CHARACTER*(LNFORM) FORME 42 : CHARACTER*(LNTYPL) TYPEL 43 : CHARACTER*(LNESP) ESPACE 44 : INTEGER DEGRE 45 : REAL*8 XCONOD(NDIMEL,NBNOD) 46 : INTEGER NPQUAF(NBDDL) 47 : INTEGER NUMCMP(NBDDL) 48 : INTEGER QUENOD(NBDDL) 49 : INTEGER ORDDER(NDIMEL,NBDDL) 50 : POINTEUR MBPOLY.POLYNS 51 : ENDSEGMENT 52 : SEGMENT ELREFS 53 : POINTEUR LISEL(0).ELREF 54 : ENDSEGMENT 55 : CENDINCLUDE SELREF 56 : POINTEUR MYLRFS.ELREFS 57 : POINTEUR MYLRF.ELREF 58 : CBEGININCLUDE SPOGAU 59 : SEGMENT POGAU 60 : CHARACTER*(LNNPG) NOMPG 61 : CHARACTER*(LNTPG) TYPMPG 62 : CHARACTER*(LNFPG) FORLPG 63 : INTEGER NORDPG 64 : REAL*8 XCOPG(NDLPG,NBPG) 65 : REAL*8 XPOPG(NBPG) 66 : ENDSEGMENT 67 : SEGMENT POGAUS 68 : POINTEUR LISPG(0).POGAU 69 : ENDSEGMENT 70 : CENDINCLUDE SPOGAU 71 : POINTEUR MYPGS.POGAUS 72 : POINTEUR MYPG.POGAU 73 : CBEGININCLUDE SMCHAEL 74 : SEGMENT MCHAEL 75 : POINTEUR IMACHE(N1).MELEME 76 : POINTEUR ICHEVA(N1).MCHEVA 77 : ENDSEGMENT 78 : SEGMENT MCHEVA 79 : REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM) 80 : ENDSEGMENT 81 : SEGMENT LCHEVA 82 : POINTEUR LISCHE(NBCHE).MCHEVA 83 : ENDSEGMENT 84 : CENDINCLUDE SMCHAEL 85 : INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM 86 : POINTEUR JCOOR.MCHEVA 87 : POINTEUR JMAJAC.MCHEVA 88 : POINTEUR JMIJAC.MCHEVA 89 : POINTEUR JDTJAC.MCHEVA 90 : POINTEUR FFGPG.MCHEVA 91 : POINTEUR DFFGPG.MCHEVA 92 : * 93 : INTEGER IMPR,IRET 94 : * 95 : REAL*8 JXX(2,2) 96 : REAL*8 X1,X2,X3,X4 97 : REAL*8 Y1,Y2,Y3,Y4 98 : REAL*8 KSIPG,ETAPG,A0,A1,A2,DETJXX,DETA,JXXA,XERR 99 : REAL*8 XERREL 100 : INTEGER NPG,INPG,ICOL,ILIG,INBTES 101 : LOGICAL LERJ 102 : * 103 : * Executable statements 104 : * 105 : WRITE(IOIMP,*) 'Entrée dans tesja1' 106 : * 107 : * On teste sur un élément quadrilatéral à quatre noeuds (D&T p.54-55) 108 : * avec douze points de Gauss, dans le cas rectangulaire, et dans le cas 109 : * où l'élément est illicite (les côtés opposés se croisent). 110 : * Dimension des espaces de référence et réels : 2, 2 111 : * 111 : CALL FILRF('H1D1QU4',MYLRFS,MYLRF,IMPR,IRET) 112 : IF (IRET.NE.0) GOTO 9999 112 : CALL FIPG('GAC2-7-12A',MYPGS,MYPG,IMPR,IRET) 113 : IF (IRET.NE.0) GOTO 9999 113 : CALL KFNREF(MYLRF,MYPG,FFGPG,DFFGPG,IMPR,IRET) 114 : IF (IRET.NE.0) GOTO 9999 115 : DO 3 INBTES=1,3 116 : IF (INBTES.EQ.1) THEN 117 : X1=0.9D0 118 : Y1=1.1D0 119 : X2=1.75D0 120 : Y2=2.25D0 121 : X3=3.4D0 122 : Y3=3.9D0 123 : X4=2.1D0 124 : Y4=6.2D0 125 : ELSEIF (INBTES.EQ.2) THEN 126 : X1=0.9D0 127 : Y1=1.1D0 128 : X2=1.75D0 129 : Y2=Y1 130 : X3=X2 131 : Y3=3.9D0 132 : X4=X1 133 : Y4=Y3 134 : ELSEIF (INBTES.EQ.3) THEN 135 : X1=0.9D0 136 : Y1=1.1D0 137 : X3=1.75D0 138 : Y3=2.25D0 139 : X2=3.4D0 140 : Y2=3.9D0 141 : X4=2.1D0 142 : Y4=6.2D0 143 : ELSE 144 : WRITE(IOIMP,*) 'Erreur dans le nombre de tests' 145 : GOTO 9999 146 : ENDIF 147 : NBELM=1 148 : NBPOI=1 149 : N2COL=2 150 : N2LIG=1 151 : NBCOL=4 152 : NBLIG=1 153 : SEGINI JCOOR 154 : JCOOR.VELCHE(1,1,1,1,1,1)=X1 155 : JCOOR.VELCHE(1,1,1,2,1,1)=Y1 156 : JCOOR.VELCHE(1,2,1,1,1,1)=X2 157 : JCOOR.VELCHE(1,2,1,2,1,1)=Y2 158 : JCOOR.VELCHE(1,3,1,1,1,1)=X3 159 : JCOOR.VELCHE(1,3,1,2,1,1)=Y3 160 : JCOOR.VELCHE(1,4,1,1,1,1)=X4 161 : JCOOR.VELCHE(1,4,1,2,1,1)=Y4 162 : SEGDES JCOOR 163 : IF (IMPR.GT.3) THEN 164 : WRITE(IOIMP,*) 'JCOOR' 164 : CALL PRCHVA(JCOOR,IMPR,IRET) 165 : IF (IRET.NE.0) GOTO 9999 166 : ENDIF 167 : LERJ=.FALSE. 167 : CALL GEOLIN(DFFGPG,JCOOR,NBELM, 168 : $ JMAJAC,JMIJAC,JDTJAC,LERJ, 169 : $ IMPR,IRET) 170 : IF (INBTES.NE.3) THEN 171 : IF (IRET.NE.0) GOTO 9999 172 : * 173 : * Test sur les valeurs de la matrice jacobienne et de son 174 : * déterminant aux points de Gauss... 175 : * 176 : SEGACT MYPG 177 : SEGACT JMAJAC 178 : SEGACT JDTJAC 179 : NPG=MYPG.XCOPG(/2) 180 : DO 32 INPG=1,NPG 181 : KSIPG=MYPG.XCOPG(1,INPG) 182 : ETAPG=MYPG.XCOPG(2,INPG) 183 : JXX(1,1)=0.25D0*((-X1+X2+X3-X4) 184 : $ +(ETAPG*(X1-X2+X3-X4))) 185 : JXX(1,2)=0.25D0*((-X1-X2+X3+X4) 186 : $ +(KSIPG*(X1-X2+X3-X4))) 187 : JXX(2,1)=0.25D0*((-Y1+Y2+Y3-Y4) 188 : $ +(ETAPG*(Y1-Y2+Y3-Y4))) 189 : JXX(2,2)=0.25D0*((-Y1-Y2+Y3+Y4) 190 : $ +(KSIPG*(Y1-Y2+Y3-Y4))) 191 : A0=0.125D0*(((Y4-Y2)*(X3-X1)) 192 : $ -((Y3-Y1)*(X4-X2))) 193 : A1=0.125D0*(((Y3-Y4)*(X2-X1)) 194 : $ -((Y2-Y1)*(X3-X4))) 195 : A2=0.125D0*(((Y4-Y1)*(X3-X2)) 196 : $ -((Y3-Y2)*(X4-X1))) 197 : DETJXX=A0+(A1*KSIPG)+(A2*ETAPG) 198 : DO 322 ILIG=1,2 199 : DO 3222 ICOL=1,2 200 : JXXA=JMAJAC.VELCHE(1,1,ILIG,ICOL,INPG,1) 201 : XERR=XERREL(JXX(ILIG,ICOL),JXXA) 202 : IF (XERR.GT.10.D0*XZPREC) THEN 203 : WRITE(IOIMP,*) 'Erreur calcul mat. jac.' 204 : WRITE(IOIMP,*) 'XERR=',XERR, ' XZPREC=',XZPREC 205 : WRITE(IOIMP,*) 'JXX=',JXX(ILIG,ICOL) 206 : WRITE(IOIMP,*) 'JXXA=',JXXA 207 : GOTO 9999 208 : ENDIF 209 : 3222 CONTINUE 210 : 322 CONTINUE 211 : DETA=JDTJAC.VELCHE(1,1,1,1,INPG,1) 212 : XERR=XERREL(DETJXX,DETA) 213 : IF (XERR.GT.10.D0*XZPREC) THEN 214 : WRITE(IOIMP,*) 'Erreur calcul det. mat. jac.' 215 : WRITE(IOIMP,*) 'XERR=',XERR, ' XZPREC=',XZPREC 216 : WRITE(IOIMP,*) 'DETJXX=',DETJXX 217 : WRITE(IOIMP,*) 'DETA=',DETA 218 : GOTO 9999 219 : ENDIF 220 : 32 CONTINUE 221 : SEGDES MYPG 222 : SEGDES JDTJAC 223 : SEGDES JMAJAC 224 : IF (IMPR.GT.3) THEN 225 : WRITE(IOIMP,*) 'JMAJAC' 225 : CALL PRCHVA(JMAJAC,IMPR,IRET) 226 : IF (IRET.NE.0) GOTO 9999 227 : WRITE(IOIMP,*) 'JMIJAC' 227 : CALL PRCHVA(JMIJAC,IMPR,IRET) 228 : IF (IRET.NE.0) GOTO 9999 229 : WRITE(IOIMP,*) 'JDTJAC' 229 : CALL PRCHVA(JDTJAC,IMPR,IRET) 230 : IF (IRET.NE.0) GOTO 9999 231 : ENDIF 232 : WRITE(IOIMP,*) 'Test',INBTES,' successful' 233 : ELSE 234 : IF (IRET.EQ.0) THEN 235 : WRITE(IOIMP,*) 'Le test aurait dû planter...' 236 : GOTO 9999 237 : ELSE 238 : WRITE(IOIMP,*) 'Le test a planté comme prévu...' 239 : ENDIF 240 : ENDIF 241 : 3 CONTINUE 242 : * 243 : * Normal termination 244 : * 245 : IRET=0 246 : RETURN 247 : * 248 : * Format handling 249 : * 250 : * 251 : * Error handling 252 : * 253 : 9999 CONTINUE 254 : IRET=1 255 : WRITE(IOIMP,*) 'An error was detected in subroutine tesja1' 256 : RETURN 257 : * 258 : * End of subroutine tesja1 259 : * 260 : END 261 : 262 :
© Cast3M 2003 - All rights reserved.
Disclaimer