  MODULE msh
!=========================== DOCUMENTATION ============================= 
!DOCDEBm*    UTI/msh  *
!
!    NOM
! msh
!    DESCRIPTION
!*
! Ensemble de routines necessaires au maillage
! Set of routines required for meshing
!*
! 
!    SOURCE
! Ce logiciel est regi par la licence [CeCILL-B]
! This software is governed by the [CeCILL-B] license
!=========================== DEBUT DES DECLARATIONS ==================== 
!.1-----  Implicit, Use
  use geo
  IMPLICIT NONE
!=========================== DEBUT DU CODE EXECUTABLE ==================

!----- Debut type ------------------------------------------------------
 TYPE vertex
    real(kind=kind(0.d0))      :: x, y, z
    integer(kind=16) :: morton ! code morton du vertex
    integer  :: num    ! numero du vertex
 END TYPE vertex
!----- Fin type --------------------------------------------------------
!----- Debut type ------------------------------------------------------
 TYPE mface
    integer(kind=16), dimension(4)  :: mv       ! codes morton des 2/3/4 vertex composant la face
    integer(kind=16)                :: morton   ! code morton de la face
    integer                         :: num      ! numero de la face
    integer                         :: qe1,qe2  ! codes quadtree/Octree des deux elements voisins si qe2 <= 0 ---> condition aux limites
    integer                         :: be1,be2  ! numero de bloc des deux elements voisins si be2 <= 0 ---> condition aux limites
    integer                         :: nr1,nr2  ! niveau de raffinement des deux elements voisins 
    integer         , dimension(2)  :: fac_loc  ! numeros locaux des deux elements connectes
    real(kind=kind(0.d0)), dimension(3)      :: vn       ! normale 
 END TYPE mface
!----- Fin type --------------------------------------------------------
!----- Debut type ------------------------------------------------------
 TYPE chaine_MF
    type(mface)              :: f
    type(chaine_MF), pointer :: prec,suiv
 END TYPE chaine_MF
!----- Fin type --------------------------------------------------------

!----- Debut type ------------------------------------------------------
 TYPE noe
    type(point)       :: coord  ! Barycentre
    integer(kind=16) :: morton ! Code morton du noeud
    integer(kind=16)  :: num    ! Numero du noeud
 END TYPE noe
!----- Fin type --------------------------------------------------------



 CONTAINS

! Routines necessaires a read_gmsh
 subroutine ini_mface(f)
     type(mface)   :: f
     f%mv=0       
     f%morton=0   
     f%num   =0   
     f%qe1=0 ; f%qe2=0  
     f%be1=0 ; f%be2=0  
     f%nr1=0 ; f%nr2=0  
     f%fac_loc=0  
 end subroutine ini_mface
!-------------------------------------------------------
   subroutine range_mface(f,lcf,ier)
    implicit none
    Type(chaine_MF), Pointer        :: lcf
    Type(mface)                     :: f
    integer                         :: ier
    logical                         :: droite
    Type(chaine_MF), Pointer        :: chloc,temp

    ! interdit de rentrer dans cette routine si la chaine pointe sur null
    droite=.false.
    if (f%morton>lcf%f%morton) droite=.true.
    do while (associated(lcf))
          temp=>lcf
          if (f%morton>lcf%f%morton) then
             if (droite) then
                lcf=>lcf%suiv
             else
                call insertdroite_MF(f,lcf,ier)
                return
             end if
          else
             if (f%morton==lcf%f%morton) then
                ier=1
!---- on actualise la face qui existe deja
                if(lcf%f%be1.le.0)then
                  lcf%f%qe1=f%qe1
                  lcf%f%be1=f%be1
                  lcf%f%nr1=f%nr1
                  lcf%f%mv =f%mv
                  lcf%f%vn =f%vn
                  lcf%f%fac_loc(1) =f%fac_loc(1)
                else
                  lcf%f%qe2=lcf%f%qe1
                  lcf%f%be2=lcf%f%be1
                  lcf%f%nr2=lcf%f%nr1
                  lcf%f%fac_loc(2) =lcf%f%fac_loc(1)
                  lcf%f%qe1=f%qe1
                  lcf%f%be1=f%be1
                  lcf%f%nr1=f%nr1
                  lcf%f%mv =f%mv
                  lcf%f%vn =f%vn
                  lcf%f%fac_loc(1) =f%fac_loc(1)
                endif
