1 | !> \file bilan_flux_output_mod.f90 |
---|
2 | !! Clacul du bilan des flux : sorties initMIP |
---|
3 | !< |
---|
4 | |
---|
5 | |
---|
6 | !> MODULE: bilan_flux |
---|
7 | !! Calcule le bilan des flux : sortie initMIP |
---|
8 | !! \author C. Dumas |
---|
9 | !! \date 10/10/2017 |
---|
10 | !! @note Calcul des flux output initMIP |
---|
11 | !> |
---|
12 | |
---|
13 | module bilan_flux_mod |
---|
14 | |
---|
15 | |
---|
16 | use module3D_phy |
---|
17 | implicit none |
---|
18 | |
---|
19 | real :: dtoutflux ! pas de temps des sorties flux |
---|
20 | real :: dtsum ! somme des pas de temps (pour verification) |
---|
21 | !real :: alpha_flot ! ro/row |
---|
22 | |
---|
23 | ! variables flux : |
---|
24 | real,dimension(nx,ny) :: bm_flux ! surface smb on ice mask |
---|
25 | real,dimension(nx,ny) :: bmeltgr_flux ! bmelt on ice grounded mask |
---|
26 | real,dimension(nx,ny) :: bmeltfl_flux ! bmelt on floating ice mask |
---|
27 | double precision,dimension(nx,ny) :: dhdt_flux ! dhdt on ice mask |
---|
28 | real,dimension(nx,ny) :: calving_flux ! calving on ice mask |
---|
29 | real,dimension(nx,ny) :: grline_flux ! grounding line flux |
---|
30 | real,dimension(nx,ny) :: calv_abl_bord_flux ! calving + abl_bord sur points en contact avec la mer |
---|
31 | |
---|
32 | double precision,dimension(nx,ny) :: vieuxH_flux ! H at previous time step |
---|
33 | !real,dimension(nx,ny) :: archimtab ! point pose si > 0 |
---|
34 | !integer,dimension(nx,ny) :: gr_line ! point grounding line (pose sous niveau de la mer avec point flottant a cote) |
---|
35 | |
---|
36 | real,parameter :: ice_density=910. !< densite de la glace pour conversion en masse |
---|
37 | |
---|
38 | contains |
---|
39 | |
---|
40 | subroutine init_bilan_flux |
---|
41 | ! initialisation des variables |
---|
42 | ! diff_H=0. |
---|
43 | ! sum_H_old = sum(H(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1) |
---|
44 | ! tot_water(:,:)=0. |
---|
45 | ! bm_dt(:,:)=0. |
---|
46 | ! bmelt_dt(:,:)=0. |
---|
47 | |
---|
48 | dtoutflux=1. ! initMIP sorties flux tous les 1 an |
---|
49 | dtsum=0. |
---|
50 | bm_flux(:,:)=0. |
---|
51 | bmeltgr_flux(:,:)=0. |
---|
52 | bmeltfl_flux(:,:)=0. |
---|
53 | dhdt_flux(:,:)=0. |
---|
54 | calving_flux(:,:)=0. |
---|
55 | grline_flux(:,:)=0. |
---|
56 | calv_abl_bord_flux(:,:)=0. |
---|
57 | vieuxH_flux(:,:) = H(:,:) |
---|
58 | |
---|
59 | end subroutine init_bilan_flux |
---|
60 | |
---|
61 | |
---|
62 | |
---|
63 | |
---|
64 | subroutine bilan_flux_output |
---|
65 | |
---|
66 | dtsum = dtsum + dt |
---|
67 | |
---|
68 | |
---|
69 | ! points englacés |
---|
70 | where (ice(2:nx-1,2:ny-1).eq.1) |
---|
71 | bm_flux(2:nx-1,2:ny-1) = bm_flux(2:nx-1,2:ny-1) + bm(2:nx-1,2:ny-1)*dt ! somme Bm sur dt |
---|
72 | where (flot(2:nx-1,2:ny-1)) |
---|
73 | bmeltfl_flux(2:nx-1,2:ny-1) = bmeltfl_flux(2:nx-1,2:ny-1) + bmelt(2:nx-1,2:ny-1)*dt ! + ablbord(2:nx-1,2:ny-1) ! somme bmelt flot sur dt |
---|
74 | elsewhere |
---|
75 | bmeltgr_flux(2:nx-1,2:ny-1) = bmeltgr_flux(2:nx-1,2:ny-1) + bmelt(2:nx-1,2:ny-1)*dt ! + ablbord(2:nx-1,2:ny-1) ! somme bmelt grounded sur dt |
---|
76 | endwhere |
---|
77 | |
---|
78 | dhdt_flux(2:nx-1,2:ny-1) = dhdt_flux(2:nx-1,2:ny-1) + (H(2:nx-1,2:ny-1)-vieuxH_flux(2:nx-1,2:ny-1)) ! somme dhdt sur dt |
---|
79 | endwhere |
---|
80 | |
---|
81 | ! points non englacés et flottant |
---|
82 | where (ice(2:nx-1,2:ny-1).eq.0.and.flot(2:nx-1,2:ny-1)) |
---|
83 | calv_abl_bord_flux(2:nx-1,2:ny-1) = calv_abl_bord_flux(2:nx-1,2:ny-1) - ablbord(2:nx-1,2:ny-1) |
---|
84 | endwhere |
---|
85 | |
---|
86 | calving_flux(2:nx-1,2:ny-1) = calving_flux(2:nx-1,2:ny-1) + calv(2:nx-1,2:ny-1) |
---|
87 | |
---|
88 | where (gr_line(:,:).eq.1) |
---|
89 | !~ grline_flux(:,:)= (((uxbar(:,:)+eoshift(uxbar(:,:),shift=1,boundary=0.,dim=1))**2+ & |
---|
90 | !~ (uybar(:,:)+eoshift(uybar(:,:),shift=1,boundary=0.,dim=2))**2)**0.5)*0.5 & |
---|
91 | !~ *H(:,:) + grline_flux(:,:) |
---|
92 | grline_flux(:,:)= - sqrt( & |
---|
93 | ( (uxbar(:,:)+eoshift(uxbar(:,:),shift=1,boundary=0.,dim=1))*dy/2. )**2 + & |
---|
94 | ( (uybar(:,:)+eoshift(uybar(:,:),shift=1,boundary=0.,dim=2))*dx/2. )**2 ) & |
---|
95 | * H(:,:) * dt /(dx*dy) + grline_flux(:,:) |
---|
96 | endwhere |
---|
97 | |
---|
98 | |
---|
99 | ! calcul de la moyenne des flux sur dtoutflux |
---|
100 | !if ((mod(abs(time),dtoutflux).lt.dtmin).and.(abs(dtsum-dtoutflux).lt.dtmin)) then |
---|
101 | if (mod(abs(time),dtoutflux).lt.dtmin) then |
---|
102 | bm_flux(:,:)=(bm_flux(:,:)/dtsum)*ice_density/secyear |
---|
103 | bmeltgr_flux(:,:)=(bmeltgr_flux(:,:)/dtsum)*ice_density/secyear |
---|
104 | bmeltfl_flux(:,:)=(bmeltfl_flux(:,:)/dtsum)*ice_density/secyear |
---|
105 | dhdt_flux(:,:)=(dhdt_flux(:,:)/dtsum)/secyear |
---|
106 | calving_flux(:,:)=(calving_flux(:,:)/dtsum)*ice_density/secyear |
---|
107 | grline_flux(:,:)= (grline_flux(:,:)/dtsum)*ice_density/secyear |
---|
108 | calv_abl_bord_flux(:,:)=(calv_abl_bord_flux(:,:)/dtsum)*ice_density/secyear |
---|
109 | |
---|
110 | ! sorties dans debug_3d: |
---|
111 | debug_3d(:,:,106) = bm_flux(:,:) ! acabf |
---|
112 | debug_3d(:,:,107) = -bmeltgr_flux(:,:) ! libmassbfgr |
---|
113 | debug_3d(:,:,108) = dhdt_flux(:,:) ! dlithkdt |
---|
114 | debug_3d(:,:,109) = calving_flux(:,:) ! licalvf |
---|
115 | debug_3d(:,:,110) = grline_flux(:,:) ! ligroundf |
---|
116 | debug_3d(:,:,121) = -bmeltfl_flux(:,:) ! libmassbffl |
---|
117 | debug_3d(:,:,123) = calv_abl_bord_flux(:,:) + calving_flux(:,:) ! lifmassbf |
---|
118 | ! remise a 0 |
---|
119 | dtsum = 0. |
---|
120 | bm_flux(:,:) = 0. |
---|
121 | bmeltgr_flux(:,:) = 0. |
---|
122 | dhdt_flux(:,:) = 0. |
---|
123 | calving_flux(:,:) = 0. |
---|
124 | grline_flux(:,:) = 0. |
---|
125 | bmeltfl_flux(:,:) = 0. |
---|
126 | calv_abl_bord_flux(:,:) = 0. |
---|
127 | |
---|
128 | !~ else if (abs(dtsum-dtoutflux).lt.dtmin) then |
---|
129 | !~ print*,'pb sorties flux dtoutflux : time=', time |
---|
130 | !~ print*,'dtsum=',dtsum,' dtoutflux', dtoutflux |
---|
131 | !~ else if (mod(abs(time),dtoutflux).lt.dtmin) then |
---|
132 | !~ print*,'pb sorties flux dtsum : time=', time |
---|
133 | !~ print*,'dtsum=',dtsum,' dtoutflux', dtoutflux |
---|
134 | endif |
---|
135 | |
---|
136 | ! reinitialisation des variables |
---|
137 | vieuxH_flux(:,:)=H(:,:) |
---|
138 | ablbord(:,:)=0. |
---|
139 | |
---|
140 | end subroutine bilan_flux_output |
---|
141 | |
---|
142 | end module bilan_flux_mod |
---|