vc11aa2.f

Go to the documentation of this file.

!                   ******************
                    SUBROUTINE vc11aa2
!                   ******************
!
     &(xmul,sf,sg,sh,f,g,h,xel,yel,
     & ikle1,ikle2,ikle3,nelem,nelmax,w1,w2,w3,icoord)
!
!***********************************************************************
! BIEF   V7P0                                   21/08/2010
!***********************************************************************
!
!brief    COMPUTES THE FOLLOWING TERMS:
!code
!+   (EXAMPLE OF THE X COMPONENT, WHICH CORRESPONDS TO ICOORD=1)
!+
!+                       /                  DF
!+    VEC(I)  =  XMUL   /  ( G . H .  P  *( --  )) D(OMEGA)
!+                     /OMEGA          I    DX
!+
!+    P   IS A LINEAR BASE
!+     I
!+
!+    F IS A VECTOR OF TYPE P1
!+    G IS A VECTOR OF TYPE DISCONTINUOUS P1
!+    H IS A VECTOR OF TYPE P0
!
!warning  THE JACOBIAN MUST BE POSITIVE
!warning  THE RESULT IS IN W IN NOT ASSEMBLED FORM
!
!history  J-M HERVOUET (LNH)    ; F LEPEINTRE (LNH)
!+        09/12/94
!+        V5P1
!+
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        13/07/2010
!+        V6P0
!+   Translation of French comments within the FORTRAN sources into
!+   English comments
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        21/08/2010
!+        V6P0
!+   Creation of DOXYGEN tags for automated documentation and
!+   cross-referencing of the FORTRAN sources
!
!history  J-M HERVOUET (EDF LAB, LNHE)
!+        12/05/2014
!+        V7P0
!+   Discontinuous elements better treated: new types 15, 16 and 17 for
!+   discontinuous linear, quasi-bubble, and quadratic, rather than
!+   using component DIMDISC=11, 12 or 13.
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| F              |-->| FUNCTION USED IN THE VECTOR FORMULA
!| G              |-->| FUNCTION USED IN THE VECTOR FORMULA
!| H              |-->| FUNCTION USED IN THE VECTOR FORMULA
!| ICOORD         |-->| 1: DERIVATIVE ALONG X, 2: ALONG Y
!| IKLE1          |-->| FIRST POINT OF TRIANGLES
!| IKLE2          |-->| SECOND POINT OF TRIANGLES
!| IKLE3          |-->| THIRD POINT OF TRIANGLES
!| NELEM          |-->| NUMBER OF ELEMENTS
!| NELMAX         |-->| MAXIMUM NUMBER OF ELEMENTS
!| SF             |-->| BIEF_OBJ STRUCTURE OF F
!| SG             |-->| BIEF_OBJ STRUCTURE OF G
!| SH             |-->| BIEF_OBJ STRUCTURE OF H
!| SURFAC         |-->| AREA OF TRIANGLES
!| W1             |<--| RESULT IN NON ASSEMBLED FORM
!| W2             |<--| RESULT IN NON ASSEMBLED FORM
!| W3             |<--| RESULT IN NON ASSEMBLED FORM
!| XEL            |-->| ABSCISSAE OF POINTS IN THE MESH, PER ELEMENT
!| XMUL           |-->| MULTIPLICATION COEFFICIENT
!| YEL            |-->| ORDINATES OF POINTS IN THE MESH, PER ELEMENT
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief, ex_vc11aa2 => vc11aa2
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN) :: NELEM,NELMAX,ICOORD
      INTEGER, INTENT(IN) :: IKLE1(nelmax),IKLE2(nelmax),IKLE3(nelmax)
!
      DOUBLE PRECISION, INTENT(IN) :: XEL(nelmax,*),YEL(nelmax,*)
      DOUBLE PRECISION,INTENT(INOUT)::W1(nelmax),W2(nelmax),W3(nelmax)
      DOUBLE PRECISION, INTENT(IN) :: XMUL
!
!     STRUCTURES OF F, G AND REAL DATA
!
      TYPE(bief_obj), INTENT(IN) :: SF,SG,SH
      DOUBLE PRECISION, INTENT(IN) :: F(*),G(*),H(*)
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER IELEM,IELMF,IELMG,IELMH
      DOUBLE PRECISION XSUR24 ,F1,F2,F3,G1,G2,G3,X2,X3,Y2,Y3,KSAT
!
!-----------------------------------------------------------------------
!
      xsur24= xmul / 24.d0
