mt06oc.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE mt06oc
!                   *****************
!
     &(a11,a12,a13,a22,a23,a33,
     & xmul,sf,f,lgseg,ikle1,ikle2,ikle3,nbor,nelem,nelmax)
!
!***********************************************************************
! BIEF   V6P1                                   21/08/2010
!***********************************************************************
!
!brief    COMPUTES THE COEFFICIENTS OF THE FOLLOWING MATRIX:
!code
!+                              /
!+                    A    =   /  F (P *P )*J(X,Y) DX
!+                     I J    /L      I  J
!+
!+     BY ELEMENTARY CELL; THE ELEMENT IS THE P2 SEGMENT
!+
!+     J(X,Y): JACOBIAN OF THE ISOPARAMETRIC TRANSFORMATION
!
!warning  THE JACOBIAN MUST BE POSITIVE
!
!history  A FROEHLY (MATMECA)
!+        01/07/08
!+        V5P9
!+
!
!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
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| A11            |<--| ELEMENTS OF MATRIX
!| A12            |<--| ELEMENTS OF MATRIX
!| A13            |<--| ELEMENTS OF MATRIX
!| A22            |<--| ELEMENTS OF MATRIX
!| A23            |<--| ELEMENTS OF MATRIX
!| A33            |<--| ELEMENTS OF MATRIX
!| F              |-->| FUNCTION F USED IN THE FORMULA
!| IKLE1          |-->| FIRST POINTS OF SEGMENTS
!| IKLE2          |-->| SECOND POINTS OF SEGMENTS
!| IKLE3          |-->| THIRD POINTS OF SEGMENTS (QUADRATIC)
!| LGSEG          |-->| LENGTH OF SEGMENTS
!| NBOR           |-->| GLOBAL NUMBER OF BOUNDARY POINTS
!| NELEM          |-->| NUMBER OF ELEMENTS
!| NELMAX         |-->| MAXIMUM NUMBER OF ELEMENTS
!| SF             |-->| BIEF_OBJ STRUCTURE OF F
!| SURFAC         |-->| AREA OF TRIANGLES
!| XMUL           |-->| MULTIPLICATION FACTOR
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief!, EX_MT06OC => MT06OC
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN) :: NELEM,NELMAX,NBOR(nelmax,*)
      INTEGER, INTENT(IN) :: IKLE1(*),IKLE2(*),IKLE3(*)
!
      DOUBLE PRECISION, INTENT(IN) :: XMUL
!
      DOUBLE PRECISION, INTENT(IN) :: F(*)
!
!     STRUCTURE OF F
      TYPE(bief_obj), INTENT(IN) :: SF
!
      DOUBLE PRECISION, INTENT(IN)    :: LGSEG(nelmax)
      DOUBLE PRECISION, INTENT(INOUT) :: A11(nelmax),A12(nelmax)
      DOUBLE PRECISION, INTENT(INOUT) :: A13(nelmax),A22(nelmax)
      DOUBLE PRECISION, INTENT(INOUT) :: A23(nelmax),A33(nelmax)
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER IELEM,IELMF
      DOUBLE PRECISION SUR30,SUR60,SUR420,DET1,F1,F2,F3
!
!-----------------------------------------------------------------------
!
      sur30  = xmul/30.d0
      sur60  = xmul/60.d0
      sur420  = xmul/420.d0
!
!-----------------------------------------------------------------------
!
      ielmf = sf%ELM
!
!     F CONSTANT BY SEGMENT, IN A BOUNDARY ARRAY
!
      IF(ielmf.EQ.0) THEN
!
      DO ielem = 1 , nelem
      f1 = f(ielem)
      det1 = lgseg(ielem) * sur30
!
      a11(ielem) = det1 * (4.d0*f1)
      a12(ielem) = det1 * (-f1)
      a13(ielem) = det1 * (2.d0*f1)
      a22(ielem) = a11(ielem)
      a23(ielem) = a13(ielem)
      a33(ielem) = det1 * (16.d0*f1)
!
      ENDDO ! IELEM
