!> \file sliding_Bindshadler_mod.f90
!! Module pour definir la loi de glissement Binshadler 
!<

!> \namespace sliding_Bindschadler
!! Module pour definir la loi de glissement Binshadler 
!! \author CatRitz et al.
!! \date 2001
!! @note Used module
!! @note   - use module3D_phy
!! @note   - use param_phy_mod
!<

module sliding_Bindschadler
! loi de glissement Binshadler modifie Ritz (Ritz et al. 2001)

USE module3D_phy
USE param_phy_mod

implicit none
real :: coefneff=1
real ::  kweert                       ! coefficient de la loi de Weertman simple 

contains
!-------------------------------------------------------------------------------
subroutine init_sliding

implicit none

namelist/slid_bindsh/kweert,loigliss,coefbmax

! loigliss est declare dans 3D-... pourrait etre dans les divers modules de sliding
! coefbmax est declare dans 3D-.. facteur de normalisation pour l'influence de l'eau


if (itracebug.eq.1)  call tracebug(' Bindschadler entree dans routine init_sliding')

rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
read(num_param,slid_bindsh)

write(num_rep_42,*)'!___________________________________________________________' 
write(num_rep_42,*)'! glissement                   module  sliding_Bindschadler '
write(num_rep_42,*)
write(num_rep_42,slid_bindsh)
write(num_rep_42,*)'! kweert : coefficent, loigliss le type de loi'
write(num_rep_42,*)'! coefbmax : facteur de normalisation pour influence eau'
write(num_rep_42,*)'!___________________________________________________________' 


coefbmax=min(coefbmax,hwatermax)

    
! tableaux  pour le glissement
mslid_mx(:,:)=1
mslid_my(:,:)=1
slid_mx(:,:)=1.
slid_my(:,:)=1.
 
    
return
end subroutine init_sliding


!________________________________________________________________________________
subroutine sliding

implicit none
real :: neffmin

! cette routine calcule le terme ddbx et ddby dans le cas Bindshadler modifie
! loigliss=2 

!      -------------------------------------------------
!      GLISSEMENT
!      DDBX et DDBY termes dus au glissement 
!      relation avec la vitesse de glissement UXB et UYB 
!      UXB=-DDBX*SDX et UYB=-DDBY*SDY
!      -------------------------------------------------


if (itracebug.eq.1)  call tracebug(' Entree dans routine sliding Bindshadler')

!coefmax=1.
ddby(:,:) = 0.
ddbx(:,:) = 0.

slidy: do j=2,NY
   do I=2,NX-1

!      write(6,*) 'boucle y, i,j',i,j

pose_y:  if (.not.flgzmy(i,j)) then

finey:   if (HMY(I,J).lt.1.) then        ! glace tres fine
      ddby(i,j)=0.

   else

!          le glissement depend de l'eau basale. 
!          ------------------------------------- 
!        la pression effective a ete calculee dans Neffect
!        Neffmy sur la demi maille y
!        coefmybmelt est calculee dans neffect et correspond a 
!        la hauteur d'eau sur les demi mailles. 


      sec_y: if ((COEFMYBMELT(i,j).le.0.).and.(.not.gzmy(i,j))) then
            DDBY(I,J)=0.
!            write(6,*) '2',i,j,COEFMYBMELT(i,j)
         else 
            neffmin=max(Neffmy(i,j),0.01)
            neffmin=max(Neffmy(i,j),1.e6)
            coefneff=min(abs(rog*hmy(i,j)/neffmin),1/0.0001)
            coefneff=min(coefneff,50.)
            coefneff=min(coefneff,5.)

            if ((Neffmy(i,j).gt.rog*2000.).and. &
                 (COEFMYBMELT(i,j).lt.coefbmax)) then

!            introduction d'un coefficient supplementaire
!            pour avoir une limitation liee a la fusion basale

               COEFNeff=coefNeff*min(1.,COEFMYBMELT(i,j)/COEFBMAX)
!               write(6,*) '3',i,j,coefneff
            endif
       
            ddby(i,j)=kweert*coefneff*(rog*hmy(i,j))**2.*sqrt(slope2my(i,j))

         end if sec_y
      end if finey
    end if pose_y

   end do
end do slidy



slidx: do j=2,NY-1
   do I=2,NX

!      write(6,*) 'boucle x, i,j',i,j

pose_x:  if (.not.flgzmx(i,j)) then

finex:   if (HMX(I,J).lt.1.) then        ! glace tres fine
      ddbx(i,j)=0.

   else

!          le glissement depend de l'eau basale. 
!          ------------------------------------- 
!        la pression effective a ete calculee dans Neffect
!        Neffmy sur la demi maille y
!        coefmybmelt est calculee dans neffect et correspond a 
!        la hauteur d'eau sur les demi mailles. 


      sec_x: if ((COEFMXBMELT(i,j).le.0.).and.(.not.gzmx(i,j))) then
            DDBX(I,J)=0.

         else 

            neffmin=max(Neffmx(i,j),0.01)
            neffmin=max(Neffmx(i,j),1.e6)

            coefneff=min(abs(rog*hmx(i,j)/neffmin),1/0.0001)
            coefneff=min(coefneff,50.)
            coefneff=min(coefneff,5.)

            if ((Neffmx(i,j).gt.rog*2000.).and. &
                 (COEFMXBMELT(i,j).lt.coefbmax)) then

!            introduction d'un coefficient supplementaire
!            pour avoir une limitation liee a la fusion basale

               COEFNeff=coefNeff*min(1.,COEFMXBMELT(i,j)/COEFBMAX)

            end if
        
            ddbx(i,j)=kweert*coefneff*(rog*hmx(i,j))**2.*sqrt(slope2mx(i,j))

         end if sec_x
      end if finex

   end if pose_x

   end do
end do slidx



return
end subroutine sliding

end module sliding_bindschadler