Changeset 102


Ignore:
Timestamp:
05/05/17 16:44:01 (7 years ago)
Author:
dumas
Message:

Update closed water cycle | H = 1m on sea suppressed

Location:
trunk/SOURCES
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • trunk/SOURCES/3D-physique-gen_mod.f90

    r99 r102  
    222222  real :: diff_H 
    223223 
     224        real,dimension(nx,ny) :: H_beau_old,diff_H_2D,diff_H_water_bilan 
    224225 
    225226  !    ************* common des tableaux a deux dimensions ****** 
  • trunk/SOURCES/Greeneem_files/output_greeneem_mod-0.4.f90

    r99 r102  
    4848mask_cal(:,:,:)=.false. 
    4949mask_cal(:,:,1)=.true. 
     50sum_H_old = sum(H(:,:),mask=ice(:,:)==1) 
    5051do j=1,ny 
    5152   do i=1,nx 
     
    167168  enddo 
    168169   
     170         
     171   
    169172 
    170173!  write(num_ritz,905)   time,((isvol(kk),inp(kk),isvolf(kk),isvol_shelf(kk),inf(kk),isacc(kk),isabl(kk),    &   
    171174!       isablbord(kk),ablatot(kk),iscalv(kk),isbm(kk),Tjulymean(kk),Hmean(kk),Hmax(kk),nabl(kk),nacc(kk)),kk=1,nbregion) 
    172175 
    173   write(num_ritz,906) time,diff_H,water_bilan,sum(calv_dtt(:,:))/dtt, & 
    174                          sum(ablbord_dtt(:,:))/dtt, & 
    175                          sum(Bm_dtt(:,:))/dtt,sum(bmelt_dtt(:,:))/dtt   
     176  write(num_ritz,906) time,diff_H,water_bilan,sum(calv_dtt(2:nx-1,2:ny-1))/dtt, & 
     177                         sum(ablbord_dtt(2:nx-1,2:ny-1))/dtt, & 
     178                         sum(Bm_dtt(2:nx-1,2:ny-1))/dtt,sum(bmelt_dtt(2:nx-1,2:ny-1))/dtt,sum(ice(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1)  
    176179 
    177180905  format(f10.1,1x,2(e11.4, 1x, i5, 1x, e11.4, 1x, e11.4, 1x, i5, 9(1x, e12.5),2(1x,i5))) 
    178 906  format (f0.2,6(1x,e15.8)) 
     181906  format (f0.2,6(1x,e15.8),1x,i7) 
    179182          
    180183! pour verifier la zone groenland sur carte  
  • trunk/SOURCES/Temperature-routines/icetemp_mod.f90

    r73 r102  
    189189          Call temp_col(Nz,Nzm,Nn,T3d_new(i,j,:),Ct(i,j,:),& 
    190190               flot(i,j),H(i,j),Tbmer(i,j),Ibase(i,j),T(i,j,:),Tpmp(i,j,:),Aa,Bb,Cc,Rr,Hh,Ncond,Dee,Cm,& 
    191                Dzm,Phid(i,j),Ghf(i,j),Ifail1,Tbdot(i,j),DTT,Dou,Ts(i,j),Bmelt(i,j)) 
     191               Dzm,Phid(i,j),Ghf(i,j),Ifail1,Tbdot(i,j),DTT,Dou,Ts(i,j)) 
    192192 
    193193          ! Test Permettant D'Ecrire Les Variables Quand Il Y A Une Erreur 
  • trunk/SOURCES/Temperature-routines/temp_col.f90

    r73 r102  
    1515 
    1616Subroutine  Temp_col(Nz,Nzm,Nn,Newtempcol,Ctij,Flotij,Hij,Tbmerij,Ibase,Tij,Tpmpij,Aa,Bb,Cc,Rr,Hh,Ncond,Dee,Cm, & 
    17                                 Dzm,Phidij,Ghfij,Ifail1,Tbdotij,DTT,Dou,Tsij,Bmeltij) 
     17                                Dzm,Phidij,Ghfij,Ifail1,Tbdotij,DTT,Dou,Tsij) 
    1818 
    1919!  Use Icetemp_declar 
     
    4747  Real, intent(inout) :: Dou 
    4848  real, intent(in) :: Tsij        !< surface ice temperature  'o' 
    49   real, intent(out) :: Bmeltij    !< fonte basale 
    5049   
    5150   
     
    220219 
    221220     Tbdotij=(Newtempcol(Nz)-Tij(Nz))/Dtt 
    222      Bmeltij=0. 
    223221     Ibase=5 
    224222     Phidij=0. 
  • trunk/SOURCES/ablation_bord.f90

    r99 r102  
    3232! mise a jour du masque "ice" 
    3333where (flot(:,:))        ! points flottants, sera éventuellement réévalué dans flottab 
    34    where(H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:))) 
     34!cdc 1m  where(H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:))) 
     35   where(H(:,:).gt.0.) 
    3536      ice(:,:)=1 
    3637   elsewhere 
    3738      ice(:,:)=0 
    38       H(:,:)=max(1.,min(0.,(sealevel - Bsoc(:,:))*row/ro-0.01)) 
     39!cdc 1m      H(:,:)=max(1.,min(0.,(sealevel - Bsoc(:,:))*row/ro-0.01)) 
     40                        H(:,:)=0. 
    3941   endwhere 
    4042elsewhere                ! points posés 
     
    4749 
    4850 
     51 
     52 
     53 
    4954where (ice(:,:).eq.1)    !  quand glace : abl=bm-acc  
    5055   ablbord(:,:)=0. 
    5156endwhere 
     57!~ where (calv(:,:).LT.0.) ! il ne peut pas y avoir calving et ablbord sur un meme point 
     58!~       ablbord(:,:)=0. 
     59!~ endwhere 
    5260 
    5361! pour être au bord, il suffit que la somme des ice voisins soit > 0 
     
    6270         somm=ice(ip1,j)+ice(im1,j)+ice(i,jp1)+ice(i,jm1) 
    6371         if (somm.gt.0.and.(ablbord(i,j).LT.0.)) then                               ! voisins englaces 
    64 !         if (somm.gt.0.and.((ablbord(i,j).LT.0.).OR.(flot(ip1,j).OR.flot(im1,j).OR.flot(i,jp1).OR.flot(i,jm1)))) then   ! voisins englaces 
    65 !            ablbord(i,j)= ablbord(i,j) - (Acc(i,j) * dt)     !  Ablvrai=bmpdd-Acc-Hcalc/dt 
    66             ablbord(i,j)= ablbord(i,j) + ((Bmelt(i,j) - Bm(i,j)) * dt)  ! version Tof : ablbord = Hcons_mass - Bm - Bmelt 
     72!~ !         if (somm.gt.0.and.((ablbord(i,j).LT.0.).OR.(flot(ip1,j).OR.flot(im1,j).OR.flot(i,jp1).OR.flot(i,jm1)))) then   ! voisins englaces 
     73!~ !            ablbord(i,j)= ablbord(i,j) - (Acc(i,j) * dt)     !  Ablvrai=bmpdd-Acc-Hcalc/dt 
     74!~           if (ablbord(i,j).LT.0..and.front(i,j).LE.4.and.front(i,j).GE.0) then 
     75            ablbord(i,j)= max(0.,ablbord(i,j) + ((Bmelt(i,j) - Bm(i,j)) * dt))  ! version Tof : ablbord = Hcons_mass - Bm - Bmelt 
    6776 
    68          else 
    69 ! version cat            ablbord(i,j)=hdot(i,j)-acc(i,j)   ! on ne fond que ce qui est tombe + existant 
    70                                                 ablbord(i,j)=0. 
    71          endif 
    72  
     77          else 
     78!~ ! version cat            ablbord(i,j)=hdot(i,j)-acc(i,j)   ! on ne fond que ce qui est tombe + existant 
     79                                                ablbord(i,j)=0. 
     80          endif 
     81        else 
     82                      ablbord(i,j)=0. 
    7383      end if 
    7484   end do 
  • trunk/SOURCES/bilan_eau.f90

    r99 r102  
    2525real :: sum_H 
    2626 
     27real,dimension(nx,ny) :: bm_dt,bmelt_dt 
    2728 
    2829tot_water(:,:)=0. 
    2930 
     31bm_dt(:,:)=0. 
     32bmelt_dt(:,:)=0. 
     33 
     34bm_dt(2:nx-1,2:ny-1)=bm(2:nx-1,2:ny-1)*dt 
     35bmelt_dt(2:nx-1,2:ny-1)=bmelt(2:nx-1,2:ny-1)*dt 
     36!~ where (Bm(:,:).lt.0.) 
     37!~      bm_dt(:,:)=max(Bm(:,:)*dt,-H(:,:)) 
     38!~ elsewhere 
     39!~      bm_dt(:,:)=bm(:,:)*dt 
     40!~ endwhere 
     41 
     42!~ where (Bmelt(:,:).lt.0.) 
     43!~      bmelt_dt(:,:)=max(Bmelt(:,:)*dt,-H(:,:)) 
     44!~ elsewhere 
     45!~      bmelt_dt(:,:)=bmelt(:,:)*dt 
     46!~ endwhere 
    3047 
    3148 
    32 sum_H = sum(H(:,:),mask=ice(:,:)==1) 
     49! on fait la somme sur la grille en excluant les bords : 
     50 
     51sum_H = sum(H(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1) 
    3352diff_H = diff_H + (sum_H - sum_H_old) ! on calcul la variation de volume à chaque dt 
    3453 
    35 where (ice(:,:).eq.1) 
    36         Bm_dtt(:,:) = Bm_dtt(:,:) + Bm(:,:) * dt ! somme Bm sur dt 
    37         bmelt_dtt(:,:) = bmelt_dtt(:,:) + bmelt(:,:) * dt ! somme bmelt sur dt 
     54diff_H_2D(2:nx-1,2:ny-1)=H(2:nx-1,2:ny-1)-H_beau_old(2:nx-1,2:ny-1) 
     55 
     56where (ice(2:nx-1,2:ny-1).eq.1) 
     57        Bm_dtt(2:nx-1,2:ny-1) = Bm_dtt(2:nx-1,2:ny-1) + Bm_dt(2:nx-1,2:ny-1) !* dt ! somme Bm sur dt 
     58        bmelt_dtt(2:nx-1,2:ny-1) = bmelt_dtt(2:nx-1,2:ny-1) + bmelt_dt(2:nx-1,2:ny-1) ! * dt ! somme bmelt sur dt 
    3859endwhere 
    3960 
    4061if (isynchro.eq.1) then 
    4162! on raisonne en bilan annuel pour simplifier : 
    42 where (ice(:,:).eq.1) 
    43         tot_water(:,:) = (Bm_dtt(:,:) - Bmelt_dtt(:,:) + calv_dtt(:,:) - ablbord_dtt(:,:)) / dtt  
    44         elsewhere 
    45         tot_water(:,:) = (calv_dtt(:,:) - ablbord_dtt(:,:)) / dtt 
    46 endwhere 
    47  
     63!~      where (ice(2:nx-1,2:ny-1).eq.1) 
     64!~              tot_water(2:nx-1,2:ny-1) = (Bm_dtt(2:nx-1,2:ny-1) - Bmelt_dtt(2:nx-1,2:ny-1) + calv_dtt(2:nx-1,2:ny-1) - ablbord_dtt(2:nx-1,2:ny-1)) / dtt  
     65!~              elsewhere 
     66!~              tot_water(2:nx-1,2:ny-1) = (calv_dtt(2:nx-1,2:ny-1) - ablbord_dtt(2:nx-1,2:ny-1)) / dtt 
     67!~      endwhere 
     68!cdc pas besoin du test sur ice ici, il a ete fait avant (et le masque ice varie a chaque dt)    
     69  tot_water(2:nx-1,2:ny-1) = (Bm_dtt(2:nx-1,2:ny-1) - Bmelt_dtt(2:nx-1,2:ny-1) + calv_dtt(2:nx-1,2:ny-1) - ablbord_dtt(2:nx-1,2:ny-1)) / dtt 
     70    
    4871! bilan d'eau sur la grille : 
    49 water_bilan=sum(tot_water(:,:)) 
    50 diff_H = diff_H/dtt 
     72        water_bilan=sum(tot_water(:,:)) 
     73        diff_H = diff_H/dtt 
    5174 
    5275!999 format(f0.2,1x,e15.8,1x,i10,8(1x,e15.8)) 
    5376!       write(6,999),time,sum_H,count(ice(:,:)==1),diff_H,water_bilan,sum(calv_dtt(:,:))/dtt,sum(ablbord_dtt(:,:))/dtt,sum(bmelt_dtt(:,:),mask=ice(:,:)==1)/dtt,sum(bm(:,:),mask=ice(:,:)==1),sum(Bm_dtt(:,:))/dtt,sum(bmelt_dtt(:,:))/dtt 
     77        diff_H_water_bilan(2:nx-1,2:ny-1)=tot_water(2:nx-1,2:ny-1)-diff_H_2D(2:nx-1,2:ny-1) 
    5478 
    5579endif 
    5680sum_H_old=sum_H 
    5781 
     82H_beau_old(:,:)=H(:,:) 
     83 
    5884end subroutine bilan_eau 
  • trunk/SOURCES/calving_frange.f90

    r100 r102  
    397397                     calv(i,j)=-h(i,j)  
    398398!cdc                     H(i,j)=1.  
    399                      H(i,j)=min(1.,max(0.,(sealevel - Bsoc(i,j))*row/ro-0.01)) 
     399!cdc 1m                     H(i,j)=min(1.,max(0.,(sealevel - Bsoc(i,j))*row/ro-0.01)) 
     400                     H(i,j)=0. 
    400401                   endif   
    401402 
  • trunk/SOURCES/conserv-mass-adv-diff_sept2009_mod.f90

    r99 r102  
    153153   dt = dtmax 
    154154end if 
    155  
     155!~ print*,'testdiag/Maxdia, ',time, dt   
    156156! next_time ajuste le pas de temps pour faciliter la synchronisation  
    157157! avec le pas de temps long (dtt) 
     
    354354   if (testdiag/aux.lt.dt) then 
    355355      time=time-dt 
    356       dt=testdiag/aux*4. 
    357 if (itracebug.eq.1) call tracebug("aux tres petit,deuxieme appel a next_time dt*4") 
     356      dt=min(testdiag/aux,dtmax)   !cdc *4. 
     357!~       print*,'aux*4, ',time, dt, aux, testdiag/aux 
     358                        if (itracebug.eq.1) call tracebug("aux tres petit,deuxieme appel a next_time dt*4") 
    358359      call next_time(time,dt,dtt,dtmax,dtmin,isynchro,itracebug,num_tracebug) 
     360   else 
     361!~       print*,'testdiag/Maxdia, ',time, dt   
    359362   end if 
    360363end if 
     
    462465!!$endwhere 
    463466 
    464  
     467    !$OMP WORKSHARE 
     468!cdc    where (flot(:,:)) 
     469!cdc 1m       where (H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:))) 
     470!cdc      where (H(:,:).gt.0.) 
     471!cdc          ice(:,:)=1 
     472!cdc       elsewhere 
     473!cdc          ice(:,:)=0 
     474!cdc       end where 
     475!cdc    elsewhere 
     476       where (H(:,:).gt.0.) 
     477          ice(:,:)=1 
     478       elsewhere 
     479          ice(:,:)=0 
     480       end where 
     481!cdc    end where 
     482    !$OMP END WORKSHARE 
     483     
    465484! Appel a la routine de resolution resol_adv_diff_2D_vect 
    466485!------------------------------------------------------------ 
     
    484503! on met à 0 l'épaisseur sur les bords si nécessaire 
    485504do i=1,nx 
    486         if (H(i,1).gt.max(Hmin,Hmin+BM(i,1)-Bmelt(i,1))) then 
     505!cdc 1m if (H(i,1).gt.max(Hmin,Hmin+BM(i,1)-Bmelt(i,1))) then 
     506        if (H(i,1).gt.0.) then 
    487507                ice(i,1)=1 
    488508  else 
     
    490510    H(i,1)=0. 
    491511  endif   
    492         if (H(i,ny).gt.max(Hmin,Hmin+BM(i,ny)-Bmelt(i,ny))) then 
     512!cdc 1m if (H(i,ny).gt.max(Hmin,Hmin+BM(i,ny)-Bmelt(i,ny))) then 
     513        if (H(i,ny).gt.0.) then 
    493514                ice(i,ny)=1 
    494515  else 
     
    500521!$OMP DO 
    501522do j=1,ny 
    502         if (H(1,j).gt.max(Hmin,Hmin+BM(1,j)-Bmelt(1,j))) then 
     523!cdc 1m if (H(1,j).gt.max(Hmin,Hmin+BM(1,j)-Bmelt(1,j))) then 
     524        if (H(1,j).gt.0.) then 
    503525                ice(1,j)=1 
    504526  else 
     
    506528    H(1,j)=0. 
    507529  endif   
    508         if (H(nx,j).gt.max(Hmin,Hmin+BM(nx,j)-Bmelt(nx,j))) then 
     530!cdc 1m if (H(nx,j).gt.max(Hmin,Hmin+BM(nx,j)-Bmelt(nx,j))) then 
     531        if (H(nx,j).gt.0.) then 
    509532                ice(nx,j)=1 
    510533  else 
     
    549572!$OMP END WORKSHARE 
    550573 
    551 if (igrdline.ne.3) then 
    552    !$OMP WORKSHARE 
    553    Hdot(:,:)=min(Hdot(:,:),10.) 
    554    Hdot(:,:)=max(Hdot(:,:),-10.) 
    555    !$OMP END WORKSHARE 
    556 endif 
     574!~ if (igrdline.ne.3) then 
     575!~    !$OMP WORKSHARE 
     576!~    Hdot(:,:)=min(Hdot(:,:),10.) 
     577!~    Hdot(:,:)=max(Hdot(:,:),-10.) 
     578!~    !$OMP END WORKSHARE 
     579!~ endif 
    557580 
    558581!$OMP WORKSHARE 
     
    566589! en attendant un masque plus propre 
    567590 
    568 where(mk0(:,:).eq.-1) 
    569    H(:,:)=vieuxH(:,:) 
    570    Hdot(:,:)=0. 
    571 end where 
     591!~ where(mk0(:,:).eq.-1) 
     592!~    H(:,:)=vieuxH(:,:) 
     593!~    Hdot(:,:)=0. 
     594!~ end where 
    572595!$OMP END WORKSHARE 
    573596!$OMP END PARALLEL 
  • trunk/SOURCES/flottab2-0.7.f90

    r99 r102  
    224224!cdc             ice(i,j)=0 
    225225!cdc         H(i,j)=1. 
    226              H(i,j)=min(1.,max(0.,(sealevel - Bsoc(i,j))*row/ro-0.01)) 
     226!cdc  1m           H(i,j)=min(1.,max(0.,(sealevel - Bsoc(i,j))*row/ro-0.01)) 
     227             H(i,j)=0. 
    227228             surnet=H(i,j)*(1.-ro/row) 
    228229             S(i,j)=surnet+sealevel 
     
    257258!!$    end do 
    258259!!$ end do 
    259  
    260260 
    261261 
     
    592592    !$OMP WORKSHARE 
    593593    where (flot(:,:)) 
    594        where (H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:)))  
     594!cdc 1m       where (H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:))) 
     595                         where (H(:,:).gt.0.) 
    595596          ice(:,:)=1 
    596597       elsewhere  
     
    598599       end where 
    599600    elsewhere 
    600        where (H(:,:).gt.(0.))  
     601       where (H(:,:).gt.0.)  
    601602          ice(:,:)=1 
    602603       elsewhere  
  • trunk/SOURCES/steps_time_loop.f90

    r99 r102  
    5555 
    5656     call icethick3 
    57       
     57 
    5858     call flottab 
    5959 
    6060     call calving 
    61       
     61 
    6262     call ablation_bord 
    63       
     63 
    6464     call bilan_eau 
    6565      
     
    112112 
    113113     call shortoutput() 
    114 !cdc initialisation des variables bilan eau 
    115                  diff_H = 0. 
    116                  Bm_dtt(:,:) = 0. 
    117                  bmelt_dtt(:,:) = 0. 
    118114 
    119115  endif 
     
    171167   
    172168  ! end of outputs 
    173   !=================================================================== 
    174         ! remise a 0 de ablbord_dtt et calving_dtt 
    175  
    176         if (isynchro.eq.1) then 
    177                 calv_dtt(:,:)=0. 
    178                 ablbord_dtt(:,:)=0. 
    179         endif 
    180  
     169  if (isynchro.eq.1) then 
     170                 diff_H = 0. 
     171                 Bm_dtt(:,:) = 0. 
     172                 bmelt_dtt(:,:) = 0. 
     173           calv_dtt(:,:)=0. 
     174                 ablbord_dtt(:,:)=0. 
     175                 diff_H_2D(:,:)=0. 
     176        endif     
    181177  !==================================================================== 
    182178  !  
     
    245241  !------------------------------------------- 
    246242  call masque() 
    247   call flottab() 
     243 
    248244  if (iter_beta.gt.0)   call init_dragging 
    249245 
     
    256252 
    257253  isync_1: if (ISYNCHRO.eq.1) then          ! debut bloc appels tous les dtt     
    258  
    259254 
    260255     ! climatic forcing 
     
    461456              call  mix_SIA_L1   
    462457              call strain_rate() 
    463               call flottab()                           ! pour iteration dragging : flottab appelle dragging 
    464 !              write(444,*) time,test_iter_diag,m 
    465 !              if (test_iter_diag.lt.1.e-4) exit 
    466 !              write(555,*) time,test_iter_diag,m 
    467 !               write(222,*) time,test_iter_diag,m 
     458              if (iter_beta.eq.0) then 
     459                call dragging           
     460              endif 
    468461              if (test_iter_diag.lt.1.e-2) exit 
    469 !              write(333,*) time,test_iter_diag,m 
    470 !              if (test_iter_diag.lt.1.e-3) exit 
    471  
    472462           end do 
    473463        endif isync_2 
    474464 
    475  
    476465     call  mix_SIA_L1   
    477466 
Note: See TracChangeset for help on using the changeset viewer.