Sources Esope | Cast3M

jo4loc
C JO4LOC    SOURCE    CB215821  17/11/30    21:16:30     9639           
      SUBROUTINE JO4LOC(XE,SHPTOT,NBNO,XEL,BPSS,NOQUAL)
C=======================================================================
C
C    -TEST DE VOISINAGE DES NOEUDS D'UN ELEMENT JOI4
C    -TEST DE PLANEITE DES FACES DE L'ELEMENT
C    -CALCUL DE LA MATRICE DE PASSAGE BPSS
C    -CALCUL DES COORDONNEES LOCALES XEL
C          ROUTINE FORTRAN PUR
C          CODE S. FELIX MAI 92
C=======================================================================
C  INPUT
C     XE     = COORDONNEES  DE L ELEMENT
C     SHPTOT = FONCTIONS DE FORME
C            = SHPTOT(1,...) = FONCTIONS DE FORME
C            = SHPTOT(2,...) = DERIVEE PAR RAPPORT A QSI
C            = SHPTOT(3,...) = DERIVEE PAR RAPPORT A ETA
C     NBNO   = NOMBRE DE NOEUDS DE L'ELEMENT
C     IVRF   = DEMANDE DE VERIFICATION DE L 'ELEMENT
C  OUTPUT
C     XEL    = COORDONNEES LOCALES
C     BPSS   = MATRICE DE PASAGE REPERE GLOBAL/REPERE LOCAL
C     NOQUAL = INDICE DE QUALITE
C            = 0 SI OK
C            = 1 SI NOEUD TROP VOISINS
C            = 2 SI NOEUD NON COPLANAIRES
C
C  REMARQUE : ATTENTION : DANS LES CAS CONTRAINTES PLANES, DEFO. PLANES
C             AXISYMETRIQUE, LA MATRICE TETA SERA UNE MATRICE DE
C             DIMENSION (2X2), ET SERA CONSTITUEE PAR LES VECTEURS
C             S1 ET SN. LES CAS CONT.PLANES,DEF.PLANES ET AXISYMETRIQUE
C             SERONT DONC SIMILAIRES AU CAS D'UN JOINT BIDIMENSIONNEL
C
C=======================================================================
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
      INTEGER IND4(0:5)
      DIMENSION XE(3,8),XEL(3,8),BPSS(3,3),SHPTOT(6,NBNO,*)
      DIMENSION U1(3),V1(3),XD(3,8),V2(3)
      DIMENSION S1(3),S2(3),SN(3)
      DATA IND4/4,1,2,3,4,1/
C
      NOQUAL = 0
      C1 = 0.0D0
      C2 = 0.0D0
      C3 = 0.0D0
C
C----------     VERIFICATION DE LA DISTANCE MINIMALE ENTRE LES POINTS
C----------     PAR COMPARAISON AVEC LE POURTOUR DU JOINT
C
      PP = 0.0D0
      DO 1 I=1,4
       II = I+1
       IF (II.EQ.5) II=1
       C1 = ABS(XE(1,I)-XE(1,II))
       C2 = ABS(XE(2,I)-XE(2,II))
       C3 = ABS(XE(3,I)-XE(3,II))
       C1 = C1*C1 + C2*C2 + C3*C3
       PP = PP + SQRT(C1)
    1 CONTINUE
C
      DMIN=PP/50.0D0
      DO 2 I=1,3
       I1 = I+1
       DO 2 N=I1,4
        IF ( ABS(XE(1,I)-XE(1,N)).LE.DMIN.AND.
     &       ABS(XE(2,I)-XE(2,N)).LE.DMIN.AND.
     &       ABS(XE(3,I)-XE(3,N)).LE.DMIN ) THEN
C                          NOEUDS TROP VOISINS
            NOQUAL=1
        ENDIF
    2 CONTINUE
C
C----------     CALCUL DE LA MATRICE DE  PASSAGE
C
      DO 6 I=1,3
         S1(I)=0.0D0
         S2(I)=0.0D0
         SN(I)=0.0D0
         V2(I)=0.0D0
    6 CONTINUE
C
      DO 7 I=1,NBNO
C
C-------------------TANGENTE AU POINT DE GAUSS 1 SELON QSI
C
         S1(1) = S1(1) + ( SHPTOT(2,I,1)*XE(1,I) )
         S1(2) = S1(2) + ( SHPTOT(2,I,1)*XE(2,I) )
         S1(3) = S1(3) + ( SHPTOT(2,I,1)*XE(3,I) )
C
C-------------------TANGENTE AU POINT DE GAUSS 1 SELON ETA
C
         V2(1) = V2(1) + ( SHPTOT(3,I,1)*XE(1,I) )
         V2(2) = V2(2) + ( SHPTOT(3,I,1)*XE(2,I) )
         V2(3) = V2(3) + ( SHPTOT(3,I,1)*XE(3,I) )
    7  CONTINUE