!----
                return
             else
                if (droite) then
                   call insertgauche_MF(f,lcf,ier)
                   return
                else
                   lcf=>lcf%prec
                end if
             end if
          end if
       end do
!cas du bout de liste droite ou gauche
       allocate(chloc)
       chloc%f=f
       if (droite) then
!ajout a droite de temp
          chloc%prec=>temp
          nullify(chloc%suiv)
          temp%suiv=>chloc
       else
!ajout a gauche  de temp
          chloc%suiv=>temp
          nullify(chloc%prec)
          temp%prec=>chloc
       end if
       lcf=>chloc!temp
       ier=0
       return
  end subroutine range_mface
!-------------------------------------------------------

 subroutine mortonface(f,lnf,nn,ndim)
  use geo_typ
  use zorder
  implicit none
  integer        :: nn,ndim,i,mortonx,mortony,mortonz
  type(mface)             :: f
  type(noe), dimension(8) :: lnf
  real(kind=kind(0.d0))   :: x,y,z
  integer(kind=16)       :: nmorton

  x= 0.d0 ; y= 0.d0 ; z= 0.d0
  do i=1,nn
    f%mv(i)=lnf(i)%morton
    x=x+lnf(i)%coord%x
    y=y+lnf(i)%coord%y
    z=z+lnf(i)%coord%z
  enddo
  x=x/nn ; y=y/nn ;z=z/nn
  mortonx=coord_real2int(x,x_min_morton,petit_morton)
  mortony=coord_real2int(y,y_min_morton,petit_morton)
  if (ndim == 2) then
    call zordercurve2d(mortonx,mortony,nmorton)
  else
    mortonz=coord_real2int(z,z_min_morton,petit_morton)
    call zordercurve3d(mortonx,mortony,mortonz,nmorton)
  endif
  f%morton=nmorton
  select case (nn)
      case(2) ! L2
        f%vn(1) = -(lnf(2)%coord%y-lnf(1)%coord%y)
        f%vn(2) =  (lnf(2)%coord%x-lnf(1)%coord%x)
        f%vn(3) = 0.d0
      case(3) ! T3
        f%vn(1) = ((lnf(2)%coord%y-lnf(1)%coord%y)*(lnf(3)%coord%z-lnf(1)%coord%z)-&
        (lnf(2)%coord%z-lnf(1)%coord%z)*(lnf(3)%coord%y-lnf(1)%coord%y))*0.5d0
        f%vn(2) = ((lnf(2)%coord%z-lnf(1)%coord%z)*(lnf(3)%coord%x-lnf(1)%coord%x)-&
        (lnf(2)%coord%x-lnf(1)%coord%x)*(lnf(3)%coord%z-lnf(1)%coord%z))*0.5d0
        f%vn(3) = ((lnf(2)%coord%x-lnf(1)%coord%x)*(lnf(3)%coord%y-lnf(1)%coord%y)-&
        (lnf(2)%coord%y-lnf(1)%coord%y)*(lnf(3)%coord%x-lnf(1)%coord%x))*0.5d0
      case(4) ! Q4
        f%vn(1) = ((lnf(2)%coord%y-lnf(1)%coord%y)*(lnf(4)%coord%z-lnf(1)%coord%z)-&
        (lnf(2)%coord%z-lnf(1)%coord%z)*(lnf(4)%coord%y-lnf(1)%coord%y))
        f%vn(2) = ((lnf(2)%coord%z-lnf(1)%coord%z)*(lnf(4)%coord%x-lnf(1)%coord%x)-&
        (lnf(2)%coord%x-lnf(1)%coord%x)*(lnf(4)%coord%z-lnf(1)%coord%z))
        f%vn(3) = ((lnf(2)%coord%x-lnf(1)%coord%x)*(lnf(4)%coord%y-lnf(1)%coord%y)-&
        (lnf(2)%coord%y-lnf(1)%coord%y)*(lnf(4)%coord%x-lnf(1)%coord%x))
 end select

 end subroutine mortonface
