Sources Esope | Cast3M

angqua
C ANGQUA    SOURCE    PV        07/11/23    21:15:23     5978
      SUBROUTINE ANGQUA (NOD1,NOD2,NOD3,NOD4,  ANG1,ANG2,ANG3,ANG4)
      implicit real*8(a-h,o-z)
************************************************************************
*
*                             A N G Q U A
*                             -----------
*
* FONCTION:
* ---------
*
*     CALCUL DES ANGLES AUX COINS D'UN QUADRANGLE.
*
* MODULES UTILISES:
* -----------------
*
      IMPLICIT INTEGER(I-N)
-INC CCREEL
 
-INC PPARAM
-INC CCOPTIO
-INC SMCOORD
*
* PARAMETRES:   (E)=ENTREE   (S)=SORTIE   (+ = CONTENU DANS UN COMMUN)
* -----------
*
*    +XCOOR   (E)  VOIR LE COMMUN "COPTIO" ET LE SEGMENT "MCOORD".
*                  (SEGMENT SUPPOSE ACTIF)
*    +IDIM    (E)  VOIR LE COMMUN "COPTIO".
*     NOD?    (E)  NUMERO DU NOEUD "?", "?"=1,4.
*     ANG?    (S)  ANGLE AU NOEUD "?" ENTRE 0 ET 2*PI.
*
      REAL *8 ANG1,ANG2,ANG3,ANG4
*
* VARIABLES:
* ----------
*
      INTEGER IND4(0:5)
      REAL*8 XE(3,4),V1(3),V2(3),ANG(4)
*
* AUTEUR, DATE DE CREATION:
* -------------------------
*
*     PASCAL MANIGOT     5 OCTOBRE 1990
*
************************************************************************
*
*     INDICAGE CIRCULAIRE:
      DATA IND4/4,1,2,3,4,1/
*
      IDIMP1 = IDIM + 1
      XE(1,1) = XCOOR((NOD1-1)*IDIMP1+1)
      XE(1,2) = XCOOR((NOD2-1)*IDIMP1+1)
      XE(1,3) = XCOOR((NOD3-1)*IDIMP1+1)
      XE(1,4) = XCOOR((NOD4-1)*IDIMP1+1)
      XE(2,1) = XCOOR((NOD1-1)*IDIMP1+2)
      XE(2,2) = XCOOR((NOD2-1)*IDIMP1+2)
      XE(2,3) = XCOOR((NOD3-1)*IDIMP1+2)
      XE(2,4) = XCOOR((NOD4-1)*IDIMP1+2)
      IF (IDIM.EQ.3) THEN
         XE(3,1) = XCOOR((NOD1-1)*IDIMP1+3)
         XE(3,2) = XCOOR((NOD2-1)*IDIMP1+3)
         XE(3,3) = XCOOR((NOD3-1)*IDIMP1+3)
         XE(3,4) = XCOOR((NOD4-1)*IDIMP1+3)
      ELSE
         XE(3,1) = 0.
         XE(3,2) = 0.
         XE(3,3) = 0.
         XE(3,4) = 0.
      END IF
*
*     CALCUL DE LA NORMALE "EXTERIEURE" MOYENNE:
*
      XN = (XE(2,3)-XE(2,1))*(XE(3,4)-XE(3,2))
     &     - (XE(3,3)-XE(3,1))*(XE(2,4)-XE(2,2))
      YN = (XE(3,3)-XE(3,1))*(XE(1,4)-XE(1,2))
     &     - (XE(1,3)-XE(1,1))*(XE(3,4)-XE(3,2))
      ZN = (XE(1,3)-XE(1,1))*(XE(2,4)-XE(2,2))
     &     - (XE(2,3)-XE(2,1))*(XE(1,4)-XE(1,2))
      XNORM1 = SQRT(XN**2 + YN**2 + ZN**2)
      XN = XN / XNORM1
      YN = YN / XNORM1
      ZN = ZN / XNORM1
*
*     CALCUL DES ANGLES:
*
      DO 100 K=1,4
         KP1 = IND4(K+1)
         KM1 = IND4(K-1)
         DO 110 J=1,3
            V1(J) = XE(J,KP1) - XE(J,K)
            V2(J) = XE(J,KM1) - XE(J,K)
 110     CONTINUE
*        END DO
         XNORM1 = SQRT(V1(1)**2 + V1(2)**2 + V1(3)**2)
     &        * SQRT(V2(1)**2 + V2(2)**2 + V2(3)**2)
         COSINU = (V1(1)*V2(1) + V1(2)*V2(2) + V1(3)*V2(3)) / XNORM1
         ANG(K) = ACOS(COSINU)
         XPVECT = V1(2)*V2(3) - V1(3)*V2(2)
         YPVECT = V1(3)*V2(1) - V1(1)*V2(3)
         ZPVECT = V1(1)*V2(2) - V1(2)*V2(1)
         S = XPVECT*XN + YPVECT*YN + ZPVECT*ZN
         IF (S .LT. 0.) THEN
            ANG(K) = (2. * REAL(XPI)) - ANG(K)
         END IF
 100  CONTINUE
*     END DO
*
      ANG1 = ANG(1)
      ANG2 = ANG(2)
      ANG3 = ANG(3)
      ANG4 = ANG(4)
*
      END