vc01ft.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE vc01ft
!                   *****************
!
     &( xmul,sf,f,x,y,z,
     &  ikle1,ikle2,ikle3,nbor,nelem,nelmax,w1,w2,w3)
!
!***********************************************************************
! BIEF   V6P1                                   21/08/2010
!***********************************************************************
!
!brief    COMPUTES THE FOLLOWING VECTOR IN FINITE ELEMENTS:
!code
!+                    /
!+    VEC(I) = XMUL  /    PSI(I) * F  D(OMEGA)
!+                  /OMEGA
!+
!+    PSI(I) IS A BASE OF TYPE P1 SEGMENT BUT ON A MESH WITH
!+    VERTICAL TRIANGLES IN THE X,Y,Z SPACE
!+
!+    F IS A VECTOR OF TYPE IELMF
!
!warning  THE JACOBIAN MUST BE POSITIVE
!warning  THE RESULT IS IN W IN NOT ASSEMBLED FORM
!
!history  J-M HERVOUET (LNH)    ; F LEPEINTRE (LNH)
!+        26/04/04
!+        V5P5
!+
!
!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
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| F              |-->| FUNCTION USED IN THE VECTOR FORMULA
!| IKLE1          |-->| FIRST POINT OF TRIANGLES
!| IKLE2          |-->| SECOND POINT OF TRIANGLES
!| IKLE3          |-->| THIRD POINT OF TRIANGLES
!| NBOR           |-->| GLOBAL NUMBER OF BOUNDARY POINTS
!| NELEM          |-->| NUMBER OF ELEMENTS
!| NELMAX         |-->| MAXIMUM NUMBER OF ELEMENTS
!| SF             |-->| BIEF_OBJ STRUCTURE OF F
!| W1             |<--| RESULT IN NON ASSEMBLED FORM
!| W2             |<--| RESULT IN NON ASSEMBLED FORM
!| W3             |<--| RESULT IN NON ASSEMBLED FORM
!| X              |-->| ABSCISSAE OF POINTS IN THE MESH
!| XMUL           |-->| MULTIPLICATION COEFFICIENT
!| Y              |-->| ORDINATES OF POINTS IN THE MESH
!| Z              |-->| ELEVATIONS OF POINTS IN THE MESH
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief, ex_vc01ft => vc01ft
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN) :: NELEM,NELMAX
      INTEGER, INTENT(IN) :: NBOR(*)
      INTEGER, INTENT(IN) :: IKLE1(nelmax),IKLE2(nelmax),IKLE3(nelmax)
!
      DOUBLE PRECISION, INTENT(IN) :: X(*),Y(*),Z(*)
      DOUBLE PRECISION, INTENT(INOUT) ::W1(nelmax),W2(nelmax),W3(nelmax)
      DOUBLE PRECISION, INTENT(IN) :: XMUL
!
!     STRUCTURE OF F AND REAL DATA
!
      TYPE(bief_obj), INTENT(IN) :: SF
      DOUBLE PRECISION, INTENT(IN) :: F(*)
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER IELEM,IELMF,I1,I2,I3
      DOUBLE PRECISION XSUR12,COEF,F123,F1,F2,F3,X1,X2,X3,Y1,Y2
      DOUBLE PRECISION Y3,Z1,Z2,Z3,S
!
      INTRINSIC sqrt
!
!***********************************************************************
!
      ielmf=sf%ELM
!
!-----------------------------------------------------------------------
!
!     F IS LINEAR BY BOUNDARY SIDE
!
      IF(ielmf.EQ.61.OR.ielmf.EQ.81) THEN
!
        xsur12 = xmul/12.d0
!
!   LOOP ON THE BOUNDARY SIDES
!
        DO ielem = 1,nelem
!
!         GLOBAL NUMBERING OF THE SIDE NODES
!
          i1 = nbor(ikle1(ielem))
          i2 = nbor(ikle2(ielem))
          i3 = nbor(ikle3(ielem))
!
          x1 = x(i1)
          y1 = y(i1)
          z1 = z(i1)
!
          x2 = x(i2)-x1
          x3 = x(i3)-x1
          y2 = y(i2)-y1
          y3 = y(i3)-y1
          z2 = z(i2)-z1
          z3 = z(i3)-z1
!
          f1 = f(ikle1(ielem))
          f2 = f(ikle2(ielem))
          f3 = f(ikle3(ielem))
          f123  = f1 + f2 + f3
!
!         COMPUTES THE AREA OF THE TRIANGLE (BY VECTOR PRODUCT)
!
          s=0.5d0*sqrt(  (y2*z3-y3*z2)**2
     &                  +(x3*z2-x2*z3)**2
     &                  +(x2*y3-x3*y2)**2  )
!
          coef=xsur12*s
!
          w1(ielem) = coef * ( f123 + f1 )
          w2(ielem) = coef * ( f123 + f2 )
          w3(ielem) = coef * ( f123 + f3 )
!
        ENDDO ! IELEM
!
!-----------------------------------------------------------------------
!
      ELSE
!
!-----------------------------------------------------------------------
!
        WRITE(lu,101) ielmf,sf%NAME
101     FORMAT(1x,'VC01FT (BIEF) :',/,
     &         1x,'DISCRETIZATION OF F NOT AVAILABLE:',1i6,
     &         1x,'REAL NAME: ',a6)
        CALL plante(1)
        stop
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!     NOTE: FOR A PLANE TRIANGLE (VC01AA):
!
!     XSUR12 = XMUL / 12.D0
!
!     DO IELEM = 1 , NELEM
!
!       F1  = F(IKLE1(IELEM))
!       F2  = F(IKLE2(IELEM))
!       F3  = F(IKLE3(IELEM))
!       F123  = F1 + F2 + F3
!
!       COEF = XSUR12 * SURFAC(IELEM)
!
!       W1(IELEM) = COEF * ( F123 + F1 )
!       W2(IELEM) = COEF * ( F123 + F2 )
!       W3(IELEM) = COEF * ( F123 + F3 )
!
!     ENDDO ! IELEM
!
!-----------------------------------------------------------------------
!
      RETURN
      END