!-------------------------------------------------------
 subroutine permnoe(ln,n)
  use geo_typ
  type(noe), dimension(8)  :: ln,temp
  integer(kind=16)        :: imin
  integer                  :: i,j,n
  imin=huge(imin)
  temp=ln
  j = 1
  do i=1,n
      if (ln(i)%morton <= imin) then
          imin=ln(i)%morton
          j=i
      endif
  enddo  
  ln(:)%coord%x= 0.d0 ;   ln(:)%coord%y= 0.d0 ;   ln(:)%coord%z= 0.d0
  ln(:)%morton= 0
  ln(:)%num   = 0
  
  do i=1,n
    if (j+i-1 <= n) then
      ln(i)=temp(j+i-1)
    else
      ln(i)=temp(j+i-n-1)
    endif
  enddo
 end subroutine permnoe

!-------------------------------------------------------
 subroutine perm4(array,n)
  integer(kind=16), dimension(4) :: array,temp
  integer(kind=16)               :: imin
  integer                         :: i,j,n
  imin=huge(imin)
  temp=array
  if (n==2)return
  do i=1,n
      if (array(i) <= imin) then
          imin=array(i)
          j=i
      endif
  enddo  
  array= 0
  do i=1,n
    if (j+i-1 <= n) then
      array(i)=temp(j+i-1)
    else
      array(i)=temp(j+i-n-1)
    endif
  enddo
 end subroutine perm4
!-------------------------------------------------------

 subroutine noe2bloc(ln,b,ns,ndim,nmorton)
  use geo_typ
  use zorder
  implicit none
  integer        :: ns,ndim,i,mortonx,mortony,mortonz
  type(bloc)              :: b
  type(noe), dimension(8) :: ln
  real(kind=kind(0.d0))   :: x,y,z
  integer(kind=16)       :: nmorton
      x=0.d0 ; y=0.d0 ;z=0.d0

    do i=1,ns
        b%s(i)%x=ln(i)%coord%x
        b%s(i)%y=ln(i)%coord%y
        b%s(i)%z=ln(i)%coord%z
        x=x+ln(i)%coord%x
        y=y+ln(i)%coord%y
        z=y+ln(i)%coord%z
    enddo
    x=x/ns ; y=y/ns ;z=z/ns
    mortonx=coord_real2int(x,x_min_morton,petit_morton)
    mortony=coord_real2int(y,y_min_morton,petit_morton)
    if (ndim == 2) then
      call zordercurve2d(mortonx,mortony,nmorton)
    else
      mortonz=coord_real2int(z,z_min_morton,petit_morton)
      call zordercurve3d(mortonx,mortony,mortonz,nmorton)
    endif
    !b%morton  = nmorton
  end subroutine noe2bloc
!-------------------------------------------------------
 subroutine insertdroite_MF(f,lcf,ier)
    implicit none
    integer                 :: ier
    Type(chaine_MF), Pointer     :: lcf
    Type(mface)                  :: f
    Type(chaine_MF), Pointer     :: temp, locch

!cas du bout de liste droite ou gauche
      temp=>lcf%suiv
       allocate (locch)
       locch%f=f
       locch%suiv=>temp
       temp%prec=>locch
       locch%prec=>lcf
       lcf%suiv=>locch
!modif
       lcf=>locch
       ier=0
       return
  end subroutine insertdroite_MF
!-------------------------------------------------------
 subroutine insertgauche_MF(f,lcf,ier)
    implicit none
    integer                 :: ier
    Type(chaine_MF), Pointer     :: lcf
    Type(mface)                  :: f
    Type(chaine_MF), Pointer     :: temp, locch
!cas du bout de liste droite ou gauche
    temp=>lcf%prec
       allocate (locch)
       locch%f=f
       locch%suiv=>lcf
       lcf%prec=>locch
       locch%prec=>temp
       temp%suiv=>locch
!modif
       lcf=>locch
       ier=0
       return
  end subroutine insertgauche_MF
!===========================   FIN DU MODULE   =========================
 END MODULE msh
!DOCFIN 