!
!-----------------------------------------------------------------------
!
      ielmf=sf%ELM
      ielmg=sg%ELM
      ielmh=sh%ELM
!
!-----------------------------------------------------------------------
!
!     F IS P1, G IS DISCONTINUOUS P1, H IS P0
!
      IF(ielmg.EQ.15.AND.ielmf.EQ.11.AND.ielmh.EQ.10) THEN
!
!     X COORDINATE
!
      IF(icoord.EQ.1) THEN
!
      DO ielem = 1 , nelem
!
        ksat=h(ielem)
        g1 = g(ielem)
        g2 = g(ielem+nelem)
        g3 = g(ielem+2*nelem)
        f1 = f(ikle1(ielem))
        f2 = f(ikle2(ielem))
        f3 = f(ikle3(ielem))
        y2 = yel(ielem,2)
        y3 = yel(ielem,3)
!
        w1(ielem)=(y2*(-g3*f3+g3*f1-g2*f3+g2*f1-2*g1*f3+2*g1*f1)+y3*(
     &             g3*f2-g3*f1+g2*f2-g2*f1+2*g1*f2-2*g1*f1))* xsur24
        w2(ielem)=(y2*(-g3*f3+g3*f1-2*g2*f3+2*g2*f1-g1*f3+g1*f1)+y3*(
     &             g3*f2-g3*f1+2*g2*f2-2*g2*f1+g1*f2-g1*f1))* xsur24
        w3(ielem)=(y2*(-2*g3*f3+2*g3*f1-g2*f3+g2*f1-g1*f3+g1*f1)+y3*(
     &             2*g3*f2-2*g3*f1+g2*f2-g2*f1+g1*f2-g1*f1))* xsur24
!
        w1(ielem)=ksat*w1(ielem)
        w2(ielem)=ksat*w2(ielem)
        w3(ielem)=ksat*w3(ielem)
!
      ENDDO
!
      ELSEIF(icoord.EQ.2) THEN
!
!     Y COORDINATE
!
      DO ielem = 1 , nelem
!
        ksat=h(ielem)
        g1 = g(ielem)
        g2 = g(ielem+nelem)
        g3 = g(ielem+2*nelem)
        f1 = f(ikle1(ielem))
        f2 = f(ikle2(ielem))
        f3 = f(ikle3(ielem))
        x2 = xel(ielem,2)
        x3 = xel(ielem,3)
!
        w1(ielem)=(x2*(g3*f3-g3*f1+g2*f3-g2*f1+2*g1*f3-2*g1*f1)+x3*(-
     &             g3*f2+g3*f1-g2*f2+g2*f1-2*g1*f2+2*g1*f1)) * xsur24
        w2(ielem)=(x2*(g3*f3-g3*f1+2*g2*f3-2*g2*f1+g1*f3-g1*f1)+x3*(-
     &             g3*f2+g3*f1-2*g2*f2+2*g2*f1-g1*f2+g1*f1)) * xsur24
        w3(ielem)=(x2*(2*g3*f3-2*g3*f1+g2*f3-g2*f1+g1*f3-g1*f1)+x3*(-
     &             2*g3*f2+2*g3*f1-g2*f2+g2*f1-g1*f2+g1*f1)) * xsur24
!
        w1(ielem)=ksat*w1(ielem)
        w2(ielem)=ksat*w2(ielem)
        w3(ielem)=ksat*w3(ielem)
!
      ENDDO
!
      ELSE
!
          WRITE(lu,21) icoord
21        FORMAT(1x,'VC11AA2 (BIEF) : IMPOSSIBLE COMPONENT ',
     &              1i6,' CHECK ICOORD')
          CALL plante(0)
          stop
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!
!-----------------------------------------------------------------------
!
      ELSE
!
!-----------------------------------------------------------------------
!
        WRITE(lu,101) ielmf,sf%NAME
        WRITE(lu,201) ielmg,sg%NAME
        WRITE(lu,301)
101     FORMAT(1x,'VC11AA2 (BIEF) :',/,
     &         1x,'DISCRETIZATION OF F:',1i6,
     &         1x,'REAL NAME: ',a6)
201     FORMAT(1x,'DISCRETIZATION OF G:',1i6,
     &         1x,'REAL NAME: ',a6)
301     FORMAT(1x,'CASE NOT IMPLEMENTED')
        CALL plante(1)
        stop
!
      ENDIF
!
!-----------------------------------------------------------------------
!
      RETURN
      END