Changeset 122 for branches/iLoveclim/SOURCES/calving_frange.f90

Timestamp:

06/29/17 16:27:48 (7 years ago)

Author:

aquiquet

Message:

Updating iLOVECLIM branch (calving / flottab / outputs) , working version but will be merged to trunk soon for water conservation

File:

• branches/iLoveclim/SOURCES/calving_frange.f90 (modified) (7 diffs)

branches/iLoveclim/SOURCES/calving_frange.f90

@@ -21 +21 @@
 real, dimension (nx,ny) :: hmhc  ! hauteur au dessus de la coupure
 real, dimension (nx,ny) :: bil_tot ! bilan surface et fond (bm-bmelt)
-real :: hcoup ! epaiseur de coupure au temps time
-real :: hcoup_max  ! epaisseur max
-real :: hcoup_min=100.  ! epaisseur min
+real, dimension (nx,ny) :: hcoup ! epaiseur de coupure au temps time
 real :: mass_total
+real :: hcoup_plateau  ! coupure points peu profonds
+real :: hcoup_abysses  ! coupure points ocean profond
+real :: prof_plateau  ! profondeur max des points peu profonds
+real :: prof_abysses  ! profondeur min des points ocean profond
 integer :: meth_hcoup  ! pour avoir hcoup dépendant du coefbmshelf
 integer :: ifrange        !  0 pas de traitement particulier pres du bord, 1 -> franges
@@ -39 +41 @@
 subroutine init_calving
 
-namelist/calving/Hcoup,ifrange,meth_hcoup
+namelist/calving/Hcoup_plateau,Hcoup_abysses,prof_plateau,prof_abysses,ifrange,meth_hcoup
 
 
@@ -59 +61 @@
 write(num_rep_42,428) '&calving              ! nom du bloc calving méthode Vincent '
 write(num_rep_42,*)
-write(num_rep_42,*)   'Hcoup      = ', Hcoup
+write(num_rep_42,*)   'Hcoup_plateau  = ', Hcoup_plateau
+write(num_rep_42,*)   'Hcoup_abysses  = ', Hcoup_abysses
+write(num_rep_42,*)   'prof_plateau  = ', prof_plateau
+write(num_rep_42,*)   'prof_abysses  = ', prof_abysses
 write(num_rep_42,*)   'ifrange    = ',ifrange
 write(num_rep_42,*)   'meth_hcoup = ', meth_hcoup
 write(num_rep_42,*)'/'
-                            
-write(num_rep_42,428) '! Hcoup epaisseur de coupure, le max si attache a coefbmshelf '
+
+write(num_rep_42,428) '! Hcoup epaisseurs de coupure pour les zones peu prodondes et profondes'
+write(num_rep_42,428) '! Hcoup_plateau<Hcoup_abysses && prof_plateau<prof_abysses'
+write(num_rep_42,428) '! prof profondeur delimitant les zones peu prodondes et profondes'
 write(num_rep_42,428) '! ifrange=0 -> pas de traitement particulier sur les bords'
 write(num_rep_42,428) '! ifrange=1 -> traitement de Vincent avec ice shelves frangeants partout'
 write(num_rep_42,428) '! ifrange=2 -> ice shelves frangeants seulement si bm-bmelt positif'
 write(num_rep_42,*)   '! meth_hcoup pour faire eventuellement varier Hcoup avec le climat'
-write(num_rep_42,*)   '! Hcoup_min=',Hcoup_min
+write(num_rep_42,*)   '! Hcoup_plateau=',Hcoup_plateau
+write(num_rep_42,*)   '! Hcoup_abysses=',Hcoup_abysses
+write(num_rep_42,*)   '! prof_plateau=',prof_plateau
+write(num_rep_42,*)   '! prof_abysses=',prof_abysses
 write(num_rep_42,*)
 
-Hcoup_max=Hcoup
+
+! afq -- coupure depend de la profondeur:
+!
+!  hcoup                       prof_abysses
+!   ^                               v 
+!   |                                _______ hcoup_abysses
+!   |                               /
+!   |                              /
+!   |                             /
+!   |                            /
+!   |      hcoup_plateau _______/
+!                               ^
+!                          prof_plateau     
+
+Hcoup(:,:) = min ( max(                                         &
+     (-(Bsoc0(:,:)-sealevel) - prof_plateau)/(prof_abysses-prof_plateau)    &
+     *(hcoup_abysses-hcoup_plateau)+hcoup_plateau               &
+     , hcoup_plateau), hcoup_abysses )
+
+if (meth_hcoup.eq.1) then
+   Hcoup(:,:)=coefbmshelf*Hcoup(:,:)
+   Hcoup(:,:)=min( max(Hcoup (:,:),Hcoup_plateau),Hcoup_abysses)
+else if (meth_hcoup.eq.2) then
+   Hcoup(:,:)=coefbmshelf*Hcoup(:,:)
+   Hcoup(:,:)=max(Hcoup(:,:),Hcoup_plateau)
+endif
+
 
 end subroutine init_calving
@@ -93 +129 @@
 
 ! coupure
-if (meth_hcoup.eq.0) then
-   Hcoup=Hcoup_max
-else if (meth_hcoup.eq.1) then
-   Hcoup=coefbmshelf*Hcoup_max
-   Hcoup=min(Hcoup,Hcoup_max)
-   Hcoup=max(Hcoup,Hcoup_min)
-
+Hcoup(:,:) = min ( max(                                         &
+     (-(Bsoc(:,:)-sealevel) - prof_plateau)/(prof_abysses-prof_plateau)    &
+     *(hcoup_abysses-hcoup_plateau)+hcoup_plateau               &
+     , hcoup_plateau), hcoup_abysses )
+
+if (meth_hcoup.eq.1) then
+   Hcoup(:,:)=coefbmshelf*Hcoup(:,:)
+   Hcoup(:,:)=min( max(Hcoup (:,:),Hcoup_plateau),Hcoup_abysses)
 else if (meth_hcoup.eq.2) then
-   Hcoup=coefbmshelf*Hcoup_max
-   Hcoup=max(Hcoup,Hcoup_min)
+   Hcoup(:,:)=coefbmshelf*Hcoup(:,:)
+   Hcoup(:,:)=max(Hcoup(:,:),Hcoup_plateau)
 endif
 
 ! hauteur au dessus de la coupure
-hmhc(:,:)=H(:,:)-Hcoup
+hmhc(:,:)=H(:,:)-Hcoup(:,:)
 
 ! coupure de l'ice shelf
@@ -117 +154 @@
 
 
-ifext:  if((flot(i,j)).and.(h(i,j).le.hcoup))  then
+ifext:  if((flot(i,j)).and.(h(i,j).le.hcoup(i,j)))  then
 ! ifext: pour les noeuds flottants avec h < hcoup 
 
@@ -373 +410 @@
    if (testipj) testH=testH + abs(uybar(i,j+1))*H(i,j+1)
    
-   if (testH.gt.Hcoup) then
+   if (testH.gt.Hcoup(i,j)) then
       testmij=.true.
    else
@@ -446 +483 @@
 ! mab: copy the calving in the right format for coupling                     
 ! calvin_GRIS =  volume flux per year (m3/a)
-                     calvin_GRIS(i,j)=calvin_GRIS(i,j)+(-calv(i,j)*dx*dy)
+                     calvingGRIS(i,j)=calvingGRIS(i,j)+(-calv(i,j)*dx*dy)
 !                     calvin_GRIS(i,j)=-calv(i,j)*40000.*40000.
 !                     write(15,*) 'calvheight:',h(i,j),calv(i,j),i,j,dtt,dt,calvin_GRIS(i,j)