!> \file dragging_heino_mod.f90 !! Module qui defini la loi de glissement dans les experiences Heino !< !> \namespace sliding_heino !! Definition de la loi de glissement dans les experiences Heino !! \author Catherine !! \date ... !! @note Used module !! @note - use module3D_phy !< module sliding_heino ! Défini la loi de glissement dans les expériences Heino USE module3D_phy implicit none integer,dimension(nx,ny) :: mksedim !< masque des regions sediment-hard rock integer,dimension(nx,ny) :: mkxsedim !< décalle sur les demi mailles integer,dimension(nx,ny) :: mkysedim !< décalle sur les demi mailles real :: Cr !< coefficient glissement pour hard rock real :: Cs !< coefficient glissement pour sediment real :: coefmax real :: coefslid real :: coefdrag logical,dimension(nx,ny) :: gzmx_heino !< pour transporter la valeur dans flottab logical,dimension(nx,ny) :: gzmy_heino real :: longslope contains !------------------------------------------------------------------------------- !> SUBROUTINE: init_sliding !! Cette routine fait l'initialisation du glissement et est appelée juste apres !! la lecture du masque !< subroutine init_sliding if (itracebug.eq.1) call tracebug(' Heino: entree dans routine init_sliding') ! Cette routine fait l'initialisation du glissement et est appelée juste apres ! la lecture du masque ! ecriture dans le fichier parametres write(num_rep_42,*)'glissement Heino - initialisation masque' write(num_rep_42,*)'loi glissement=',loigliss write(num_rep_42,*)'----------------------------------------------------------' ! Le masque est égal au mk0 lu du début, mksedim(:,:)=mk0(:,:) where (mk0>1) mk0=1 !do j=ny, 1, -1 ! write(6,'(81i1)') (mk0(i,j), i=1,nx) !end do !write(6,*) mksedim(41,41) !do j=ny, 1, -1 ! write(6,'(81i1)') (mksedim(i,j), i=1,nx) !end do ! calcul du masque sur les demi mailles : frontiere -> sedim do j=2,ny-1 do i=2,nx-1 mkxsedim(i,j)=max(mksedim(i,j),mksedim(i-1,j)) mkysedim(i,j)=max(mksedim(i,j),mksedim(i,j-1)) end do end do Cr=1.e5 ! en a-1 coefficient hard rock (loi en puissance 3) Cs=500. ! en a-1 coefficient sediment (lineaire) write(num_rep_42,*)'coefficients: Cr=',Cr,' Cs=',Cs,' (a-1)' write(num_rep_42,*)'glissement seulement si les 2 voisins sont au point de fusion' if ((loigliss.eq.5).and.(Cs.gt.1.)) then coefdrag=rog/Cs Cs=0. write(num_rep_42,*)'passage en stream si zone sediment + fusion' write(num_rep_42,*)'coefficient de frottement dans la zone sediment=',coefdrag endif return end subroutine init_sliding !________________________________________________________________________________ !>SUBROUTINE: sliding !! Cette routine calcule le terme ddbx et ddby dans le cas Heino !< subroutine sliding ! cette routine calcule le terme ddbx et ddby dans le cas Heino ! ------------------------------------------------- ! 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(' Heino: entree dans routine sliding') coefmax=1. slidy: do j=2,NY do I=2,NX-1 pose_y: if (.not.flgzmy(i,j)) then if (HMY(I,J).lt.1.) then ! glace tres fine ddby(i,j)=0. ! else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).and.(T(i,j-1,nz).lt.Tpmp(i,j-1,nz))) then ! ddby(i,j)=0. ! pas de glissement seulement ! ! si les 2 voisins sont froids else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i,j-1,nz).lt.Tpmp(i,j-1,nz))) then ddby(i,j)=0. ! pas de glissement si au moins un seul ! des 2 voisins est froid else ! glissement coefslid=(hwater(i,j)+hwater(i,j-1))*0.5 coefslid=max(coefslid,0.) coefslid=min(coefslid/coefmax,1.) coefslid=1. if (mkysedim(i,j).eq.1) then ! loi hard rock ddby(i,j)=Cr*slope2my(i,j)*Hmy(i,j)*coefslid else if (mkysedim(i,j).eq.2) then ! loi sediment ddby(i,j)=Cs*Hmy(i,j)*coefslid ! if (abs(sdy(i,j)).gt.1.e-8) then ! j_1=max(j-1,1) ! j1=min(j+1,ny) ! longslope=(sdy(i,j_1)+sdy(i,j)+sdy(i,j1))/3. ! ddby(i,j)=ddby(i,j)/sdy(i,j)*longslope ! endif else ddby(i,j)=0. endif endif endif pose_y end do end do slidy slidx: do j=2,NY-1 do I=2,NX pose_x: if (.not.flgzmx(i,j)) then if (HMx(I,J).lt.1.) then ! glace tres fine ddbx(i,j)=0. ! else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).and.(T(i-1,j,nz).lt.Tpmp(i-1,j,nz))) then ! ddbx(i,j)=0. ! pas de glissement seulement ! ! si les 2 voisins sont froids else if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i-1,j,nz).lt.Tpmp(i-1,j,nz))) then ddbx(i,j)=0. ! pas de glissement si au moins un seul ! des 2 voisins est froid else ! glissement coefslid=(hwater(i,j)+hwater(i-1,j))*0.5 coefslid=max(coefslid,0.) coefslid=min(coefslid/coefmax,1.) coefslid=1. if (mkxsedim(i,j).eq.1) then ! loi hard rock ddbx(i,j)=Cr*slope2mx(i,j)*Hmx(i,j)*coefslid else if (mkxsedim(i,j).eq.2) then ! loi sediment ddbx(i,j)=Cs*Hmx(i,j)*coefslid ! if (abs(sdx(i,j)).gt.1.e-8) then ! i_1=max(i-1,1) ! i1=min(i+1,nx) ! longslope=(sdx(i_1,j)+sdx(i,j)+sdx(i1,j))/3. ! ddbx(i,j)=ddbx(i,j)/sdx(i,j)*longslope ! endif else ddbx(i,j)=0. endif endif endif pose_x end do end do slidx return end subroutine sliding !------------------------------------------------------------------------- !> SUBROUTINE: dragging !! Defini le basal dragging !< subroutine dragging ! dans la zone sediment : appele si loigliss=5, ! defini le basal drag flgzmx(:,:)=.false. flgzmy(:,:)=.false. betamx(:,:)=1.e5 betamy(:,:)=1.e5 froty: do j=2,ny do i=2,nx-1 if (mkysedim(i,j).eq.2) then ! loi sediment ! seulement si au dessus du point de fusion if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i,j-1,nz).lt.Tpmp(i,j-1,nz))) then ! base froide else betamy(:,:)=coefdrag flgzmy(i,j)=.true. endif endif end do end do froty frotx: do j=2,ny-1 do i=2,nx if (mkxsedim(i,j).eq.2) then ! loi sediment ! seulement si au dessus du point de fusion if ((T(i,j,nz).lt.Tpmp(i,j,nz)).or.(T(i-1,j,nz).lt.Tpmp(i-1,j,nz))) then ! base froide else betamx(:,:)=coefdrag flgzmx(i,j)=.true. endif endif end do end do frotx gzmx(:,:)=flgzmx(:,:) gzmy(:,:)=flgzmy(:,:) gzmx_heino(:,:)=gzmx(:,:) gzmy_heino(:,:)=gzmy(:,:) return end subroutine dragging END MODULE sliding_heino