Changeset 461
- Timestamp:
- 01/12/24 16:26:19 (4 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/GRISLIv3/SOURCES/bmelt-beckmann-gcm_mod.f90
r449 r461 35 35 real*8 :: K_t 36 36 real*8 :: n_tour 37 real*8, dimension (:,:,:), allocatable :: temp_ocean , mask_oce2, temp_ocean_2!< thermal forcing , input37 real*8, dimension (:,:,:), allocatable :: temp_ocean !< thermal forcing , input 38 38 real*8, dimension (:,:,:), allocatable :: salinity_ocean !< thermal forcing , input 39 39 … … 41 41 42 42 ! integer, dimension (nx,ny) :: profondeur 43 real*8, dimension (:), allocatable :: mean_TF !< mean thermal forcing for a given basin (non-local formula) 44 real*8, dimension (:), allocatable :: IS_area !< ice shelf area for a given basin (non-local formula) 45 real*8, dimension(nx,ny) :: deltaT_basin 43 46 integer, dimension(nx,ny) :: bassin 44 47 real*8, dimension (:), allocatable :: zoc !< depth of oceanic levels , input … … 47 50 character(len=200) :: TO_file, SO_file, Bassin_file ! fichiers de forcage 48 51 logical, dimension(:), allocatable :: mask_bassin_vide ! mask true if bassin without any GCM data 52 integer :: nonlocal_bool !< non local ISMIP parametrisation? (1=yes) 53 integer :: grlmelt_bool !< do we apply ocean melt at the grounded grounding line points (1=yes) 49 54 50 55 contains … … 63 68 character(len=100) :: SO_file !< files read 64 69 character(len=100) :: Bassin_file !< files read 70 character(len=100) :: deltaT_file !< files read, spatial correction for Toc 71 integer :: deltaT_bool !< Do we apply a spatial correction for Toc (1=yes) 65 72 integer :: i,j,k 66 73 integer :: n,imin,imax,jmin,jmax … … 74 81 integer :: nmax ! to avoid inifite loop 75 82 integer :: err 76 namelist/bmelt_beckmann_gcm_mod/K_t,TO_file,SO_file,Bassin_file,nbassins,bmelt_empty 83 84 namelist/bmelt_beckmann_gcm_mod/K_t,TO_file,SO_file,Bassin_file,nbassins,deltaT_bool,deltaT_file,bmelt_empty,nonlocal_bool,grlmelt_bool 77 85 78 86 rewind(num_param) ! loop back at the beginning of the param_list.dat file … … 112 120 end if 113 121 end if 114 if (.not.allocated(mask_oce2)) then115 allocate(mask_oce2(nx,ny,nzoc),stat=err)116 if (err/=0) then117 print *,"erreur a l'allocation du tableau mask_oce2",err118 stop 4119 end if120 end if121 if (.not.allocated(temp_ocean_2)) then122 allocate(temp_ocean_2(nx,ny,nzoc),stat=err)123 if (err/=0) then124 print *,"erreur a l'allocation du tableau temp_ocean_2",err125 stop 4126 end if127 end if128 122 if (.not.allocated(salinity_ocean)) then 129 123 allocate(salinity_ocean(nx,ny,nzoc),stat=err) … … 151 145 if (err/=0) then 152 146 print *,"erreur a l'allocation du tableau max_depth_bassin",err 147 stop 4 148 end if 149 end if 150 if (.not.allocated(mean_TF)) then 151 allocate(mean_TF(nbassins),stat=err) 152 if (err/=0) then 153 print *,"erreur a l'allocation du tableau mean_TF",err 154 stop 4 155 end if 156 end if 157 if (.not.allocated(IS_area)) then 158 allocate(IS_area(nbassins),stat=err) 159 if (err/=0) then 160 print *,"erreur a l'allocation du tableau IS_area",err 153 161 stop 4 154 162 end if … … 163 171 call Read_Ncdf_var('so',file_inputs,tab3d) 164 172 salinity_ocean(:,:,:) = tab3d(:,:,:) 165 173 174 if (deltaT_bool.eq.1) then 175 file_inputs=TRIM(DIRNAMEINP)//'Snapshots-GCM/'//TRIM(deltaT_file) 176 call Read_Ncdf_var('deltaT_basin',file_inputs,tab2d) 177 deltaT_basin(:,:) = tab2d(:,:) 178 else 179 deltaT_basin(:,:) = 0. 180 endif 181 166 182 file_inputs=TRIM(DIRNAMEINP)//'Snapshots-GCM/'//TRIM(Bassin_file) 167 183 call Read_Ncdf_var('mask',file_inputs,tab2d) … … 226 242 where (temp_ocean(:,:,:).gt.1.e10) 227 243 mask_gcm_init(:,:,:) = 0 228 244 endwhere 229 245 230 246 ! mask_gcm = 1 si valeur dans gcm, 0 si pas de valeur => tableau update lorsqu'on rempli la grille … … 232 248 where (temp_ocean(:,:,:).gt.1.e10) 233 249 mask_gcm(:,:,:) = 0 234 250 endwhere 235 251 !$OMP END WORKSHARE 236 252 … … 352 368 !$OMP WORKSHARE 353 369 ice_draft(:,:) = S(:,:)-H(:,:)-sealevel_2d(:,:) 370 mean_TF(:) = 0.d0 371 IS_area(:) = 0.d0 354 372 !$OMP END WORKSHARE 355 373 … … 370 388 kinf = ksup - 1 371 389 if ( ice_draft(i,j) .gt. zoc(1) ) then 372 TO_draft(i,j) = temp_ocean(i,j,1)373 SO_draft(i,j) = salinity_ocean(i,j,1)374 T_freez(i,j) = 0.0939 - 0.057 * salinity_ocean(i,j,1) + 0.000764 * zoc(1)390 TO_draft(i,j) = temp_ocean(i,j,1) 391 SO_draft(i,j) = salinity_ocean(i,j,1) 392 T_freez(i,j) = 0.0939 - 0.057 * salinity_ocean(i,j,1) + 0.000764 * zoc(1) 375 393 elseif ( ice_draft(i,j) .lt. zoc(nzoc) ) then 376 TO_draft(i,j) = temp_ocean(i,j,nzoc)377 SO_draft(i,j) = salinity_ocean(i,j,nzoc)378 T_freez(i,j) = 0.0939 - 0.057 * salinity_ocean(i,j,nzoc) + 0.000764 * zoc(nzoc)394 TO_draft(i,j) = temp_ocean(i,j,nzoc) 395 SO_draft(i,j) = salinity_ocean(i,j,nzoc) 396 T_freez(i,j) = 0.0939 - 0.057 * salinity_ocean(i,j,nzoc) + 0.000764 * zoc(nzoc) 379 397 else 380 398 !TO_draft(i,j) = ( (zoc(ksup)-ice_draft(i,j)) * temp_ocean(i,j,kinf) & … … 382 400 !SO_draft(i,j) = ( (zoc(ksup)-ice_draft(i,j)) * salinity_ocean(i,j,kinf) & 383 401 ! & + (ice_draft(i,j)-zoc(kinf)) * salinity_ocean(i,j,ksup) ) / (zoc(ksup)-zoc(kinf)) 384 TO_draft(i,j) = temp_ocean(i,j,ksup) + ((ice_draft(i,j)-zoc(ksup))*(temp_ocean(i,j,kinf)-temp_ocean(i,j,ksup))/(zoc(kinf)-zoc(ksup))) 385 SO_draft(i,j) = salinity_ocean(i,j,ksup) + ((ice_draft(i,j)-zoc(ksup))*(salinity_ocean(i,j,kinf)-salinity_ocean(i,j,ksup))/(zoc(kinf)-zoc(ksup))) 386 T_freez(i,j) = 0.0939 - 0.057 * SO_draft(i,j) + 0.000764 * ice_draft(i,j) 402 TO_draft(i,j) = temp_ocean(i,j,ksup) + ((ice_draft(i,j)-zoc(ksup))*(temp_ocean(i,j,kinf)-temp_ocean(i,j,ksup))/(zoc(kinf)-zoc(ksup))) 403 SO_draft(i,j) = salinity_ocean(i,j,ksup) + ((ice_draft(i,j)-zoc(ksup))*(salinity_ocean(i,j,kinf)-salinity_ocean(i,j,ksup))/(zoc(kinf)-zoc(ksup))) 404 T_freez(i,j) = 0.0939 - 0.057 * SO_draft(i,j) + 0.000764 * ice_draft(i,j) 405 endif 406 407 if (nonlocal_bool.eq.1) then 408 mean_TF( bassin(i,j) ) = mean_TF( bassin(i,j) ) + mesh_area(i,j) * (TO_draft(i,j)-T_freez(i,j)) 409 IS_area( bassin(i,j) ) = IS_area( bassin(i,j) ) + mesh_area(i,j) 387 410 endif 388 411 … … 401 424 enddo 402 425 !$OMP END DO 426 427 where ( IS_area(:).gt.0.d0) ! for all basins that have ice shelves 428 mean_TF(:) = mean_TF(:) / IS_area(:) 429 elsewhere ! we have no floating points over the whole basin 430 mean_TF(:) = -9999.d0 431 endwhere 403 432 404 433 ! 3 - Calculation of melting rate : … … 409 438 do j=1,ny 410 439 do i=1,nx 411 if ( TO_draft(i,j) .gt. -9000.d0 ) then ! floating points 440 !if ( TO_draft(i,j) .gt. -9000.d0 ) then ! floating points 441 if ( flot(i,j) ) then ! floating points 412 442 !if (mask_bassin_vide(bassin(i,j)) .or. (bassin(i,j).eq.1)) then ! bassins without any GCM data or bassin 1 413 443 !if (mask_bassin_vide(bassin(i,j)) .or. (bassin(i,j).eq.1)) then !bassins without any GCM data … … 415 445 bmelt(i,j) = bmelt_empty 416 446 else ! Bassins with GCM data 417 bmelt(i,j) = coef_OM * ( TO_draft(i,j) - T_freez(i,j) )* abs( TO_draft(i,j) - T_freez(i,j) ) 447 if (nonlocal_bool.eq.1) then 448 bmelt(i,j) = coef_OM * ( TO_draft(i,j) + deltaT_basin(i,j) - T_freez(i,j) )* abs( mean_TF(bassin(i,j)) + deltaT_basin(i,j) ) 449 else 450 bmelt(i,j) = coef_OM * ( TO_draft(i,j) + deltaT_basin(i,j) - T_freez(i,j) )* abs( TO_draft(i,j) + deltaT_basin(i,j) - T_freez(i,j) ) 451 endif 418 452 endif 419 453 else … … 427 461 do j=2,ny-1 428 462 do i=2,nx-1 429 if ( TO_draft(i,j) .le. -9000.d0 ) then !grounded points 430 ngr=0 431 bmloc=0.d0 432 !if (TO_draft(i+1,j) .le. -9000.d0 ) then 433 if (flot(i+1,j)) then 434 ngr=ngr+1 435 bmloc=bmloc+bmelt(i+1,j) 463 !if ( TO_draft(i,j) .le. -9000.d0 ) then !grounded points 464 if ( .not. flot(i,j) ) then 465 if (grlmelt_bool.eq.1) then 466 ngr=0 467 bmloc=0.d0 468 !if (TO_draft(i+1,j) .le. -9000.d0 ) then 469 if (flot(i+1,j)) then 470 ngr=ngr+1 471 bmloc=bmloc+bmelt(i+1,j) 472 endif 473 !if (TO_draft(i-1,j) .le. -9000.d0 ) then 474 if (flot(i-1,j)) then 475 ngr=ngr+1 476 bmloc=bmloc+bmelt(i-1,j) 477 endif 478 !if (TO_draft(i,j+1) .le. -9000.d0 ) then 479 if (flot(i,j+1)) then 480 ngr=ngr+1 481 bmloc=bmloc+bmelt(i,j+1) 482 endif 483 !if (TO_draft(i,j-1) .le. -9000.d0 ) then 484 if (flot(i,j-1)) then 485 ngr=ngr+1 486 bmloc=bmloc+bmelt(i,j-1) 487 endif 488 else 489 ngr=0 ! afq -- I do not allow subgrid ocean melt 436 490 endif 437 !if (TO_draft(i-1,j) .le. -9000.d0 ) then438 if (flot(i-1,j)) then439 ngr=ngr+1440 bmloc=bmloc+bmelt(i-1,j)441 endif442 !if (TO_draft(i,j+1) .le. -9000.d0 ) then443 if (flot(i,j+1)) then444 ngr=ngr+1445 bmloc=bmloc+bmelt(i,j+1)446 endif447 !if (TO_draft(i,j-1) .le. -9000.d0 ) then448 if (flot(i,j-1)) then449 ngr=ngr+1450 bmloc=bmloc+bmelt(i,j-1)451 endif452 ngr=0 ! afq -- I do not allow subgrid ocean melt453 491 bmelt(i,j)=(ngr/4.)*bmloc+(1-ngr/4.)*bmelt(i,j) 454 492 endif
Note: See TracChangeset
for help on using the changeset viewer.