C TESTFA SOURCE CHAT 05/01/13 03:36:24 5004 SUBROUTINE TESTFA(EPSILO,NDIM,ITYG $ ,XREEL,IZSH,JFACE,PT1,PT2,PT3,PT4) ************************************************************************* *** SP 'TESTFA' : permet de tester % element considéré appartenance pt à *** l'une des faces associées, renvoie alors n° local face ou se trouve pt. *** *** APPELES 1 = aucun *** APPELES 2 = 'LIEUP2' *** *** E = 'EPSILO' marge relative acceptée position pt % element *** 'NDIM' dimension de l'espace *** 'ITYG' entier caracterisant la geometrie de l'element *** 'XREEL' coordonnees reelles du pt considéré *** 'IZSH' segment content coords reelles des noeuds de l'elemt considere *** *** S = 'JFACE' n° local face ou se trouve particule, -1 sinon *** 'PT1', 'PT2', 'PT3', 'PT4' noeuds appartenant à la face *** *** Rq : developpee seulement pour ITYPEL = 4,8,14,16 et 23 *** *** *** Auteur Cyril Nou ************************************************************************* IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) SEGMENT IZSH REAL*8 SHP(6,MNO9),SHY(12,MNO9),XYZL(3,MNO9) ENDSEGMENT DIMENSION XREEL(3) DIMENSION PT1(3),PT2(3),PT3(3),PT4(3) DIMENSION ITRIAN(5),ICARRE(6),ICUBE(4,6),IPRISM(4,5),ITETRA(4,4) *** donnees ordonnées specifiquement pour chaque type d'element afin de *** parcourir les faces dans l'ordre croissant dans les diverses boucles DATA ITRIAN/1,2,3,1,2/ DATA ICARRE/1,2,3,4,1,2/ DATA ICUBE/1,2,3,5, 5,6,7,1, 1,2,6,3, 2,3,7,1, 3,4,8,1, 1,4,8,2/ DATA IPRISM/1,2,3,4, 4,5,6,1, 1,2,5,3, 2,3,6,1, 1,3,6,2/ DATA ITETRA/1,2,4,3, 1,2,3,4, 2,3,4,1, 1,3,4,2/ *** on postule initialement non appartenance à face JFACE=-1 *** cas TRI3 (triangles), boucle sur les différentes faces IF (ITYG.EQ.4) THEN DO 10 I=1,3 *** recuperation des pts definissant la face DO 20 J=1,NDIM PT1(J)=XYZL(J,ITRIAN(I)) PT2(J)=XYZL(J,ITRIAN(I+1)) PT3(J)=XYZL(J,ITRIAN(I+2)) 20 CONTINUE CALL LIEUP2(NDIM,PT1,PT2,PT3,PT4,XREEL,TEST) *** 'TEST'=0. si appartenance à face considérée IF (ABS(TEST).LE.EPSILO) THEN JFACE=I RETURN ENDIF 10 CONTINUE *** cas QUA4 (quadrangles) ELSEIF (ITYG.EQ.8) THEN DO 30 I=1,4 DO 40 J=1,NDIM PT1(J)=XYZL(J,ICARRE(I)) PT2(J)=XYZL(J,ICARRE(I+1)) PT3(J)=XYZL(J,ICARRE(I+2)) 40 CONTINUE CALL LIEUP2(NDIM,PT1,PT2,PT3,PT4,XREEL,TEST) IF (ABS(TEST).LE.EPSILO) THEN JFACE=I RETURN ENDIF 30 CONTINUE *** cas CUB8 (cubes) ELSEIF (ITYG.EQ.14) THEN DO 50 I=1,6 DO 60 J=1,NDIM PT1(J)=XYZL(J,ICUBE(1,I)) PT2(J)=XYZL(J,ICUBE(2,I)) PT3(J)=XYZL(J,ICUBE(3,I)) PT4(J)=XYZL(J,ICUBE(4,I)) 60 CONTINUE CALL LIEUP2(NDIM,PT1,PT2,PT3,PT4,XREEL,TEST) IF (ABS(TEST).LE.EPSILO) THEN JFACE=I RETURN ENDIF 50 CONTINUE *** cas PRI6 (prismes) ELSEIF (ITYG.EQ.16) THEN DO 70 I=1,5 DO 80 J=1,NDIM PT1(J)=XYZL(J,IPRISM(1,I)) PT2(J)=XYZL(J,IPRISM(2,I)) PT3(J)=XYZL(J,IPRISM(3,I)) PT4(J)=XYZL(J,IPRISM(4,I)) 80 CONTINUE CALL LIEUP2(NDIM,PT1,PT2,PT3,PT4,XREEL,TEST) IF (ABS(TEST).LE.EPSILO) THEN JFACE=I RETURN ENDIF 70 CONTINUE *** cas TET4 (tetraedres) ELSEIF (ITYG.EQ.23) THEN DO 90 I=1,4 DO 100 J=1,NDIM PT1(J)=XYZL(J,ITETRA(1,I)) PT2(J)=XYZL(J,ITETRA(2,I)) PT3(J)=XYZL(J,ITETRA(3,I)) PT4(J)=XYZL(J,ITETRA(4,I)) 100 CONTINUE CALL LIEUP2(NDIM,PT1,PT2,PT3,PT4,XREEL,TEST) IF (ABS(TEST).LE.EPSILO) THEN JFACE=I RETURN ENDIF 90 CONTINUE ENDIF RETURN END