C TESJA1 SOURCE GOUNAND 05/12/21 21:37:04 5281 SUBROUTINE TESJA1(MYLRFS,MYPGS, $ IMPR,IRET) IMPLICIT REAL*8 (A-H,O-Z) IMPLICIT INTEGER (I-N) C*********************************************************************** C NOM : TESJA1 C PROJET : Noyau linéaire NLIN C DESCRIPTION : On vérifie le calcul de la matrice jacobienne et de son C déterminant sur un carré... C C LANGAGE : ESOPE C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF) C mél : gounand@semt2.smts.cea.fr C*********************************************************************** C APPELES : C APPELES (E/S) : C APPELES (BLAS) : C APPELES (CALCUL) : C APPELE PAR : TESTJA C*********************************************************************** C SYNTAXE GIBIANE : C ENTREES : C ENTREES/SORTIES : C SORTIES : C*********************************************************************** C VERSION : v1, 16/08/99, version initiale C HISTORIQUE : v1, 16/08/99, création C HISTORIQUE : C HISTORIQUE : C*********************************************************************** C Prière de PRENDRE LE TEMPS de compléter les commentaires C en cas de modification de ce sous-programme afin de faciliter C la maintenance ! C*********************************************************************** -INC CCOPTIO -INC CCREEL CBEGININCLUDE SELREF SEGMENT ELREF CHARACTER*(LNNOM) NOMLRF CHARACTER*(LNFORM) FORME CHARACTER*(LNTYPL) TYPEL CHARACTER*(LNESP) ESPACE INTEGER DEGRE REAL*8 XCONOD(NDIMEL,NBNOD) INTEGER NPQUAF(NBDDL) INTEGER NUMCMP(NBDDL) INTEGER QUENOD(NBDDL) INTEGER ORDDER(NDIMEL,NBDDL) POINTEUR MBPOLY.POLYNS ENDSEGMENT SEGMENT ELREFS POINTEUR LISEL(0).ELREF ENDSEGMENT CENDINCLUDE SELREF POINTEUR MYLRFS.ELREFS POINTEUR MYLRF.ELREF CBEGININCLUDE SPOGAU SEGMENT POGAU CHARACTER*(LNNPG) NOMPG CHARACTER*(LNTPG) TYPMPG CHARACTER*(LNFPG) FORLPG INTEGER NORDPG REAL*8 XCOPG(NDLPG,NBPG) REAL*8 XPOPG(NBPG) ENDSEGMENT SEGMENT POGAUS POINTEUR LISPG(0).POGAU ENDSEGMENT CENDINCLUDE SPOGAU POINTEUR MYPGS.POGAUS POINTEUR MYPG.POGAU CBEGININCLUDE SMCHAEL SEGMENT MCHAEL POINTEUR IMACHE(N1).MELEME POINTEUR ICHEVA(N1).MCHEVA ENDSEGMENT SEGMENT MCHEVA REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM) ENDSEGMENT SEGMENT LCHEVA POINTEUR LISCHE(NBCHE).MCHEVA ENDSEGMENT CENDINCLUDE SMCHAEL INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM POINTEUR JCOOR.MCHEVA POINTEUR JMAJAC.MCHEVA POINTEUR JMIJAC.MCHEVA POINTEUR JDTJAC.MCHEVA POINTEUR FFGPG.MCHEVA POINTEUR DFFGPG.MCHEVA * INTEGER IMPR,IRET * REAL*8 JXX(2,2) REAL*8 X1,X2,X3,X4 REAL*8 Y1,Y2,Y3,Y4 REAL*8 KSIPG,ETAPG,A0,A1,A2,DETJXX,DETA,JXXA,XERR REAL*8 XERREL INTEGER NPG,INPG,ICOL,ILIG,INBTES LOGICAL LERJ * * Executable statements * WRITE(IOIMP,*) 'Entrée dans tesja1' * * On teste sur un élément quadrilatéral à quatre noeuds (D&T p.54-55) * avec douze points de Gauss, dans le cas rectangulaire, et dans le cas * où l'élément est illicite (les côtés opposés se croisent). * Dimension des espaces de référence et réels : 2, 2 * CALL FILRF('H1D1QU4',MYLRFS,MYLRF,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 CALL FIPG('GAC2-7-12A',MYPGS,MYPG,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 CALL KFNREF(MYLRF,MYPG,FFGPG,DFFGPG,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 DO 3 INBTES=1,3 IF (INBTES.EQ.1) THEN X1=0.9D0 Y1=1.1D0 X2=1.75D0 Y2=2.25D0 X3=3.4D0 Y3=3.9D0 X4=2.1D0 Y4=6.2D0 ELSEIF (INBTES.EQ.2) THEN X1=0.9D0 Y1=1.1D0 X2=1.75D0 Y2=Y1 X3=X2 Y3=3.9D0 X4=X1 Y4=Y3 ELSEIF (INBTES.EQ.3) THEN X1=0.9D0 Y1=1.1D0 X3=1.75D0 Y3=2.25D0 X2=3.4D0 Y2=3.9D0 X4=2.1D0 Y4=6.2D0 ELSE WRITE(IOIMP,*) 'Erreur dans le nombre de tests' GOTO 9999 ENDIF NBELM=1 NBPOI=1 N2COL=2 N2LIG=1 NBCOL=4 NBLIG=1 SEGINI JCOOR JCOOR.VELCHE(1,1,1,1,1,1)=X1 JCOOR.VELCHE(1,1,1,2,1,1)=Y1 JCOOR.VELCHE(1,2,1,1,1,1)=X2 JCOOR.VELCHE(1,2,1,2,1,1)=Y2 JCOOR.VELCHE(1,3,1,1,1,1)=X3 JCOOR.VELCHE(1,3,1,2,1,1)=Y3 JCOOR.VELCHE(1,4,1,1,1,1)=X4 JCOOR.VELCHE(1,4,1,2,1,1)=Y4 SEGDES JCOOR IF (IMPR.GT.3) THEN WRITE(IOIMP,*) 'JCOOR' CALL PRCHVA(JCOOR,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 ENDIF LERJ=.FALSE. CALL GEOLIN(DFFGPG,JCOOR,NBELM, $ JMAJAC,JMIJAC,JDTJAC,LERJ, $ IMPR,IRET) IF (INBTES.NE.3) THEN IF (IRET.NE.0) GOTO 9999 * * Test sur les valeurs de la matrice jacobienne et de son * déterminant aux points de Gauss... * SEGACT MYPG SEGACT JMAJAC SEGACT JDTJAC NPG=MYPG.XCOPG(/2) DO 32 INPG=1,NPG KSIPG=MYPG.XCOPG(1,INPG) ETAPG=MYPG.XCOPG(2,INPG) JXX(1,1)=0.25D0*((-X1+X2+X3-X4) $ +(ETAPG*(X1-X2+X3-X4))) JXX(1,2)=0.25D0*((-X1-X2+X3+X4) $ +(KSIPG*(X1-X2+X3-X4))) JXX(2,1)=0.25D0*((-Y1+Y2+Y3-Y4) $ +(ETAPG*(Y1-Y2+Y3-Y4))) JXX(2,2)=0.25D0*((-Y1-Y2+Y3+Y4) $ +(KSIPG*(Y1-Y2+Y3-Y4))) A0=0.125D0*(((Y4-Y2)*(X3-X1)) $ -((Y3-Y1)*(X4-X2))) A1=0.125D0*(((Y3-Y4)*(X2-X1)) $ -((Y2-Y1)*(X3-X4))) A2=0.125D0*(((Y4-Y1)*(X3-X2)) $ -((Y3-Y2)*(X4-X1))) DETJXX=A0+(A1*KSIPG)+(A2*ETAPG) DO 322 ILIG=1,2 DO 3222 ICOL=1,2 JXXA=JMAJAC.VELCHE(1,1,ILIG,ICOL,INPG,1) XERR=XERREL(JXX(ILIG,ICOL),JXXA) IF (XERR.GT.10.D0*XZPREC) THEN WRITE(IOIMP,*) 'Erreur calcul mat. jac.' WRITE(IOIMP,*) 'XERR=',XERR, ' XZPREC=',XZPREC WRITE(IOIMP,*) 'JXX=',JXX(ILIG,ICOL) WRITE(IOIMP,*) 'JXXA=',JXXA GOTO 9999 ENDIF 3222 CONTINUE 322 CONTINUE DETA=JDTJAC.VELCHE(1,1,1,1,INPG,1) XERR=XERREL(DETJXX,DETA) IF (XERR.GT.10.D0*XZPREC) THEN WRITE(IOIMP,*) 'Erreur calcul det. mat. jac.' WRITE(IOIMP,*) 'XERR=',XERR, ' XZPREC=',XZPREC WRITE(IOIMP,*) 'DETJXX=',DETJXX WRITE(IOIMP,*) 'DETA=',DETA GOTO 9999 ENDIF 32 CONTINUE SEGDES MYPG SEGDES JDTJAC SEGDES JMAJAC IF (IMPR.GT.3) THEN WRITE(IOIMP,*) 'JMAJAC' CALL PRCHVA(JMAJAC,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 WRITE(IOIMP,*) 'JMIJAC' CALL PRCHVA(JMIJAC,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 WRITE(IOIMP,*) 'JDTJAC' CALL PRCHVA(JDTJAC,IMPR,IRET) IF (IRET.NE.0) GOTO 9999 ENDIF WRITE(IOIMP,*) 'Test',INBTES,' successful' ELSE IF (IRET.EQ.0) THEN WRITE(IOIMP,*) 'Le test aurait dû planter...' GOTO 9999 ELSE WRITE(IOIMP,*) 'Le test a planté comme prévu...' ENDIF ENDIF 3 CONTINUE * * Normal termination * IRET=0 RETURN * * Format handling * * * Error handling * 9999 CONTINUE IRET=1 WRITE(IOIMP,*) 'An error was detected in subroutine tesja1' RETURN * * End of subroutine tesja1 * END