source: trunk/SOURCES/Ant40_files/dragging_stream_impose_vitbil_mod.f90 @ 93

!> \file dragging_stream_impose_vitbil_mod.f90
!! Module qui defini les zones de stream d'apres les vitesses de bilan
!<

!> \namespace dragging_vitbil
!! Module qui defini les zones de stream d'apres les vitesses de bilan
!! \author ...
!! \date ...
!! @note Used module
!! @note   - use module3D_phy
!<
module dragging_vitbil
! Défini les zones de stream d'apres les vitesses de bilan

USE module3D_phy

implicit none
character(len=80) :: fil_vitbil
integer,dimension(nx,ny) :: uxbil
integer,dimension(nx,ny) :: uybil
integer :: icfcote
real :: cfcote

contains

!-------------------------------------------------------------------------------
subroutine init_dragging      ! Cette routine fait l'initialisation du dragging.

implicit none

integer :: ii,jj
integer,parameter :: ntotbil=50
integer :: kux,kuy,ncolbil
integer,dimension(ntotbil) :: index
real,dimension(ntotbil) :: variables
real :: seuil_stream

if (itracebug.eq.1)  call tracebug(' Dragging avec vitbil')

fil_vitbil='./Balance_velocities/depart-Camb-004+c01/vitbil_iter_100000'
!fil_vitbil='file_vitbil-impose'