!
!     F LINEAR BY SEGMENT, IN A BOUNDARY ARRAY
!     NOTE: IKLE IS HERE A BOUNDARY IKLE
!
      ELSEIF(ielmf.EQ.1) THEN
!
      DO ielem = 1 , nelem
!
      f1 = f(ikle1(ielem))
      f2 = f(ikle2(ielem))
!
      det1 = lgseg(ielem) * sur60
!
      a11(ielem) = det1 * (7.d0*f1+f2)
      a12(ielem) = det1 * (-f1-f2)
      a13(ielem) = det1 * (4.d0*f1)
      a22(ielem) = det1 * (f1+7.d0*f2)
      a23(ielem) = det1 * (4.d0*f2)
      a33(ielem) = det1 * 16.d0 * (f1+f2)
!
      ENDDO ! IELEM
!
!     F LINEAR, IN AN ARRAY DEFINED ON THE DOMAIN
!
      ELSEIF(ielmf.EQ.11.OR.ielmf.EQ.21) THEN
!
      DO ielem = 1 , nelem
!
      f1 = f(nbor(ielem,1))
      f2 = f(nbor(ielem,2))
!
      det1 = lgseg(ielem) * sur60
!
      a11(ielem) = det1 * (7.d0*f1+f2)
      a12(ielem) = det1 * (-f1-f2)
      a13(ielem) = det1 * (4.d0*f1)
      a22(ielem) = det1 * (f1+7.d0*f2)
      a23(ielem) = det1 * (4.d0*f2)
      a33(ielem) = det1 * 16.d0 * (f1+f2)
!
      ENDDO ! IELEM
!
!     F QUADRATIC BY SEGMENT, IN A BOUNDARY ARRAY
!     NOTE: IKLE IS HERE A BOUNDARY IKLE
!
      ELSEIF(ielmf.EQ.2) THEN
!
      DO ielem = 1 , nelem
!
      f1 = f(ikle1(ielem))
      f2 = f(ikle2(ielem))
      f3 = f(ikle3(ielem))
      det1 = lgseg(ielem) * sur420
!
      a11(ielem) = det1 * (39.d0*f1-3.d0*f2+20.d0*f3)
      a12(ielem) = det1 * (-3.d0*f1-3.d0*f2-8.d0*f3)
      a13(ielem) = det1 * (20.d0*f1-8.d0*f2+16.d0*f3)
      a22(ielem) = det1 * (-3.d0*f1+39.d0*f2+20.d0*f3)
      a23(ielem) = det1 * (-8.d0*f1+20.d0*f2+16.d0*f3)
      a33(ielem) = det1 * 16.d0 * (f1+f2+12.d0*f3)
!
      ENDDO ! IELEM
!
!     F QUADRATIC, IN AN ARRAY DEFINED ON THE DOMAIN
!
      ELSEIF(ielmf.EQ.13) THEN
!
      DO ielem = 1 , nelem
!
      f1 = f(nbor(ielem,1))
      f2 = f(nbor(ielem,2))
      f3 = f(nbor(ielem,3))
!
      det1 = lgseg(ielem) * sur420
!
      a11(ielem) = det1 * (39.d0*f1-3.d0*f2+20.d0*f3)
      a12(ielem) = det1 * (-3.d0*f1-3.d0*f2-8.d0*f3)
      a13(ielem) = det1 * (20.d0*f1-8.d0*f2+16.d0*f3)
      a22(ielem) = det1 * (-3.d0*f1+39.d0*f2+20.d0*f3)
      a23(ielem) = det1 * (-8.d0*f1+20.d0*f2+16.d0*f3)
      a33(ielem) = det1 * 16.d0 * (f1+f2+12.d0*f3)
!
      ENDDO ! IELEM
!
!     OTHER TYPES OF DISCRETISATION OF F
!
      ELSE
!
        WRITE(lu,101) ielmf,sf%NAME
101     FORMAT(1x,'MT06OC (BIEF) :',/,
     &         1x,'DISCRETIZATION OF F NOT AVAILABLE:',1i6,
     &         1x,'REAL NAME: ',a6)
        CALL plante(1)
        stop
!
      ENDIF
!
!-----------------------------------------------------------------------
!
      RETURN
      END