[112] | 1 | !> \file bilan_eau.f90 |
---|
| 2 | !! Clacul du bilan d'eau sur la calotte |
---|
| 3 | !< |
---|
| 4 | |
---|
| 5 | |
---|
| 6 | !> MODULE: bilan_eau |
---|
| 7 | !! Calcule le bilan d'eau sur la grille GRISLI |
---|
| 8 | !! \author C. Dumas |
---|
| 9 | !! \date 20/02/2017 |
---|
| 10 | !! @note Permet de vérifier que le bilan d'eau est ok et qu'il n'y a pas de fuite. |
---|
| 11 | !> |
---|
| 12 | |
---|
| 13 | ! bilan d'eau : |
---|
| 14 | ! ce qui tombe : Acc |
---|
| 15 | ! ce qui sort : Bmelt, Abl, calving, ablbord |
---|
| 16 | ! ce qui est stocké : masse de glace |
---|
| 17 | ! le total de tous ces termes doit être constant (sur la grille, pas localement) |
---|
| 18 | |
---|
| 19 | module bilan_eau_mod |
---|
| 20 | |
---|
| 21 | |
---|
| 22 | use module3D_phy |
---|
| 23 | implicit none |
---|
| 24 | |
---|
| 25 | real :: sum_H |
---|
| 26 | |
---|
| 27 | |
---|
| 28 | real,dimension(nx,ny) :: tot_water !< bilan d'eau |
---|
| 29 | real,dimension(nx,ny) :: H_beau_old, diff_H_2D, diff_H_water_bilan |
---|
| 30 | real,dimension(nx,ny) :: Bm_dtt !< mass balance on ice points accumulated during dtt |
---|
| 31 | real,dimension(nx,ny) :: Bmelt_dtt !< basal melting on ice points accumulated during dtt |
---|
| 32 | real,dimension(nx,ny) :: calv_dtt !< calving sur dtt (pour calcul bilan d'eau) |
---|
[136] | 33 | real,dimension(nx,ny) :: archimtab ! point pose si > 0 |
---|
| 34 | real,dimension(nx,ny) :: grline_dtt ! grounding line flux during dtt |
---|
[112] | 35 | |
---|
| 36 | real :: sum_H_old |
---|
| 37 | real :: diff_H |
---|
[136] | 38 | real :: alpha_flot ! ro/row |
---|
[112] | 39 | |
---|
| 40 | real,dimension(nx,ny) :: bm_dt,bmelt_dt |
---|
| 41 | |
---|
| 42 | |
---|
| 43 | contains |
---|
| 44 | |
---|
| 45 | subroutine init_bilan_eau |
---|
| 46 | ! initialisation des variables |
---|
| 47 | diff_H=0. |
---|
| 48 | sum_H_old = sum(H(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1) |
---|
| 49 | tot_water(:,:)=0. |
---|
| 50 | bm_dt(:,:)=0. |
---|
| 51 | bmelt_dt(:,:)=0. |
---|
[136] | 52 | alpha_flot=ro/row |
---|
[112] | 53 | end subroutine init_bilan_eau |
---|
| 54 | |
---|
| 55 | |
---|
| 56 | |
---|
| 57 | |
---|
| 58 | subroutine bilan_eau |
---|
| 59 | |
---|
| 60 | tot_water(:,:)=0. |
---|
| 61 | |
---|
| 62 | bm_dt(:,:)=0. |
---|
| 63 | bmelt_dt(:,:)=0. |
---|
| 64 | |
---|
| 65 | bm_dt(2:nx-1,2:ny-1)=bm(2:nx-1,2:ny-1)*dt |
---|
| 66 | bmelt_dt(2:nx-1,2:ny-1)=bmelt(2:nx-1,2:ny-1)*dt |
---|
| 67 | !~ where (Bm(:,:).lt.0.) |
---|
| 68 | !~ bm_dt(:,:)=max(Bm(:,:)*dt,-H(:,:)) |
---|
| 69 | !~ elsewhere |
---|
| 70 | !~ bm_dt(:,:)=bm(:,:)*dt |
---|
| 71 | !~ endwhere |
---|
| 72 | |
---|
| 73 | !~ where (Bmelt(:,:).lt.0.) |
---|
| 74 | !~ bmelt_dt(:,:)=max(Bmelt(:,:)*dt,-H(:,:)) |
---|
| 75 | !~ elsewhere |
---|
| 76 | !~ bmelt_dt(:,:)=bmelt(:,:)*dt |
---|
| 77 | !~ endwhere |
---|
| 78 | |
---|
| 79 | |
---|
| 80 | ! on fait la somme sur la grille en excluant les bords : |
---|
| 81 | |
---|
| 82 | sum_H = sum(H(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1) |
---|
| 83 | |
---|
| 84 | diff_H = diff_H + (sum_H - sum_H_old) ! on calcul la variation de volume a chaque dt |
---|
| 85 | |
---|
| 86 | diff_H_2D(2:nx-1,2:ny-1)=H(2:nx-1,2:ny-1)-H_beau_old(2:nx-1,2:ny-1) |
---|
| 87 | |
---|
| 88 | where (ice(2:nx-1,2:ny-1).eq.1) |
---|
| 89 | 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 |
---|
| 90 | 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 |
---|
| 91 | endwhere |
---|
| 92 | |
---|
[136] | 93 | archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot - sealevel |
---|
| 94 | gr_line(:,:)=0 |
---|
| 95 | do j=1,ny |
---|
| 96 | do i=1,nx |
---|
[145] | 97 | !afq if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel)) then ! grounded with ice |
---|
| 98 | if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.)) then ! grounded with ice |
---|
[136] | 99 | if (archimtab(i-1,j).LT.0..and.Uxbar(i,j).LT.0..and..not.flot_marais(i-1,j)) gr_line(i,j)=1 |
---|
| 100 | if (archimtab(i+1,j).LT.0..and.Uxbar(i+1,j).GT.0..and..not.flot_marais(i+1,j)) gr_line(i,j)=1 |
---|
| 101 | if (archimtab(i,j-1).LT.0..and.Uybar(i,j).LT.0..and..not.flot_marais(i,j-1)) gr_line(i,j)=1 |
---|
| 102 | if (archimtab(i,j+1).LT.0..and.Uybar(i,j+1).GT.0..and..not.flot_marais(i,j+1)) gr_line(i,j)=1 |
---|
| 103 | endif |
---|
| 104 | enddo |
---|
| 105 | enddo |
---|
| 106 | |
---|
| 107 | where (gr_line(:,:).eq.1) |
---|
| 108 | !~ grline_dtt(:,:)= (((uxbar(:,:)+eoshift(uxbar(:,:),shift=1,boundary=0.,dim=1))**2+ & |
---|
| 109 | !~ (uybar(:,:)+eoshift(uybar(:,:),shift=1,boundary=0.,dim=2))**2)**0.5)*0.5 & |
---|
[145] | 110 | !~ *H(:,:) + grline_dtt(:,:) |
---|
| 111 | grline_dtt(:,:)= - sqrt( & |
---|
| 112 | ( (uxbar(:,:)+eoshift(uxbar(:,:),shift=1,boundary=0.,dim=1))*dy/2. )**2 + & |
---|
| 113 | ( (uybar(:,:)+eoshift(uybar(:,:),shift=1,boundary=0.,dim=2))*dx/2. )**2 ) & |
---|
| 114 | * H(:,:) * dt + grline_dtt(:,:) |
---|
[136] | 115 | endwhere |
---|
| 116 | |
---|
| 117 | |
---|
[112] | 118 | if (isynchro.eq.1) then |
---|
| 119 | ! on raisonne en bilan annuel pour simplifier : |
---|
| 120 | !~ where (ice(2:nx-1,2:ny-1).eq.1) |
---|
| 121 | !~ 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 |
---|
| 122 | !~ elsewhere |
---|
| 123 | !~ 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 |
---|
| 124 | !~ endwhere |
---|
| 125 | !cdc pas besoin du test sur ice ici, il a ete fait avant (et le masque ice varie a chaque dt) |
---|
| 126 | 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 |
---|
| 127 | |
---|
| 128 | ! bilan d'eau sur la grille : |
---|
| 129 | water_bilan=sum(tot_water(:,:)) |
---|
| 130 | diff_H = diff_H/dtt |
---|
| 131 | |
---|
| 132 | !999 format(f0.2,1x,e15.8,1x,i10,8(1x,e15.8)) |
---|
| 133 | ! 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 |
---|
| 134 | 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) |
---|
| 135 | |
---|
| 136 | endif |
---|
| 137 | sum_H_old=sum_H |
---|
| 138 | |
---|
| 139 | H_beau_old(:,:)=H(:,:) |
---|
| 140 | |
---|
| 141 | end subroutine bilan_eau |
---|
| 142 | |
---|
| 143 | end module bilan_eau_mod |
---|