betamax=200.   ! frottement maxi  (utilise dans remplimat-ant-0.4-40km 
              ! pour le calcul de frotmx

toblim=0.25e5 ! test pour les iles
seuil_stream=50
cf=-10.e-5    ! coefficient de la loi de frottement fonction
              ! de la vitesse et pression effective
icfcote=0     ! traitement particulier pour la cote
cfcote=0      ! valeur de cf pour la cote

!    ecriture dans le fichier parametres
     write(num_rep_42,*)'dragging vitbil-  initialisation '
     write(num_rep_42,*)'vitesses venant du fichier',fil_vitbil
     write(num_rep_42,*)'seuil stream=',seuil_stream,'(m)'
     write(num_rep_42,*)' frottement maxi, betamax=',betamax
     write(num_rep_42,fmt='(a,e10.3)') ' coefficient de frottement: cf=',cf

     if (icfcote.eq.1) then
        write(num_rep_42,*)'traitement particulier cote,  cfcote=',cfcote
     else
        write(num_rep_42,*)'pas de traitement particulier cote'
     endif


     write(num_rep_42,*)'----------------------------------------------------------'

! lecture des masques dans un fichier hz ------------------------

open(12,file=trim(fil_vitbil))
read(12,*)   ! lecture ligne 1
read(12,*)  i,ncolbil ! lecture ligne 2, 2eme position, nb de colonnes
read(12,*) (index(k),k=1,ncolbil)   

do k=1,ntotbil
   if (index(k).eq.32) kux=k   ! colonne ux
   if (index(k).eq.33) kuy=k   ! colonne uy
end do

read(12,*)   ! lecture ligne titres variables
do j=1,ny
   do i=1,nx
      read(12,*)(variables(k),k=1,ncolbil)
      uxbil(i,j)=variables(kux)
      uybil(i,j)=variables(kuy)
   end do
end do

close(12)

!-------------------------------------------------------------------
! masque stream

mstream_mx(:,:)=0
mstream_my(:,:)=0

do j=2,ny-1
   do i=2,nx-1
      if (abs(uxbil(i,j)).ge.seuil_stream) then
         mstream_mx(i,j)=1
      endif

      if (abs(uybil(i,j)).ge.seuil_stream) then
         mstream_my(i,j)=1
      endif
   end do
end do

open(122,file='masque-vitbil')
do j=1,ny
   do i=1,nx
      write(122,*) i,j,mstream_mx(i,j),mstream_my(i,j)
   end do
end do

! coefficient permettant de modifier le basal drag.
drag_mx(:,:)=1.
drag_my(:,:)=1.

! idee future : pour les points limite mer : drag_mx=0.

  
return
end subroutine init_dragging
!________________________________________________________________________________


!-------------------------------------------------------------------------
subroutine dragging   ! defini la localisation des streams et le frottement basal

implicit none
integer :: ii,jj

!         les iles n'ont pas de condition neff mais ont la condition toblim
!         (ce bloc etait avant dans flottab)

do j=2,ny
   do i=2,nx
      ilemx(i,j)=ilemx(i,j).and.(abs(tobmx(i,j)).lt.toblim)
      ilemy(i,j)=ilemy(i,j).and.(abs(tobmy(i,j)).lt.toblim)
   end do
end do


gzm_x: do j=2,ny
   do i=2,nx-1
   
! gzmx pour les hauteurs d'eau supérieures à un seuil
! attention, il semble qu'il y avait une erreur dans la precedente
! formulation car le test etait seulement sur hwater(i,j).

      gzmx(i,j)=gzmx(i,j).or. & 
           ((mstream_mx(i,j).eq.1) &
           .and..not.flotmx(i,j))


! calcul du frottement basal (ce bloc etait avant dans neffect)
      if (gzmx(i,j)) then                          ! stream
         betamx(i,j)=cf*neffmx(i,j) 

! si limite region flottante        
       if (icfcote.eq.1) then
         if ((flot(i,j).and..not.flot(i-1,j)) &                !    F x  P
             .or.(.not.flot(i,j).and.flot(i-1,j))) then        !    P x  F
            betamx(i,j)=cfcote*neffmx(i,j)
         endif
      endif

! limite sur deux points de grille
      if (icfcote.eq.2) then
         ii=max(1,i-2)
         if (flot(i,j).or.flot(i+1,j) &               !  P x F  ou P x P    F
              .or. flot(i-1,j).or.flot(ii,j)) then    !  F x P  ou F  P x P
            betamx(i,j)=cfcote*neffmx(i,j)
         end if
      end if




      else if (flotmx(i,j).or.ilemx(i,j)) then     ! flottant ou ile
         betamx(i,j)=0.
      else                                         ! grounded, SIA
         betamx(i,j)=1.e5                          ! frottement 1 bar
      endif      

   end do
end do gzm_x

gzm_y: do j=2,ny-1
   do i=2,nx
   
! gzmy pour les hauteurs d'eau supérieures à un seuil
! attention, il semble qu'il y avait une erreur dans la precedente
! formulation car le test etait seulement sur hwater(i,j).

      gzmy(i,j)=gzmy(i,j).or. & 
           ((mstream_my(i,j).eq.1) &
           .and..not.flotmy(i,j))


 ! calcul du frottement basal (ce bloc etait avant dans neffect)
      if (gzmy(i,j)) then                ! stream
         betamy(i,j)=cf*neffmy(i,j) 

! si limite region flottante        
       if (icfcote.eq.1) then
         if ((flot(i,j).and..not.flot(i,j-1)) &                !    F    P
             .or.(.not.flot(i,j).and.flot(i,j-1))) then        !    P    F
            betamy(i,j)=cfcote*neffmy(i,j) 
         endif
      endif

! limite sur deux points de grille
      if (icfcote.eq.2) then
         jj=max(1,j-2)
         if (flot(i,j).or.flot(i,j+1) &               !  P x F  ou P x P    F
              .or. flot(i,j-1).or.flot(i,jj)) then    !  F x P  ou F  P x P
            betamy(i,j)=cfcote*neffmy(i,j)
         end if
      end if


      else if (flotmy(i,j).or.ilemy(i,j)) then     ! flottant ou ile
         betamy(i,j)=0.
      else                               ! grounded, SIA
         betamy(i,j)=1.e5                ! frottement 1 bar
      endif      

   end do
end do gzm_y

return
end subroutine dragging

end module dragging_vitbil