!> \file dragging-vit_bil_CISM_gen_neem_mod.f90 !! Definition des zones de streams apres lecture d'une carte evaluee offline !< !> \namespace dragging_neem !! Definition zone stream. Input neem : critere sur vitesse de bilan et sur la 'courbure' !! \author ... !! \date ... !! @note Defini les zones de stream avec : !! @note - un critere sur la hauteur d'eau !! @note - un critere de fleuves actuels !! @note Attention il faut fonctionner en additionnant SIA et Stream pour éviter les pb !! de fleuves en pointille !! @note Used module !! @note - use module3D_phy !! @note - use param_phy_mod !< module dragging_plastic_LGM ! Defini les zones de stream avec : ! * un critere sur la hauteur d'eau ! * un critere de fleuves actuels ! attention il faut fonctionner en additionnant SIA et Stream pour éviter les pb ! de fleuves en pointille use module3d_phy use param_phy_mod use interface_input implicit none logical,dimension(nx,ny) :: fleuvemx !< fleuves sont les tableaux courants (dep. time) logical,dimension(nx,ny) :: fleuvemy logical,dimension(nx,ny) :: fleuve logical,dimension(nx,ny) :: cote real, dimension(nx,ny) :: Vcol_x !< uniquement pour compatibilite dans le spinup cat (aurel) real, dimension(nx,ny) :: Vcol_y !< uniquement pour compatibilite dans le spinup cat (aurel) real, dimension(nx,ny) :: Vsl_x !< pour compatibilite dragging beta (aurel) real, dimension(nx,ny) :: Vsl_y !< pour compatibilite dragging beta (aurel) real, dimension(nx,ny) :: work_tab ! pour lire les masques real, dimension(nx,ny) :: drag_centre ! real, dimension(nx,ny) :: uslid_centre ! pour calcul beta vitesse de glissement du pas précédent real, dimension(nx,ny) :: uslid_old ! vitesse de glissement du pas antepenultieme pour test convergence logical :: corr_def = .false. !< for deformation correction, pour compatibilite beta (aurel) character(len=100) :: masque_stream !< file masque streams autorises character(len=100) :: file_ncdf !< fichier netcdf issue des fichiers .dat real :: Tstream_lim !< limit for streaming (lower than melting point to account for spatial variability) real :: tob_prescrit !< valeur pour les ice streams real :: valmax integer :: imax,jmax integer :: i_moins1,i_plus1,j_moins1,j_plus1 integer :: lmax=20 integer :: idep,jdep,iloc,jloc integer :: nestim !< pour les estimateurs de convergence integer :: mm real :: uslid_diff real :: uslid_absdiff real :: uslid_moy real :: tostick ! pour la glace posee real :: tob_ile ! pour les iles contains !------------------------------------------------------------------------------- !> SUBROUTINE: init_dragging !! Cette routine fait l'initialisation du dragging. !< subroutine init_dragging ! Cette routine fait l'initialisation du dragging. implicit none integer :: bidon namelist/drag_plastic_LGM/Tstream_lim,tob_prescrit,betamax,masque_stream if (itracebug.eq.1) write(num_tracebug,*)' init_dragging avec dragging_plastic_LGM' ! formats pour les ecritures dans 42 428 format(A) ! lecture des parametres du run !-------------------------------------------------------------------- rewind(num_param) ! pour revenir au debut du fichier param_list.dat read(num_param,drag_plastic_LGM) write(num_rep_42,428) '!___________________________________________________________' write(num_rep_42,428) '! module dragging_plasticLGM ' write(num_rep_42,drag_plastic_LGM) write(num_rep_42,428) '! critere de temperature pour passage en stream' write(num_rep_42,428) '! si T > Tlim' write(num_rep_42,428) '! seulement pour les points cotiers' write(num_rep_42,428) '! tob_prescrit : (Pa) frottement maxi sous les streams ' write(num_rep_42,428) '! masque_stream : nom du fichier qui contient le masque de streams autorises' write(num_rep_42,*) tostick=betamax ! valeurs pour les points non flgzmx tob_ile=tob_prescrit/2. !moteurmax=1.e5 !moteurmax=0. !------------------------------------------------------------------- ! masque stream : lecture du masque masque_stream = trim(dirnameinp)//trim(masque_stream) !call lect_input(1,'z',1,work_tab,masque_stream,file_ncdf) call lect_input(1,'z',1,work_tab,masque_stream,file_ncdf) !call lect_ncfile('cry',work_tab,masque_stream) mstream(:,:) = int(work_tab(:,:)) !i = 154 !j = 123 !if (itracebug.eq.1) write(num_tracebug,*) 'init_dragging', masque_stream, work_tab(i,j) ! elargit les streams autorisés aux noeuds mineurs adjacents do j=1,ny do i=1,nx if (mstream(i,j).eq.1) then ! allowed ice streams mstream_mx(i:i+1,j) = 1 mstream_my(i,j:j+1) = 1 end if end do end do debug_3D(:,:,25) = mstream(:,:) ! coefficient permettant de modifier le basal drag. drag_centre(:,:) = 1 drag_mx(:,:) = 1. drag_my(:,:) = 1. do j=2,ny-1 do i=2,nx-1 uslid_centre(i,j) = (ux(i,j,nz)+ux(i-1,j,nz))**2+(uy(i,j,nz)+uy(i,j-1,nz))**2 uslid_centre(i,j) = uslid_centre(i,j)**0.5 uslid_centre(i,j) = max(uslid_centre(i,j),0.1) end do end do uslid_old(:,:) = uslid_centre(:,:) !!$call diffusiv() !!$call mix_SIA_L1 !!$call strain_rate() !!$ !!$do mm=1,50 !!$ write(6,*)'time, m',time, mm !!$ call diagnoshelf !!$ call mix_SIA_L1 !!$ call strain_rate() !!$end do return end subroutine init_dragging !________________________________________________________________________________ !> SUBROUTINE: dragging !! defini la localisation des streams et le frottement basal !< !------------------------------------------------------------------------- subroutine dragging ! defini la localisation des streams et le frottement basal ! 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(rog*Hmx(i,j)*sdx(i,j)).lt.toblim) ilemy(i,j)=ilemy(i,j).and.(abs(rog*Hmy(i,j)*sdy(i,j)).lt.toblim) end do end do fleuvemx(:,:)=.false. fleuvemy(:,:)=.false. fleuve(:,:)=.false. cote(:,:)=.false. cotemx(:,:)=.false. cotemy(:,:)=.false. ! detection des cotes do j=2,ny-1 do i=2,nx-1 if ((.not.flot(i,j)).and. & ((flot(i+1,j)).or.(flot(i,j+1)).or.(flot(i-1,j)).or.(flot(i,j-1)))) then cote(i,j)=.true. cotemx(i:i+1,j) = .true. cotemy(i,j:j+1) = .true. gzmx(i:i+1,j)=.true. gzmy(i,j:j+1)=.true. endif end do end do ! le calcul des fleuves se fait sur les mailles majeures ! normalement, les points cotiers sont deja tagges gzmx dans la routine flottab do j=2,ny-1 do i=2,nx-1 uslid_centre(i,j) = (ux(i,j,nz)+ux(i-1,j,nz))**2+(uy(i,j,nz)+uy(i,j-1,nz))**2 uslid_centre(i,j) = uslid_centre(i,j)**0.5 uslid_centre(i,j) = max(uslid_centre(i,j),0.1) end do end do ! calcul des estimateurs nestim = 0 uslid_diff = 0. uslid_absdiff = 0. uslid_moy = 0. do j=2,ny-1 do i=2,nx-1 if ((uslid_centre(i,j).gt..2).and.(.not.flot(i,j))) then nestim = nestim +1 uslid_diff = uslid_diff + (uslid_centre(i,j)-uslid_old(i,j)) uslid_absdiff = uslid_absdiff + abs( (uslid_centre(i,j)-uslid_old(i,j))) uslid_moy = uslid_moy + uslid_centre(i,j) end if end do end do if (nestim.gt.0) then write(6,*) 'conv uslid',nestim,uslid_diff/nestim,uslid_absdiff/nestim,uslid_moy/nestim end if uslid_old(:,:) = uslid_centre(:,:) do j=2,ny-1 do i=2,nx-1 if ( (T(i,j,nz)-Tpmp(i,j,nz).gt.Tstream_lim).and.(mstream(i,j).eq.1).and.(.not.flot(i,j))) then fleuve(i,j) = .true. drag_centre(i,j) = tob_prescrit/uslid_centre(i,j) drag_centre(i,j) = max(drag_centre(i,j),10.) fleuvemx(i:i+1,j)=.true. fleuvemy(i,j:j+1)=.true. gzmx(i:i+1,j)=.true. gzmy(i,j:j+1)=.true. else drag_centre(i,j) = tostick endif end do end do !i = 154 !j = 123 !write(num_tracebug,*)'dans drag_plastic', mstream(i,j),drag_centre(i,j),tob_prescrit,uslid_centre(i,j) where (flot(:,:)) drag_centre(:,:) = 0. end where call beta_c2stag(nx,ny,drag_mx,drag_my,drag_centre) ! redistribute on staggered grid where (cotemx(:,:)) ! point cotier forcement beta faible drag_mx(:,:)=min(drag_mx(:,:),1000.) end where where (cotemy(:,:)) ! point cotier forcement beta faible drag_my(:,:)=min(drag_my(:,:),1000.) end where betamx(:,:) = drag_mx(:,:) betamy(:,:) = drag_my(:,:) beta_centre(:,:) = drag_centre(:,:) where (fleuve(:,:)) debug_3D(:,:,1)=1 elsewhere (flot(:,:)) debug_3D(:,:,1)=2 elsewhere debug_3D(:,:,1)=0 endwhere where (cote(:,:)) debug_3D(:,:,2)=1 elsewhere debug_3D(:,:,2)=0 endwhere where (fleuvemx(:,:)) debug_3D(:,:,3)=1 elsewhere debug_3D(:,:,3)=0 endwhere where (flotmx(:,:)) debug_3D(:,:,3)=-1 endwhere !_____________________ where (fleuvemy(:,:)) debug_3D(:,:,4)=1 elsewhere debug_3D(:,:,4)=0 endwhere where (flotmy(:,:)) debug_3D(:,:,4)=-1 end where debug_3D(:,:,23)= abs(RO*G*HMX(:,:)*sdx(:,:)*1.e-5) debug_3D(:,:,24)= abs(RO*G*HMY(:,:)*sdy(:,:)*1.e-5) return end subroutine dragging !______________________________________________________________________________________ !>SUBROUTINE: beta_c2stag_min !! redistribute central values on staggered grid !< subroutine beta_c2stag(nxx,nyy,R_mx,R_my,R_centre) ! redistribute central values on staggered grid implicit none integer :: nxx,nyy ! dimension des tableaux real, dimension(nxx,nyy), intent(out) :: R_mx ! valeur du beta sur maille mx real, dimension(nxx,nyy), intent(out) :: R_my ! valeur du beta sur maille my real, dimension(nxx,nyy), intent(in) :: R_centre ! valeur du beta sur maille centree do j=2,nyy do i=2,nxx R_mx(i,j)= (R_centre(i-1,j)+R_centre(i,j))*0.5 R_my(i,j)= (R_centre(i,j-1)+R_centre(i,j))*0.5 end do end do end subroutine beta_c2stag end module dragging_plastic_LGM