source: trunk/SOURCES/Ant16_files/bmelt-ant-regions-initmip_mod.f90 @ 334

Last change on this file since 334 was 225, checked in by dumas, 5 years ago

format of the parameter file created by the code, possibility to use it in reading for Ant-40
File size: 8.2 KB
Line 
1!> \file bmelt-ant-regions-initmip_mod.f90
2!! Module qui calcule la fusion basale (grounded ou ice shelves).
3!<
4
5!> \namespace BMELT_ANT_REGIONS_INITMIP
6!! Module qui calcule la fusion basale (grounded ou ice shelves)
7!! \author CatRitz, aquiquet
8!! \date juillet
9!! @note pour les ice shelves Antarctique, tient compte des différentes régions
10!! A choisir dans le module_choix
11!! @note Used module
12!! @note   - use module3D_phy
13!<
14MODULE  BMELT_ANT_REGIONS_INITMIP         ! cat juillet 2005 - afq 2017
15
16  use module3D_phy
17  use netcdf
18  use io_netcdf_grisli
19 
20  implicit none
21
22  REAL,dimension(nx,ny) ::  bmgrz       !< fusion basale a la grounding zone
23  real,dimension(nx,ny) ::  bmshelf     !< fusion basale sous shelf
24
25  REAL,dimension(nx,ny) :: typeshelf    !< Type de shelf Ronne->1 Ross ->2 ....
26
27  integer, parameter :: nb_reg = 18
28
29  real,dimension (nb_reg) :: bmelt_regions, bmgrz_regions
30  real :: bmelt_talus,bmgrz_talus       !< Au dela du talus continental
31  real :: bmelt_coef                    !< Coef pour corriger les valeurs std de bmelt (lu dans namelist)
32  real, parameter :: depth_talus=-2500. !< Profondeur du talus continental
33
34  ! For initMIP, we need to read an anomaly of bmelt...:
35  integer               :: bmelt_time        !< pour selectionner le type de calcul de bmelt 
36  real,dimension(nx,ny) :: bmelt_anom        !< bilan de masse des anomalies indices : nx,ny
37  real,dimension(nx,ny) ::  bmgrz_0               !< fusion basale a la grounding zone
38  real,dimension(nx,ny) ::  bmshelf_0             !< fusion basale sous shelf
39  integer :: flag_dist
40  real,dimension(nx,ny) :: coef_bmelt_dist
41 
42CONTAINS
43  !-------------------------------------------------------------------------------
44  subroutine init_bmelt
45
46
47    ! Cette routine fait l'initialisation pour la fusion basale.
48    ! Elle est appelée par inputfile-vec-0.5.f90
49
50    real*8, dimension(:,:),   pointer    :: tab               !< tableau 2d real pointer
51    real :: bmelt_dist0
52    character(len=100) :: file_number_shelves ! fichier avec les zones ice-shelves
53    character(len=100) :: file_bmelt_anom
54    real,dimension(nx,ny) :: dist_talus  ! distance au talus
55    character(len=100) :: file_dist_talus ! fichier distance au talus
56   
57    namelist/bmelt_ant_reg_initmip/bmelt_regions,bmgrz_regions,  & 
58         bmelt_talus,bmgrz_talus,bmelt_coef,file_number_shelves, &
59         flag_dist,file_dist_talus,bmelt_dist0 
60    namelist/bmelt_anom_initMIP/file_bmelt_anom, bmelt_time 
61
62
63    rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
64    read(num_param,bmelt_ant_reg_initmip)
65
66    write(num_rep_42,'(A)') '!  module bmelt_ant_regions_initmip '
67    write(num_rep_42,bmelt_ant_reg_initmip)                       
68    write(num_rep_42,'(A)') '! bmelt_regions    :  fonte basale sous shelves pour les 18 regions initMIP'
69    write(num_rep_42,'(A)') '! bmgrz_regions    :  fonte basale grounding zone pour les 18 regions initMIP'
70    write(num_rep_42,'(A)') '! bmelt_talus & bmgrz_talus  :  fonte basale apres talus cont'
71    write(num_rep_42,'(A)') '! bmelt_coef                 :  coef fonte (1 pour conserver val)'
72    write(num_rep_42,'(A)') '! file_numer_ice-shelves     : fichier zones ice shelves'
73    write(num_rep_42,'(A)') '! flag_dist : flag pour bmelt fnct distance talus'
74    write(num_rep_42,'(A)') '! file_dist_talus : fichier de distance talus'
75    write(num_rep_42,'(A)') '! bmelt_dist0 : coef de bmelt au talus'
76    write(num_rep_42,'(A)')'!_______________________________________________________________________'
77
78
79    ! lecture du fichier contenant les types de shelves
80    ! 18 regions sur la base des bassins fournis par le projet ISMIP6-initmip antarctica
81    typeshelf(:,:)=100
82    file_number_shelves=TRIM(DIRNAMEINP)//trim(file_number_shelves)
83    call Read_Ncdf_var('z',file_number_shelves,tab)
84    typeshelf(:,:)  = tab(:,:)
85
86
87    if (flag_dist.eq.1) then
88       file_dist_talus=TRIM(DIRNAMEINP)//trim(file_dist_talus)
89       call Read_Ncdf_var('z',file_dist_talus,tab)
90       dist_talus(:,:)  = min(max(tab(:,:),0.),400.)
91       coef_bmelt_dist(:,:)= 1. / ((1/bmelt_dist0)*(1+(bmelt_dist0-1)*(dist_talus(:,:)/400.)**2))
92    else
93       coef_bmelt_dist(:,:)= 1.
94    endif
95
96    bms_init:    do j=1,ny
97       do i=1,nx
98
99          talus:    if (Bsoc(i,j).GT.depth_talus) then
100
101             bmshelf(i,j) = bmelt_regions(typeshelf(i,j))
102             bmgrz(i,j)   = bmgrz_regions(typeshelf(i,j))
103
104          else   ! au dela du talus continental
105             
106             bmshelf(i,j)=bmelt_talus 
107             bmgrz(i,j)=bmgrz_talus   
108
109          endif talus   ! fin du test sur dist_talus
110
111       enddo
112    enddo bms_init
113   
114    bmshelf(:,:)=bmshelf(:,:)*bmelt_coef
115    bmgrz(:,:)=bmgrz(:,:)*bmelt_coef
116
117    bmshelf_0(:,:)=bmshelf(:,:)
118    bmgrz_0(:,:)=bmgrz(:,:)
119
120
121  ! ______ Anomalies...
122 
123  rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
124  read(num_param,bmelt_anom_initMIP)
125
126  write(num_rep_42,'(A)')'!  module bmelt_ant_regions_initmip                                     '
127  write(num_rep_42,bmelt_anom_initMIP)
128  write(num_rep_42,'(A)')'! file_bmelt_anom   = fichier anomalie bmelt                            '
129  write(num_rep_42,'(A)')'! bmelt_time        = 0:fixe, 1:anomalies                               '
130  write(num_rep_42,'(A)')'!_______________________________________________________________________' 
131
132  if ( bmelt_time .eq. 1 ) then 
133     file_bmelt_anom  = trim(dirnameinp)//trim(file_bmelt_anom)
134     call Read_Ncdf_var('abmb',file_bmelt_anom,tab)
135     bmelt_anom (:,:) = Tab(:,:)
136  end if
137 
138   
139    debug_3D(:,:,33)=typeshelf(:,:)
140    debug_3D(:,:,34)=bmshelf(:,:)
141
142   
143    return
144  end subroutine init_bmelt
145
146  !________________________________________________________________________________
147  subroutine bmeltshelf
148
149
150    ! cette routine calcule la fusion basale proprement dite
151
152    integer :: ngr           ! nombre de voisins flottants
153    real :: coef_talus       ! pour ne pas changer la fusion au dessus de l'ocean profond
154
155    real             :: coefanomtime
156
157    coefanomtime = min ( real(time/40.) , 1. ) 
158    if (bmelt_time.eq.0) then
159       bmshelf(:,:) = bmshelf_0(:,:)
160       bmgrz(:,:) = bmgrz_0(:,:)
161    else if (bmelt_time.eq.1) then
162       bmshelf(:,:) = bmshelf_0(:,:)+ (coefanomtime * bmelt_anom(:,:))
163       bmgrz(:,:) = bmgrz_0(:,:)    + (coefanomtime * bmelt_anom(:,:))
164    end if
165
166   
167    do j=2,ny-1
168       do i=2,nx-1
169
170          talus_nochange : if (Bsoc(i,j).GT.depth_talus) then
171             coef_talus = coefbmshelf
172          else
173             coef_talus = 1         ! pas de changement au dela du talus continental               
174          endif talus_nochange
175
176
177          shelf:    if (flot(i,j)) then    ! partie flottante
178
179             bmelt(i,j)=coef_talus*bmshelf(i,j)*coef_bmelt_dist(i,j)
180
181             if (fbm(i,j)) then
182                bmelt(i,j)=coef_talus*bmgrz(i,j)*coef_bmelt_dist(i,j)
183             endif
184
185             ! ATTENTION le bloc suivant sert a determiner la fusion basale d'equilibre
186             ! pour les shelves stationnaires
187
188             if (igrdline.eq.1) then
189                corrbmelt(i,j)=hdot(i,j)+bmelt(i,j)   ! le bmelt d'equilibre
190                debug_3D(i,j,28)=corrbmelt(i,j)
191             endif
192
193
194          else                   ! point posé, on compte le nombre de voisins flottants
195             ngr=0
196             if (flot(i+1,j)) ngr=ngr+1
197             if (flot(i-1,j)) ngr=ngr+1
198             if (flot(i,j+1)) ngr=ngr+1
199             if (flot(i,j-1)) ngr=ngr+1
200
201             !   la fusion des points limites est une combinaison entre valeur posée et valeur flottante
202             !   en fonction du nombre de points flottants
203
204             bmelt(i,j)= (ngr/4.*bmgrz(i,j)*coef_talus)*coef_bmelt_dist(i,j)+(1.-ngr/4.)*bmelt(i,j)
205
206
207          endif shelf
208
209       end do
210    end do
211
212
213    ! les bords de la grille sont mis a dist_talus=0 de force
214    ! attention, pourrait poser des problemes avec AGRIF
215
216    bmelt(1,:)=bmelt_talus   !a gauche
217    bmelt(2,:)=bmelt_talus
218
219    bmelt(nx-1,:)=bmelt_talus !  a droite
220    bmelt(nx,:)=bmelt_talus 
221
222    bmelt(:,1)=bmelt_talus    ! en bas   
223    bmelt(:,2)=bmelt_talus
224
225    bmelt(:,ny-1)=bmelt_talus    ! en haut
226    bmelt(:,ny)=bmelt_talus
227
228
229    if (igrdline.eq.1) then
230       bmelt(:,:)=0.      ! hdot donne alors le -bmelt
231    endif
232
233    return
234  end subroutine  bmeltshelf
235
236
237END MODULE  BMELT_ANT_REGIONS_INITMIP
Note: See TracBrowser for help on using the repository browser.