C
       CALL NORMER(S1)
       CALL NORMER(V2)
C-------------------NORMALE AU PLAN DU JOINT
C
       SN(1) = (S1(2)*V2(3)) - (S1(3)*V2(2))
       SN(2) = (S1(3)*V2(1)) - (S1(1)*V2(3))
       SN(3) = (S1(1)*V2(2)) - (S1(2)*V2(1))
       CALL NORMER(SN)
C
C-------------------ORTHOGONALISATION DE S2
C
       S2(1) = (SN(2)*S1(3)) - (SN(3)*S1(2))
       S2(2) = (SN(3)*S1(1)) - (SN(1)*S1(3))
       S2(3) = (SN(1)*S1(2)) - (SN(2)*S1(1))
       CALL NORMER(S2)
C
C-------------------STOCKAGE DE LA MATRICE DE PASSAGE
C
       DO 10 I=1,3
          BPSS(1,I) = S1(I)
          BPSS(2,I) = S2(I)
          BPSS(3,I) = SN(I)
   10  CONTINUE
C
C----------     CALCUL DES COORDONNEES LOCALES DE L'ELEMENT
C
C
C-------------------CHANGEMENT D'ORIGINE ( ORIGINE AU NOEUD 1)
C
*      DO 8 J=1,NBNO
*         DO 8 I=1,3
*            XD(I,J) = XE(I,J) - XE(I,1)
*   8  CONTINUE
C
C-------------------PROJECTION SUR LE PLAN DU JOINT
C
*      DO 9 J=1,NBNO
*       XEL(1,J)=(XD(1,J)*S1(1))+(XD(2,J)*S1(2))+(XD(3,J)*S1(3))
*       XEL(2,J)=(XD(1,J)*S2(1))+(XD(2,J)*S2(2))+(XD(3,J)*S2(3))
*       XEL(3,J)=0.0D0
*   9  CONTINUE
C+PPj
C
C----------    CALCUL DES COORDONNEES GLOBALES DU PLAN MOYEN DU JOINT
C              QUE L'ON STOCKE DANS LA FIN DE XEL
C
       NBNOS2=NBNO/2
       DO J=1,NBNOS2
         DO I=1,3
           XEL(I,J+NBNOS2) = (XE(I,J) + XE(I,NBNOS2+J))/2
         ENDDO
       ENDDO
C
C-----------   CHANGEMENT D'ORIGINE DU PLAN MOYEN DU JOINT
C              (ORIGINE AU NOEUD 1)
C
       DO J=1,NBNOS2
         DO I=1,3
           XD(I,J) = XEL(I,J+NBNOS2) - XEL(I,1+NBNOS2)
         ENDDO
       ENDDO
C
C-----------   PROJECTION SUR LE PLAN DU JOINT ET STOCKAGE DANS LE
C              DEBUT DE XEL
C
       DO J=1,NBNOS2
         XEL(1,J)=(XD(1,J)*S1(1))+(XD(2,J)*S1(2))+(XD(3,J)*S1(3))
         XEL(2,J)=(XD(1,J)*S2(1))+(XD(2,J)*S2(2))+(XD(3,J)*S2(3))
         XEL(3,J)=0.0D0
       ENDDO
C
C+PPj
C
C----------     TEST DE PLANEITE
C
C     CALCUL DES 4 NORMALES LOCALES
C
      DO 3 K=1,4
       KP1 = IND4(K+1)
       KM1 = IND4(K-1)
       DO 4  J=1,3
        U1(J) = XE(J,KP1) - XE(J,K)
        V1(J) = XE(J,KM1) - XE(J,K)
    4  CONTINUE
C
       XD(1,K) = U1(2)*V1(3) - U1(3)*V1(2)
       XD(2,K) = U1(3)*V1(1) - U1(1)*V1(3)
       XD(3,K) = U1(1)*V1(2) - U1(2)*V1(1)
C
       XXD = (XD(1,K)**2) + (XD(2,K)**2) + (XD(3,K)**2)
       XXD = SQRT(XXD)
C
      IF (XXD.GT.0.0D0) THEN
        DO 5 J=1,3
         XD(J,K) = XD(J,K)/XXD
    5   CONTINUE
      ENDIF
    3 CONTINUE
C
C     CALCUL DE LA NON PLANEITE
C
      COS13 = XD(1,3)*XD(1,1) + XD(2,3)*XD(2,1) + XD(3,3)*XD(3,1)
      COS24 = XD(1,4)*XD(1,2) + XD(2,4)*XD(2,2) + XD(3,4)*XD(3,2)
      IF (MIN(COS13,COS24).LT.0.99999) THEN
C            NON PLANEITE DE 0.25 DEGRE OU PLUS
C            0.99999, QUI EQUIVAUT A 1 DEGRE, EST INSUFFISANT
       NOQUAL = 2
      ENDIF
C
      RETURN
      END