- Timestamp:
- 01/09/19 17:09:26 (5 years ago)
- Location:
- trunk/SOURCES
- Files:
-
- 51 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/SOURCES/3D-physique-gen_mod.f90
r166 r237 141 141 real :: SECYEAR !< for relation an/seconds 142 142 !real :: S22 ! ct for PDD calculation 143 real :: sealevel !< niveau des mers 143 real :: sealevel !< niveau des mers, afq: now used as eustatic sea level 144 144 ! real :: SIGMA ! variabilite Tday 145 145 real :: SURF !< … … 235 235 integer,dimension(nx,ny) :: Mk_init !< initial mask (with islands, outcrops, ... 236 236 integer,dimension(nx,ny) :: MK !< masks (ice sheet, max, above water, below water, 1) 237 integer,dimension(nx,ny) :: M_sealev !< masks for sea level calculation (present ice sheet)238 237 integer,dimension(nx,ny) :: MNEG !< masks (ice sheet, max, above water, below water, 1) 239 238 integer,dimension(nx,ny) :: IBASE !< type de base (froide, temperee) … … 245 244 real,dimension(nx,ny) :: ACQUA !< Surface des surfaces en eau 246 245 real,dimension(nx,ny) :: B !< Altitude de la base de la glace 'o' 247 real,dimension(nx,ny) :: B_sealev !< Altitude de la base de la glace pour le248 !< calcul du nivx de la mer (/present)249 246 real,dimension(nx,ny) :: BDOT !< derivee de B / t 250 247 real,dimension(nx,ny) :: B1 !< … … 310 307 311 308 312 real,dimension(nx,ny) :: H_sealev !< ice thickness for sea level calculation ! (present ice sheet)313 309 real,dimension(nx,ny) :: HDOT !< ice thickness derivee / t 314 310 real,dimension(nx,ny) :: HDOTWATER … … 334 330 real,dimension(nx,ny) :: SW !< for bedrock isostasy 335 331 real,dimension(nx,ny) :: S !< altitude of ice sheet surface 336 real,dimension(nx,ny) :: S_sealev !< altitude of ice sheet surface for sea 337 !< level calculation (present ice sheet) 332 real,dimension(nx,ny) :: sealevel_2d !< local sea surface elevation 338 333 real,dimension(nx,ny) :: SLOPE !< slope 'o' 339 334 real,dimension(nx,ny) :: SDX !< slope derivee / x '>' … … 344 339 real,dimension(nx,ny) :: SLOPE2my !< = Sdy**2 + Sdxmy**2 '^' 345 340 real,dimension(nx,ny) :: S0 !< altitude actuelle de la surface 346 real,dimension(nx,ny) :: slv !< niveau de flottaison (sealevel et lakes)341 !afq -- replaced by sealevel_2D: real,dimension(nx,ny) :: slv !< niveau de flottaison (sealevel et lakes) 347 342 real,dimension(nx,ny) :: TJULY !< Ground air temperature July 348 343 real,dimension(nx,ny) :: TANN !< Ground air temperature annual -
trunk/SOURCES/Ant15_CISM_files/output_anta_mod-0.4.f90
r4 r237 83 83 84 84 ! calcul de la hauteur au dessus de la flottaison 85 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers85 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 86 86 volf=volf+h(i,j) 87 87 else 88 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))88 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 89 89 endif 90 90 -
trunk/SOURCES/Ant16_files/climat_InitMIP_years_perturb_mod.f90
r180 r237 16 16 17 17 18 use module3d_phy,only: nx,ny,S,S0,H,Bsoc,acc,abl,BM,Tann,Tjuly,Ts,time,dt,num_param,num_rep_42,num_forc,dirforcage,dirnameinp,tafor,time,sealevel, coefbmshelf18 use module3d_phy,only: nx,ny,S,S0,H,Bsoc,acc,abl,BM,Tann,Tjuly,Ts,time,dt,num_param,num_rep_42,num_forc,dirforcage,dirnameinp,tafor,time,sealevel,sealevel_2d,coefbmshelf 19 19 use netcdf 20 20 use io_netcdf_grisli … … 285 285 100 continue 286 286 287 sealevel_2d(:,:) = sealevel 288 287 289 Tann (:,:) = Ta0 (:,:) + T_lapse_rate * (S(:,:)-S0(:,:)) +Tafor 288 290 Ts(:,:) = Tann(:,:) -
trunk/SOURCES/Ant40_files/lect-anteis_mod.f90
r227 r237 228 228 do J=1,NY 229 229 do I=1,NX 230 if ((BSOC(I,J)+H(I,J)*RO/ROW - SEALEVEL).LT.0.) then230 if ((BSOC(I,J)+H(I,J)*RO/ROW -sealevel_2d(i,j)).LT.0.) then 231 231 FLOT(I,J)=.TRUE. 232 232 else … … 260 260 ymax=ycc(nx,ny)/1000. 261 261 262 !!$ ! lecture du fichier de reference pour le calcul du niveau des mers : etat actuel263 !!$ open(88,file=TRIM(DIRNAMEINP)//'grzone56-k000.pl')264 !!$265 !!$ read(88,*)266 !!$ read(88,*)267 !!$ read(88,*)268 !!$ do j=1,ny269 !!$ do i=1,nx270 !!$ read(88,*) S_sealev(i,j),H_sealev(i,j),B_sealev(i,j),M_sealev(i,j)271 !!$ enddo272 !!$ enddo273 !!$ close(88)274 !!$275 !!$ write(num_rep_42,*) 'fichier reference pour le niveau des mers : ', &276 !!$ TRIM(DIRNAMEINP)//'grzone56-k000.pl'277 278 262 ! lecture du flux geothermique de Shapiro 279 263 open(88,file=ghf_fich) -
trunk/SOURCES/Ant40_files/output_anta40_mod-0.4.f90
r148 r237 191 191 192 192 ! calcul de la hauteur au dessus de la flottaison 193 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers193 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 194 194 volf=volf+h(i,j)*corrsurf(i,j) ! volume au-dessus de la flottaison 195 195 else 196 volf=volf+(h(i,j)-row/ro*(sealevel -b(i,j)))*corrsurf(i,j) ! volume au-dessus de la flottaison196 volf=volf+(h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)))*corrsurf(i,j) ! volume au-dessus de la flottaison 197 197 endif 198 198 -
trunk/SOURCES/Ant45_CISM_files/output_anta_mod-0.4.f90
r4 r237 84 84 85 85 ! calcul de la hauteur au dessus de la flottaison 86 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers86 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 87 87 volf=volf+h(i,j) 88 88 else 89 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))89 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 90 90 endif 91 91 -
trunk/SOURCES/Antarctique_general_files/output_anta_mod-0.4.f90
r68 r237 192 192 193 193 ! calcul de la hauteur au dessus de la flottaison 194 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers194 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 195 195 volf=volf+h(i,j) 196 196 else 197 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))197 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 198 198 endif 199 199 … … 289 289 scal_np = scal_np + 1 290 290 scal_vol = scal_vol + H(i,j) 291 scal_volbuoy = scal_volbuoy + H(i,j) - row/ro * max(0.,sealevel -b(i,j))291 scal_volbuoy = scal_volbuoy + H(i,j) - row/ro * max(0.,sealevel_2d(i,j)-b(i,j)) 292 292 scal_meanhdot = scal_meanhdot + Hdot(i,j) 293 293 scal_sigma_hdot = scal_sigma_hdot + Hdot(i,j)**2 -
trunk/SOURCES/Draggings_modules/dragging_prescr_beta_buoyency_mod.f90
r70 r237 500 500 501 501 502 where ( sealevel -B1(:,:).le.0.) ! socle au dessus du niveau des mers502 where ( sealevel_2d(:,:)-B1(:,:).le.0.) ! socle au dessus du niveau des mers 503 503 H_buoy(:,:) = H1(:,:) 504 504 505 505 elsewhere 506 H_buoy(:,:) = H1(:,:)-row/ro*(sealevel -B1(:,:))506 H_buoy(:,:) = H1(:,:)-row/ro*(sealevel_2d(:,:)-B1(:,:)) 507 507 508 508 end where -
trunk/SOURCES/Draggings_modules/dragging_prescr_beta_nolin_mod.f90
r168 r237 606 606 607 607 608 where ( sealevel -B1(:,:).le.0.) ! socle au dessus du niveau des mers608 where ( sealevel_2d(:,:)-B1(:,:).le.0.) ! socle au dessus du niveau des mers 609 609 H_buoy(:,:) = H1(:,:) 610 610 611 611 elsewhere 612 H_buoy(:,:) = H1(:,:)-row/ro*(sealevel -B1(:,:))612 H_buoy(:,:) = H1(:,:)-row/ro*(sealevel_2d(:,:)-B1(:,:)) 613 613 614 614 end where -
trunk/SOURCES/Eurasie40_files/bmelt-eurasie-depth-lake_mod.f90
r4 r237 85 85 do i=1,nx 86 86 87 if ((sealevel -BSOC(i,j)).lt.450) then87 if ((sealevel_2d(i,j)-BSOC(i,j)).lt.450) then 88 88 bmshelf(i,j)=2.00 !3.00! 2.55! 0.45! 0.65 89 89 bmgrz(i,j)=2.00 !3.00! 2.55! 0.45 ! 0.65 … … 127 127 128 128 129 if (s lv(i,j).gt.0) then129 if (sealevel_2d(i,j).gt.0) then 130 130 bmshelf(i,j)=0.05!5.0 131 131 bmgrz(i,j)=0.05!5.0 … … 135 135 136 136 137 if (s lv(i,j).gt.sealevel) then137 if (sealevel_2d(i,j).gt.sealevel) then 138 138 bmelt(i,j)=bmshelf(i,j) 139 139 if (fbm(i,j)) bmelt(i,j)=bmgrz(i,j) … … 166 166 ! en fonction du nombre de points flottants 167 167 168 if (s lv(i,j).gt.sealevel) then168 if (sealevel_2d(i,j).gt.sealevel) then 169 169 bmelt(i,j)= ngr/4.*bmgrz(i,j)+(1.-ngr/4.)*bmelt(i,j) 170 170 else -
trunk/SOURCES/Eurasie40_files/climat-forcage-eurasie_mod-0.4.f90
r4 r237 493 493 nom_table='lvflo' 494 494 call printtable_r(levelflot,nom_table) 495 495 496 sealevel_2d(:,:) = sealevel 497 496 498 ! coefbmshelf(:,:)=(1.+delTatime(:,:)/7.) 497 499 !coefbmshelf=1. ! modif pour test de sensibilite (tof 20 avril 01) -
trunk/SOURCES/Eurasie40_files/lakes-prescribed_mod-0.1.f90
r4 r237 54 54 subroutine lake_to_slv 55 55 56 ! inclut les lacs dans le tableau s lvdes niveaux de flottaison56 ! inclut les lacs dans le tableau sealevel_2d des niveaux de flottaison 57 57 58 s lv(:,:)=sealevel58 sealevel_2d(:,:)=sealevel 59 59 60 60 ! VOIR COMMENT ACTIVER-DESACTIVER LES LACS EN FONCTION DU TEMPS … … 64 64 do i_lake=1,nb_lakes 65 65 if (lake(i,j).eq.i_lake) then 66 s lv(i,j)=max(level_lakes(i_lake),Bsoc(i,j))66 sealevel_2d(i,j)=max(level_lakes(i_lake),Bsoc(i,j)) 67 67 end if 68 68 end do -
trunk/SOURCES/Eurasie40_files/output_eurasie40_mod-0.1.f90
r4 r237 85 85 86 86 ! calcul de la hauteur au dessus de la flottaison 87 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers87 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 88 88 volf=volf+h(i,j) 89 89 else 90 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))90 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 91 91 endif 92 92 -
trunk/SOURCES/GrIce2sea_files/climat_GrIce2sea_years_perturb_mod.f90
r225 r237 16 16 17 17 18 use module3d_phy,only: nx,ny,S,S0,H,Bsoc,acc,abl,BM,Tann,Tjuly,Ts,time,dt,num_param,num_rep_42,num_forc,dirforcage,dirnameinp,tafor,time,sealevel, coefbmshelf18 use module3d_phy,only: nx,ny,S,S0,H,Bsoc,acc,abl,BM,Tann,Tjuly,Ts,time,dt,num_param,num_rep_42,num_forc,dirforcage,dirnameinp,tafor,time,sealevel,sealevel_2d,coefbmshelf 19 19 !use lect_climref_Ice2sea 20 20 use netcdf … … 382 382 endif 383 383 100 continue 384 385 sealevel_2d(:,:) = sealevel 384 386 385 387 Tann (:,:) = Ta0 (:,:) + T_lapse_rate * (S(:,:)-S0(:,:)) +Tafor -
trunk/SOURCES/GrIce2sea_files/lect_GrIce2sea_gen_nc.f90
r4 r237 207 207 208 208 209 ! lecture du fichier de reference pour le calcul du niveau des mers : etat actuel210 ! pas de version correcte pour l'instant -> valeurs initiales211 S_sealev(:,:) = S0(:,:)212 H_sealev(:,:) = H0(:,:)213 B_sealev(:,:) = Bsoc(:,:)214 M_sealev(:,:) = Mk0(:,:)215 216 209 if (itracebug.eq.1) call tracebug(' apres lect entete') 217 210 -
trunk/SOURCES/GrIce2sea_files/output_Grice2sea_mod.f90
r11 r237 196 196 197 197 ! calcul de la hauteur au dessus de la flottaison 198 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers198 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 199 199 volf=volf+h(i,j) 200 200 else 201 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))201 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 202 202 endif 203 203 … … 298 298 scal_np = scal_np + 1 299 299 scal_vol = scal_vol + H(i,j) 300 scal_volbuoy = scal_volbuoy + H(i,j) - row/ro * max(0.,sealevel -b(i,j))300 scal_volbuoy = scal_volbuoy + H(i,j) - row/ro * max(0.,sealevel_2d(i,j)-b(i,j)) 301 301 scal_meanhdot = scal_meanhdot + Hdot(i,j) 302 302 scal_sigma_hdot = scal_sigma_hdot + Hdot(i,j)**2 -
trunk/SOURCES/Greeneem_files/lect-clim-act-greeneem_mar_mod.f90
r4 r237 55 55 do i=1,nx 56 56 read (20,*) S0CLIM(i,j) 57 S0CLIM(i,j)=max(S0CLIM(i,j),s lv(i,j)) ! pour etre au niveau des mers57 S0CLIM(i,j)=max(S0CLIM(i,j),sealevel_2d(i,j)) ! pour etre au niveau des mers 58 58 end do 59 59 end do … … 146 146 do j=1,ny 147 147 do i=1,nx 148 alti=max(s lv(i,j),S(I,J))148 alti=max(sealevel_2d(i,j),S(I,J)) 149 149 Tann(i,j)=TA0(i,j)+lapse_ann*(alti-S0CLIM(i,j)) 150 150 Tjuly(i,j)=TJ0(i,j)+lapse_july*(alti-S0CLIM(i,j)) -
trunk/SOURCES/Greeneem_files/lect-clim-act-greeneem_mois_lapsecouche_mod.f90
r4 r237 114 114 Tm_0(i,j,mo)=Tm(i,j,mo) ! temp sur la calotte du GCM a l'instant initial 115 115 116 zel=max(s lv(i,j),S0CLIM(i,j))116 zel=max(sealevel_2d(i,j),S0CLIM(i,j)) 117 117 ! on calcule une temperature au sol a partir de la temperature sur la calotte de ref du GCM 118 118 … … 131 131 end if 132 132 133 !zel=max(s lv(i,j),S(i,j))133 !zel=max(sealevel_2d(i,j),S(i,j)) 134 134 zel=max(0.,S(i,j)) 135 135 ! puis on calcule la temperature sur la calotte S -
trunk/SOURCES/Greeneem_files/lect-clim-act-greeneem_mois_mod.f90
r4 r237 74 74 do i=1,nx 75 75 read (20,*) S0CLIM(i,j) 76 S0CLIM(i,j)=max(S0CLIM(i,j),s lv(i,j)) ! pour etre au niveau des mers76 S0CLIM(i,j)=max(S0CLIM(i,j),sealevel_2d(i,j)) ! pour etre au niveau des mers 77 77 end do 78 78 end do … … 120 120 do i=1,nx 121 121 ! Zs est l'altitude de la surface qu'elle soit mer, glace ou lac 122 ZS(i,j)=max(s lv(i,j),S0CLIM(I,J)) !slv=sealev pour nolake122 ZS(i,j)=max(sealevel_2d(i,j),S0CLIM(I,J)) !sealevel_2d=sealev pour nolake 123 123 124 124 do mo=1,mois … … 261 261 do j=1,ny 262 262 do i=1,nx 263 ZS(i,j)=max(s lv(i,j),S(i,j))263 ZS(i,j)=max(sealevel_2d(i,j),S(i,j)) 264 264 Tann(i,j)=d_ann+gamma_ann*S(i,j)+c_ann*ylat(i,j)+kappa_ann*(360-xlong(i,j)) 265 265 Tjuly(i,j)=d_july+gamma_july*S(i,j)+c_july*ylat(i,j)+kappa_july*(360-xlong(i,j)) -
trunk/SOURCES/Greeneem_files/lect-greeneem_mod.f90
r113 r237 86 86 do i=1,nx 87 87 88 archimed = Bsoc(i,j)+H(i,j)*ro/row -sealevel 88 archimed = Bsoc(i,j)+H(i,j)*ro/row -sealevel_2d(i,j) 89 89 90 90 if ((archimed.LT.0.).and.(H(I,J).gt.1.E-3)) then ! le point flotte -
trunk/SOURCES/Greeneem_files/output_greeneem_mod-0.4.f90
r102 r237 141 141 do j=1,ny 142 142 do i=1,nx 143 if (mask_cal(i,j,k).and..not.flot(i,j).and.H(i,j).GT.1..and.((sealevel -B(i,j)).LE.0.)) then143 if (mask_cal(i,j,k).and..not.flot(i,j).and.H(i,j).GT.1..and.((sealevel_2d(i,j)-B(i,j)).LE.0.)) then 144 144 isvolf(kk) = isvolf(kk) + H(i,j) 145 145 else 146 isvolf(kk) = isvolf(kk) + H(i,j)-row/ro*(sealevel -B(i,j))146 isvolf(kk) = isvolf(kk) + H(i,j)-row/ro*(sealevel_2d(i,j)-B(i,j)) 147 147 endif 148 148 enddo -
trunk/SOURCES/Greenmint40_files/lect-greenmint_mod.f90
r4 r237 112 112 do i=1,nx 113 113 114 archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel 114 archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel_2d(i,j) 115 115 116 116 if ((ARCHIM.LT.0.).and.(H(I,J).gt.1.E-3)) then ! le point flotte -
trunk/SOURCES/Greenmint40_files/output_green_mod-0.4.f90
r4 r237 89 89 ! calcul de la hauteur au dessus de la flottaison 90 90 91 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers91 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 92 92 volf=volf+h(i,j) 93 93 else 94 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))94 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 95 95 endif 96 96 -
trunk/SOURCES/Heino_files/flottab2-0.5-heino.f90
r4 r237 135 135 uys1(i,j)=uybar(i,j) 136 136 137 archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel 137 archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel_2d(i,j) 138 138 139 139 arch: if ((ARCHIM.LT.0.).and.(H(I,J).gt.1.E-3)) then ! le point flotte … … 146 146 if (igrdline.eq.1) then ! en cas de grounding line prescrite 147 147 flot(i,j)=.false. 148 H(i,j)=(10.+sealevel -Bsoc(i,j))*row/ro ! pour avoir archim=10148 H(i,j)=(10.+sealevel_2d(i,j)-Bsoc(i,j))*row/ro ! pour avoir archim=10 149 149 new_flot_point(i,j)=.false. 150 150 endif … … 162 162 if (igrdline.eq.1) then ! en cas de grounding line prescrite 163 163 flot(i,j)=.false. 164 H(i,j)=(10.+sealevel -Bsoc(i,j))*row/ro ! pour avoir archim=10164 H(i,j)=(10.+sealevel_2d(i,j)-Bsoc(i,j))*row/ro ! pour avoir archim=10 165 165 new_flot_point(i,j)=.false. 166 166 endif … … 184 184 185 185 surnet=H(i,j)*(1.-ro/row) 186 S(i,j)=surnet+sealevel 186 S(i,j)=surnet+sealevel_2d(i,j) 187 187 B(i,j)=S(i,j)-H(i,j) 188 188 end if … … 234 234 ! sur le demi noeud condition archim > 100 m 235 235 236 archim=(Bsoc(i,j)+Bsoc(i-1,j))*0.5-sealevel +ro/row*Hmx(i,j)236 archim=(Bsoc(i,j)+Bsoc(i-1,j))*0.5-sealevel_2d(i,j)+ro/row*Hmx(i,j) 237 237 gzmx(i,j)=gzmx(i,j).and.(archim.le.100.) 238 238 … … 301 301 ! sur le demi noeud condition archim > 100 m 302 302 303 archim=(Bsoc(i,j)+Bsoc(i,j-1))*0.5-sealevel +ro/row*Hmy(i,j)303 archim=(Bsoc(i,j)+Bsoc(i,j-1))*0.5-sealevel_2d(i,j)+ro/row*Hmy(i,j) 304 304 gzmy(i,j)=gzmy(i,j).and.(archim.le.100.) 305 305 … … 554 554 !!!!!!!!! on teste si la tache n'est pas completement ice stream 555 555 if (icetrim(table_out(i,j))) THEN ! on a une ile d'ice stream !!!! TEST2 556 smax_i=0 ; smax_j=0 ; smax_=sealevel 556 smax_i=0 ; smax_j=0 ; smax_=sealevel !afq -- warning: what about local sealevel? 557 557 DO II=3,NX-2 558 558 DO JJ=3,NY-2 -
trunk/SOURCES/Hemin40_files/bmelt-hemin40-depth_mod.f90
r4 r237 85 85 do i=1,nx 86 86 87 if ((sealevel -bsoc(i,j)).lt.600) then87 if ((sealevel_2d(i,j)-bsoc(i,j)).lt.600) then 88 88 bmshelf(i,j)=0.2 ! 0.45! 0.65 89 89 bmgrz(i,j)=0.2 ! 0.45 ! 0.65 … … 129 129 shelf: if (flot(i,j)) then ! partie flottante 130 130 131 if ((sealevel -bsoc(i,j)).lt.600) then131 if ((sealevel_2d(i,j)-bsoc(i,j)).lt.600) then 132 132 if ((j.lt.75).and.(i.lt.145) ) then !Atlantique 133 133 bmelt(i,j)=(.6+.4*coefbmshelf)*bmshelf(i,j) -
trunk/SOURCES/Hudson_files/bmelt_hudson_mod.f90
r4 r237 52 52 ! bmgrz(i,j)= 10. 53 53 54 bmshelf(i,j)=(sealevel -BSOC(i,j))/200 ! jalv: eq lineaire en fonction de la profondeur de la mer54 bmshelf(i,j)=(sealevel_2d(i,j)-BSOC(i,j))/200 ! jalv: eq lineaire en fonction de la profondeur de la mer 55 55 bmgrz(i,j)=bmshelf(i,j) 56 56 ! endif -
trunk/SOURCES/Hudson_files/climat-hudson_mod.f90
r4 r237 102 102 !jalv: test pour que le bilan de masse soit 0 sur l'ocean plus profond de 600 m: 103 103 ! c'est le bon critere??: 104 where (sealevel -BSOC(:,:).gt.600) acc(:,:)=0.105 where (sealevel -BSOC(:,:).gt.600) abl(:,:)=0.104 where (sealevel_2d(:,:)-BSOC(:,:).gt.600) acc(:,:)=0. 105 where (sealevel_2d(:,:)-BSOC(:,:).gt.600) abl(:,:)=0. 106 106 107 where ((sealevel -BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2500.)) abl(:,:)=-1.108 where ((sealevel -BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2600.)) abl(:,:)=-5.109 where ((sealevel -BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2650.)) abl(:,:)=-10.110 where ((sealevel -BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2700.)) abl(:,:)=-20.107 where ((sealevel_2d(:,:)-BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2500.)) abl(:,:)=-1. 108 where ((sealevel_2d(:,:)-BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2600.)) abl(:,:)=-5. 109 where ((sealevel_2d(:,:)-BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2650.)) abl(:,:)=-10. 110 where ((sealevel_2d(:,:)-BSOC(:,:).gt.100).and.(dist_heino(:,:).gt.2700.)) abl(:,:)=-20. 111 111 112 112 bm(:,:)=acc(:,:)+abl(:,:) … … 135 135 ! c'est le bon critere??: 136 136 !where (BSOC(:,:).lt.B(:,:)) Tann(:,:)=0. 137 where (sealevel -BSOC(:,:).gt.600) Tann(:,:)=0.137 where (sealevel_2d(:,:)-BSOC(:,:).gt.600) Tann(:,:)=0. 138 138 139 139 -
trunk/SOURCES/Hudson_files/eaubasale-0.5_hudson_mod.f90
r4 r237 139 139 if (flot(i,j))then ! points flottants 140 140 klimit(i,j)=1 141 limit_hw(i,j)=(sealevel -Bsoc(i,j))*rowg/rofreshg141 limit_hw(i,j)=(sealevel_2d(i,j)-Bsoc(i,j))*rowg/rofreshg 142 142 143 143 else if (IBASE(I,J).eq.1) then ! base froide … … 182 182 pot_w(I,J)=pot_w(I,J)+rog*H(I,J) 183 183 184 pot_f(I,J)=rowg*(sealevel -S(i,j)+H(I,J)) ! pression a la base de l'ice shelf184 pot_f(I,J)=rowg*(sealevel_2d(i,j)-S(i,j)+H(I,J)) ! pression a la base de l'ice shelf 185 185 enddo 186 186 enddo … … 308 308 309 309 if (flot(i,j)) then ! if flot > hwater=hwater in ocean 310 hwater(i,j)= sealevel - bsoc(i,j)310 hwater(i,j)= sealevel_2d(i,j) - bsoc(i,j) 311 311 ! hwater(i,j)= max(0.,hwater(i,j)) 312 312 if (hwater(i,j).lt.0.) hwater(i,j)=0. … … 320 320 321 321 ! Attention le bloc suivant est pour le run 20 322 ! Zflot=row/ro*(sealevel -Bsoc(i,j))-10.322 ! Zflot=row/ro*(sealevel_2d(i,j)-Bsoc(i,j))-10. 323 323 Zflot=H(i,j)*rog/rofreshg 324 324 -
trunk/SOURCES/MISMIP3D_files/lect-mismip3d_mod.f90
r4 r237 111 111 do j=1,ny 112 112 do i=1,nx 113 if ((bsoc(i,j)+h(i,j)*ro/row -sealevel ).lt.0.) then113 if ((bsoc(i,j)+h(i,j)*ro/row -sealevel_2d(i,j)).lt.0.) then 114 114 flot(i,j)=.true. 115 115 else -
trunk/SOURCES/New-remplimat/remplimat-shelves-tabTu.f90
r180 r237 483 483 484 484 case(-2) ! bord Est 485 call vel_U_East(pvi(i,j),pvi(i-1,j),rog,H(i-1,j),rowg,s lv(i-1,j),B(i-1,j),dx,Tu(i,j,:,:),Tv(i,j,:,:),opposx(i,j))485 call vel_U_East(pvi(i,j),pvi(i-1,j),rog,H(i-1,j),rowg,sealevel_2d(i-1,j),B(i-1,j),dx,Tu(i,j,:,:),Tv(i,j,:,:),opposx(i,j)) 486 486 487 487 case(-3) ! bord nord … … 489 489 490 490 case(-4) ! bord West 491 call vel_U_West(pvi(i,j),pvi(i-1,j),rog,H(i,j),rowg,s lv(i,j),B(i,j),dx,Tu(i,j,:,:),Tv(i,j,:,:),opposx(i,j))491 call vel_U_West(pvi(i,j),pvi(i-1,j),rog,H(i,j),rowg,sealevel_2d(i,j),B(i,j),dx,Tu(i,j,:,:),Tv(i,j,:,:),opposx(i,j)) 492 492 493 493 case(-12) ! coin SE … … 527 527 528 528 case(-1) ! bord sud 529 call vel_V_South(pvi(i,j),pvi(i,j-1),rog,H(i,j),rowg,s lv(i,j),B(i,j),dx,Sv(i,j,:,:),Su(i,j,:,:),opposy(i,j))529 call vel_V_South(pvi(i,j),pvi(i,j-1),rog,H(i,j),rowg,sealevel_2d(i,j),B(i,j),dx,Sv(i,j,:,:),Su(i,j,:,:),opposy(i,j)) 530 530 531 531 case(-2) ! bord Est … … 533 533 534 534 case(-3) ! bord nord 535 call vel_V_North(pvi(i,j),pvi(i,j-1),rog,H(i-1,j),rowg,s lv(i-1,j),B(i-1,j),dx,Sv(i,j,:,:),Su(i,j,:,:),opposy(i,j))535 call vel_V_North(pvi(i,j),pvi(i,j-1),rog,H(i-1,j),rowg,sealevel_2d(i-1,j),B(i-1,j),dx,Sv(i,j,:,:),Su(i,j,:,:),opposy(i,j)) 536 536 537 537 case(-4) ! bord West -
trunk/SOURCES/Snowball_files/bmelt-snowball-depth_mod.f90
r57 r237 66 66 do i=1,nx 67 67 68 if ((sealevel -bsoc(i,j)).lt.600) then68 if ((sealevel_2d(i,j)-bsoc(i,j)).lt.600) then 69 69 bmshelf(i,j)=0.2 ! 0.45! 0.65 70 70 bmgrz(i,j)=0.2 ! 0.45 ! 0.65 … … 106 106 shelf: if (flot(i,j)) then ! partie flottante 107 107 108 if ((sealevel -bsoc(i,j)).lt.600) then108 if ((sealevel_2d(i,j)-bsoc(i,j)).lt.600) then 109 109 if ((j.lt.75).and.(i.lt.145) ) then !Atlantique 110 110 bmelt(i,j)=(.6+.4*coefbmshelf)*bmshelf(i,j) -
trunk/SOURCES/Snowball_files/output_snowball_mod-0.4.f90
r57 r237 83 83 vol=vol+h(i,j) 84 84 ! calcul de la hauteur au dessus de la flottaison 85 if (sealevel -B(i,j).le.0.) then ! socle au dessus du niveau des mers85 if (sealevel_2d(i,j)-B(i,j).le.0.) then ! socle au dessus du niveau des mers 86 86 volf=volf+h(i,j) 87 87 else 88 volf=volf+h(i,j)-row/ro*(sealevel -b(i,j))88 volf=volf+h(i,j)-row/ro*(sealevel_2d(i,j)-b(i,j)) 89 89 endif 90 90 -
trunk/SOURCES/bilan_eau_mod.f90
r192 r237 45 45 subroutine init_bilan_eau 46 46 ! initialisation des variables 47 48 49 50 51 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. 52 52 alpha_flot=ro/row 53 archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot - sealevel 53 archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot - sealevel_2d(:,:) 54 54 gr_line(:,:)=0 55 55 do j=1,ny 56 56 do i=1,nx 57 !afq if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel )) then ! grounded with ice57 !afq if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel_2d(i,j))) then ! grounded with ice 58 58 if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.)) then ! grounded with ice 59 59 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 … … 67 67 end subroutine init_bilan_eau 68 68 69 69 70 70 71 71 … … 101 101 102 102 where (ice(2:nx-1,2:ny-1).eq.1) 103 104 103 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 104 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 105 105 endwhere 106 106 107 archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot - sealevel 107 archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot - sealevel_2d(:,:) 108 108 gr_line(:,:)=0 109 109 do j=1,ny 110 110 do i=1,nx 111 !afq if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel )) then ! grounded with ice111 !afq if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel_2d(i,j))) then ! grounded with ice 112 112 if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.)) then ! grounded with ice 113 113 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 … … 141 141 142 142 ! bilan d'eau sur la grille : 143 144 143 water_bilan=sum(tot_water(:,:)) 144 diff_H = diff_H/dtt 145 145 146 146 !999 format(f0.2,1x,e15.8,1x,i10,8(1x,e15.8)) 147 147 ! 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 148 148 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) 149 149 150 150 endif -
trunk/SOURCES/calving_frange.f90
r224 r237 89 89 90 90 Hcoup(:,:) = min ( max( & 91 (-(Bsoc0(:,:)-sealevel ) - prof_plateau)/(prof_abysses-prof_plateau) &91 (-(Bsoc0(:,:)-sealevel_2d(:,:)) - prof_plateau)/(prof_abysses-prof_plateau) & 92 92 *(hcoup_abysses-hcoup_plateau)+hcoup_plateau & 93 93 , hcoup_plateau), hcoup_abysses ) … … 121 121 122 122 Hcoup(:,:) = min ( max( & 123 (-(Bsoc(:,:)-sealevel ) - prof_plateau)/(prof_abysses-prof_plateau) &123 (-(Bsoc(:,:)-sealevel_2d(:,:)) - prof_plateau)/(prof_abysses-prof_plateau) & 124 124 *(hcoup_abysses-hcoup_plateau)+hcoup_plateau & 125 125 , hcoup_plateau), hcoup_abysses ) … … 409 409 calv(i,j)=-h(i,j) 410 410 !cdc H(i,j)=1. 411 !cdc 1m H(i,j)=min(1.,max(0.,(sealevel - Bsoc(i,j))*row/ro-0.01))411 !cdc 1m H(i,j)=min(1.,max(0.,(sealevel_2d(i,j) - Bsoc(i,j))*row/ro-0.01)) 412 412 H(i,j)=0. 413 S(i,j)=H(i,j)*(1.-ro/row) + sealevel 413 S(i,j)=H(i,j)*(1.-ro/row) + sealevel_2d(i,j) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 414 414 B(i,j)=S(i,j) - H(i,j) 415 415 ! ATTENTION ne pas mettre ice=0 sinon degradation bilan d'eau (bm et bmelt non comptabilises dans ce cas) … … 439 439 ice(:,:)=0 440 440 H(:,:)=0. 441 S(:,:)=H(:,:)*(1.-ro/row) + sealevel 441 S(:,:)=H(:,:)*(1.-ro/row) + sealevel_2d(:,:) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 442 442 B(:,:)=S(:,:) - H(:,:) 443 443 endwhere -
trunk/SOURCES/calving_frange_abuk.f90
r204 r237 96 96 97 97 Hcoup(:,:) = min ( max( & 98 (-(Bsoc0(:,:)-sealevel ) - prof_plateau)/(prof_abysses-prof_plateau) &98 (-(Bsoc0(:,:)-sealevel_2d(:,:)) - prof_plateau)/(prof_abysses-prof_plateau) & 99 99 *(hcoup_abysses-hcoup_plateau)+hcoup_plateau & 100 100 , hcoup_plateau), hcoup_abysses ) … … 128 128 129 129 Hcoup(:,:) = min ( max( & 130 (-(Bsoc(:,:)-sealevel ) - prof_plateau)/(prof_abysses-prof_plateau) &130 (-(Bsoc(:,:)-sealevel_2d(:,:)) - prof_plateau)/(prof_abysses-prof_plateau) & 131 131 *(hcoup_abysses-hcoup_plateau)+hcoup_plateau & 132 132 , hcoup_plateau), hcoup_abysses ) … … 416 416 calv(i,j)=-h(i,j) 417 417 !cdc H(i,j)=1. 418 !cdc 1m H(i,j)=min(1.,max(0.,(sealevel - Bsoc(i,j))*row/ro-0.01))418 !cdc 1m H(i,j)=min(1.,max(0.,(sealevel_2d(i,j) - Bsoc(i,j))*row/ro-0.01)) 419 419 H(i,j)=0. 420 S(i,j)=H(i,j)*(1.-ro/row) + sealevel 420 S(i,j)=H(i,j)*(1.-ro/row) + sealevel_2d(i,j) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 421 421 B(i,j)=S(i,j) - H(i,j) 422 422 ! ATTENTION ne pas mettre ice=0 sinon degradation bilan d'eau (bm et bmelt non comptabilises dans ce cas) … … 446 446 !~ ice(:,:)=0 447 447 !~ H(:,:)=0. 448 !~ S(:,:)=H(:,:)*(1.-ro/row) + sealevel 448 !~ S(:,:)=H(:,:)*(1.-ro/row) + sealevel_2d(:,:) 449 449 !~ B(:,:)=S(:,:) - H(:,:) 450 450 !~ endwhere … … 455 455 ice(:,:)=0 456 456 H(:,:)=0. 457 S(:,:)=H(:,:)*(1.-ro/row) + sealevel 457 S(:,:)=H(:,:)*(1.-ro/row) + sealevel_2d(:,:) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 458 458 B(:,:)=S(:,:) - H(:,:) 459 459 endwhere -
trunk/SOURCES/climat-forcage-insolation_mod.f90
r65 r237 12 12 ! C. Dumas 07/2015 13 13 14 use module3D_phy,only:nx,ny,S,H,flot, slv,Tann,Tjuly,Tmois,acc,coefbmshelf,ro,num_param,num_rep_42,dirnameinp,time,pi,dx,xlong,ylat,bsoc14 use module3D_phy,only:nx,ny,S,H,flot,Tann,Tjuly,Tmois,acc,coefbmshelf,sealevel_2d,ro,num_param,num_rep_42,dirnameinp,time,pi,dx,xlong,ylat,bsoc 15 15 use netcdf 16 16 use io_netcdf_grisli … … 637 637 do j=1,ny 638 638 do i=1,nx 639 Zs(i,j)=max(s lv(i,j),S(i,j))639 Zs(i,j)=max(sealevel_2d(i,j),S(i,j)) 640 640 do ictr=CO2SnapMin,CO2SnapMax,1 641 641 do igtr=1,gtr -
trunk/SOURCES/climat-forcage-insolation_mod_oneway.f90
r108 r237 12 12 ! C. Dumas 06/2015 13 13 14 USE module3D_phy,only:nx,ny,S,s lv,Tann,Tjuly,Tmois,acc,coefbmshelf,ro,num_param,num_rep_42,dirnameinp,time14 USE module3D_phy,only:nx,ny,S,sealevel_2d,Tann,Tjuly,Tmois,acc,coefbmshelf,ro,num_param,num_rep_42,dirnameinp,time 15 15 !use interface_input 16 16 use netcdf … … 245 245 do j=1,ny 246 246 do i=1,nx 247 Zs(i,j)=max(s lv(i,j),S(i,j))247 Zs(i,j)=max(sealevel_2d(i,j),S(i,j)) 248 248 !Il faut mettre S0(i,j) si pas d'interpolation sinon Ssnap(i,j,ictr,igtr) 249 249 !if (Zs(i,j).ge.Ssnap(i,j,ictr,igtr)) then -
trunk/SOURCES/climat-perturb_mod-0.4.f90
r93 r237 13 13 14 14 15 use module3d_phy,only:nx,ny,S,S0,Tann,Tjuly,precip,acc,Ylat,num_forc,num_param,num_rep_42,dirforcage,dirnameinp,tafor,time,sealevel, coefbmshelf15 use module3d_phy,only:nx,ny,S,S0,Tann,Tjuly,precip,acc,Ylat,num_forc,num_param,num_rep_42,dirforcage,dirnameinp,tafor,time,sealevel,sealevel_2d,coefbmshelf 16 16 use netcdf 17 17 use io_netcdf_grisli … … 307 307 100 continue 308 308 309 sealevel_2d(:,:)=sealevel 309 310 310 311 ! coefmshelf est un coefficient qui fait varier bmgrz et bmshelf -
trunk/SOURCES/climat_forcage_mod.f90
r224 r237 6 6 ! possibilite d'utiliser autant de snapshots que l'on veut (ntr) 7 7 ! fichiers de forcage, Snapshot et valeurs de lapse rate definis dans fichier param 8 use module3d_phy,only:nx,ny,sealevel,s lv,S,S0,Tmois,Tann,Tjuly,acc,num_forc,num_param,num_rep_42,coefbmshelf,PI,ro,time,dirforcage,dirnameinp8 use module3d_phy,only:nx,ny,sealevel,sealevel_2d,S,S0,Tmois,Tann,Tjuly,acc,num_forc,num_param,num_rep_42,coefbmshelf,PI,ro,time,dirforcage,dirnameinp 9 9 use netcdf 10 10 use io_netcdf_grisli … … 279 279 100 continue 280 280 281 s lv(:,:)=sealevel281 sealevel_2d(:,:)=sealevel 282 282 283 283 if(timeloc.le.ttr(1)) then ! time en dehors des limites du fichier forcage … … 333 333 ! calcul de Tjuly et et Tann 334 334 !$OMP WORKSHARE 335 Zs(:,:)=max(s lv(:,:),S(:,:))335 Zs(:,:)=max(sealevel_2d(:,:),S(:,:)) 336 336 deltaZs(:,:)= -lapserate*(Zs(:,:)-S0(:,:)) ! difference de temperature entre Zs et S0 337 337 !$OMP END WORKSHARE -
trunk/SOURCES/eaubasale-0.5_mod.f90
r196 r237 158 158 if (flot(i,j))then ! points flottants 159 159 klimit(i,j)=1 160 limit_hw(i,j)=(sealevel -Bsoc(i,j))*rowg/rofreshg160 limit_hw(i,j)=(sealevel_2d(i,j)-Bsoc(i,j))*rowg/rofreshg 161 161 162 162 else if (IBASE(I,J).eq.1) then ! base froide … … 200 200 pot_w(I,J)=pot_w(I,J)+rog*H(I,J) 201 201 202 pot_f(I,J)=rowg*(sealevel -S(i,j)+H(I,J)) ! pression a la base de l'ice shelf202 pot_f(I,J)=rowg*(sealevel_2d(i,j)-S(i,j)+H(I,J)) ! pression a la base de l'ice shelf 203 203 enddo 204 204 enddo … … 314 314 315 315 if (flot(i,j)) then ! if flot > hwater=hwater in ocean 316 hwater(i,j)= sealevel - bsoc(i,j)316 hwater(i,j)= sealevel_2d(i,j) - bsoc(i,j) 317 317 ! hwater(i,j)= max(0.,hwater(i,j)) 318 318 if (hwater(i,j).lt.0.) hwater(i,j)=0. -
trunk/SOURCES/flottab2-0.7.f90
r223 r237 122 122 where ( mk_init(:,:).eq.1) ! pose 123 123 flot(:,:) = .False. 124 H(:,:)=max(H(:,:),(10.+sealevel -Bsoc(:,:))*row/ro) ! pour avoir archim=10124 H(:,:)=max(H(:,:),(10.+sealevel_2d(:,:)-Bsoc(:,:))*row/ro) ! pour avoir archim=10 125 125 S(:,:) = Bsoc(:,:) + H(:,:) 126 126 elsewhere … … 175 175 uys1(i,j)=uybar(i,j) 176 176 177 archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel 177 archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel_2d(i,j) 178 178 ! if ((i.eq.132).and.(j.eq.183)) print*,'archim=',archim 179 179 … … 189 189 if (igrdline.eq.1) then ! en cas de grounding line prescrite 190 190 flot(i,j)=.false. 191 H(i,j)=(10.+sealevel -Bsoc(i,j))*row/ro ! pour avoir archim=10191 H(i,j)=(10.+sealevel_2d(i,j)-Bsoc(i,j))*row/ro ! pour avoir archim=10 192 192 new_flot_point(i,j)=.false. 193 193 endif … … 202 202 if (igrdline.eq.1) then ! en cas de grounding line prescrite 203 203 flot(i,j)=.false. 204 H(i,j)=(10.+sealevel -Bsoc(i,j))*row/ro ! pour avoir archim=10204 H(i,j)=(10.+sealevel_2d(i,j)-Bsoc(i,j))*row/ro ! pour avoir archim=10 205 205 new_flot_point(i,j)=.false. 206 206 endif … … 228 228 flot(i,j)=.true. !cdc points ocean sont flot meme sans glace 229 229 H(i,j)=0. 230 S(i,j)=sealevel 230 S(i,j)=sealevel_2d(i,j) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 231 231 B(i,j)=S(i,j)-H(i,j) 232 232 endif arch … … 239 239 240 240 surnet=H(i,j)*(1.-ro/row) 241 S(i,j)=surnet+sealevel 241 S(i,j)=surnet+sealevel_2d(i,j) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 242 242 B(i,j)=S(i,j)-H(i,j) 243 243 end if … … 246 246 end do 247 247 !$OMP END DO 248 !~ print*,'flottab 2',S(132,183),H(132,183),BSOC(132,183),B(132,183),sealevel249 !~ print*,'flottab 2',flot(132,183),ice(132,183)250 !~ if (S(132,183).LT.sealevel) print*,'BUGGGGGGGGGGGGG !!!!!!!!!!!!!!!!'251 248 252 249 !!$ do i=1,nx … … 262 259 263 260 264 265 261 !----------------------------------------------------------------------- 266 262 !$OMP DO … … 298 294 ! sur le demi noeud condition archim < 100 m 299 295 300 archim=(Bsoc(i,j)+Bsoc(i-1,j))*0.5- sealevel+ro/row*Hmx(i,j)296 archim=(Bsoc(i,j)+Bsoc(i-1,j))*0.5-(sealevel_2d(i,j)+sealevel_2d(i-1,j))*0.5+ro/row*Hmx(i,j) 301 297 gzmx(i,j)=gzmx(i,j).and.(archim.le.100.) 302 298 cotemx(i,j)=gzmx(i,j) … … 337 333 ! sur le demi noeud condition archim > 100 m 338 334 339 archim=(Bsoc(i,j)+Bsoc(i,j-1))*0.5- sealevel+ro/row*Hmy(i,j)335 archim=(Bsoc(i,j)+Bsoc(i,j-1))*0.5-(sealevel_2d(i,j)+sealevel_2d(i,j-1))*0.5+ro/row*Hmy(i,j) 340 336 gzmy(i,j)=gzmy(i,j).and.(archim.le.100.) 341 337 cotemy(i,j)=gzmy(i,j) … … 605 601 ice(:,:)=0 606 602 H(:,:)=0. 607 S(:,:)=H(:,:)*(1.-ro/row) + sealevel 603 S(:,:)=H(:,:)*(1.-ro/row) + sealevel_2d(:,:) 608 604 B(:,:)=S(:,:) - H(:,:) 609 605 end where … … 790 786 flgzmy(smax_i,smax_j)=.false. ; flgzmx(smax_i,smax_j+1)=.false. 791 787 792 if (Smax_.le.sealevel) then 788 if (Smax_.le.sealevel) then ! afq -- WARNING: en fait avec le niveau marin local cest plus complique que ca! 789 !if (Smax_.le.sealevel_2d(i,j)) then ! afq --WARNING: il faudrait regarder point par point sur la tache... 793 790 write(num_tracebug,*)'Attention, une ile avec la surface sous l eau' 794 791 write(num_tracebug,*)'time=',time,' coord:',smax_i,smax_j -
trunk/SOURCES/furst_schoof_mod.f90
r228 r237 85 85 gr_line_schoof(:,:)=0 86 86 87 archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot -sealevel 87 archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot -sealevel_2d(:,:) 88 88 89 89 imx_diag_in(:,:) = imx_diag(:,:) … … 93 93 do i=3,nx-3 94 94 !archim = Bsoc(i,j)+H(i,j)*ro/row -sealevel 95 if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel )) then ! grounded and with ice95 if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel_2d(i,j))) then ! grounded and with ice 96 96 ! print*,'schoof gr line ij',i,j 97 97 ! selon x (indice i) … … 100 100 !if (flot(i-1,j).and.Uxbar(i,j).LT.0.) then ! grounding line between i-1,i CAS WEST (ecoulement vers grille West) 101 101 102 call calc_xgl4schoof(-dx,alpha_flot,H(i,j),Bsoc(i,j), H(i-1,j),Bsoc(i-1,j),Sealevel,xpos,Hgl)102 call calc_xgl4schoof(-dx,alpha_flot,H(i,j),Bsoc(i,j),sealevel_2d(i,j),H(i-1,j),Bsoc(i-1,j),sealevel_2d(i-1,j),xpos,Hgl) 103 103 xpos_tab(i,j)=xpos 104 104 Hglx_tab(i,j)=Hgl … … 194 194 !else if (flot(i+1,j).and.Uxbar(i+1,j).GT.0.) then ! grounding line between i,i+1 CAS EST (ecoulement vers grille Est) 195 195 196 call calc_xgl4schoof(dx,alpha_flot,H(i,j),Bsoc(i,j), H(i+1,j),Bsoc(i+1,j),Sealevel,xpos,Hgl)196 call calc_xgl4schoof(dx,alpha_flot,H(i,j),Bsoc(i,j),sealevel_2d(i,j),H(i+1,j),Bsoc(i+1,j),sealevel_2d(i+1,j),xpos,Hgl) 197 197 xpos_tab(i,j)=xpos 198 198 Hglx_tab(i,j)=Hgl … … 285 285 286 286 287 call calc_xgl4schoof(-dy,alpha_flot,H(i,j),Bsoc(i,j), H(i,j-1),Bsoc(i,j-1),Sealevel,ypos,Hgl)287 call calc_xgl4schoof(-dy,alpha_flot,H(i,j),Bsoc(i,j),sealevel_2d(i,j),H(i,j-1),Bsoc(i,j-1),sealevel_2d(i,j-1),ypos,Hgl) 288 288 ypos_tab(i,j)=ypos 289 289 Hgly_tab(i,j)=Hgl … … 373 373 !else if (flot(i,j+1).and.Uybar(i,j+1).GT.0.) then ! grounding line between j,j+1 CAS NORD (ecoulement vers grille Nord) 374 374 375 call calc_xgl4schoof(dy,alpha_flot,H(i,j),Bsoc(i,j), H(i,j+1),Bsoc(i,j+1),Sealevel,ypos,Hgl)375 call calc_xgl4schoof(dy,alpha_flot,H(i,j),Bsoc(i,j),sealevel_2d(i,j),H(i,j+1),Bsoc(i,j+1),sealevel_2d(i,j+1),ypos,Hgl) 376 376 ypos_tab(i,j)=ypos 377 377 Hgly_tab(i,j)=Hgl … … 507 507 !------------------------------------------------------------------------- 508 508 ! calcule la position sous maille de la ligne dechouage 509 subroutine calc_xgl4schoof(dy,alpha,H_0,B_0, H_1,B_1,SL,xpos,Hgl)509 subroutine calc_xgl4schoof(dy,alpha,H_0,B_0,SL_0,H_1,B_1,SL_1,xpos,Hgl) 510 510 511 511 … … 516 516 real,intent(in) :: H_0 !< epaisseur au point pose 517 517 real,intent(in) :: B_0 !< altitude socle au point pose 518 real,intent(in) :: SL_0 !< sea level au point pose 518 519 real,intent(in) :: H_1 !< epaisseur au point flottant 519 520 real,intent(in) :: B_1 !< altitude socle au point flottant 520 real,intent(in) :: SL !< Sea level521 real,intent(in) :: SL_1 !< sea level au point flottant 521 522 real,intent(out) :: xpos !< position de la ligne (en distance depuis le point pose) 522 real,intent(out) :: Hgl !< Epaisseur de glace a la grounding line523 real,intent(out) :: Hgl !< epaisseur de glace a la grounding line 523 524 524 525 … … 541 542 !~ endif 542 543 543 Cgl = (B_1-B_0) + alpha * (H_1-H_0) 544 Cgl = (B_1-B_0) + alpha * (H_1-H_0) - (SL_1-SL_0) 544 545 545 546 if (abs(Cgl).gt.1.e-5) then 546 xpos = (SL - (B_0 + alpha * H_0)) / Cgl547 xpos = (SL_0 - (B_0 + alpha * H_0)) / Cgl 547 548 else 548 549 xpos = 0.5 ! verifier … … 558 559 print*,'calc_xgl4schoof : xpos value error i,j=',i,j 559 560 print*, 'xpos,Cgl', xpos,Cgl 560 print*,'B_0, alpha, H_0, SL ', B_0, alpha, H_0, SL561 print*,'archim=',B_1+H_1*alpha - SL 561 print*,'B_0, alpha, H_0, SL_0', B_0, alpha, H_0, SL_0 562 print*,'archim=',B_1+H_1*alpha - SL_1 562 563 !stop 563 564 xpos = min(max(0.,xpos),1.) -
trunk/SOURCES/initial-phy-2-job.f90
r10 r237 192 192 TJFOR=0.0 193 193 SEALEVEL=0.0 194 sealevel_2d(:,:)=0. 194 195 195 196 SECYEAR=365.*24.*3600. -
trunk/SOURCES/initial-phy-2.f90
r233 r237 185 185 TJFOR=0.0 186 186 SEALEVEL=0.0 187 sealevel_2d(:,:)=0. 187 188 188 189 SECYEAR=365.*24.*3600. -
trunk/SOURCES/initial2-0.4.f90
r4 r237 71 71 do J=1,NY 72 72 do I=1,NX 73 if ((BSOC(I,J)+H(I,J)*RO/ROW -SEALEVEL ).LT.0.) then73 if ((BSOC(I,J)+H(I,J)*RO/ROW -SEALEVEL_2D(I,J)).LT.0.) then 74 74 FLOT(I,J)=.TRUE. 75 75 else -
trunk/SOURCES/isostasie_mod-0.3.f90
r181 r237 103 103 do i=1,nx 104 104 do j=1,ny 105 if (ro*h0(i,j).ge.-row*(BSOC0(i,j)-sealevel )) then105 if (ro*h0(i,j).ge.-row*(BSOC0(i,j)-sealevel_2d(i,j))) then 106 106 107 107 ! glace ou terre … … 110 110 else 111 111 ! ocean 112 charge(i,j)=-rowg*(Bsoc0(i,j)-sealevel )112 charge(i,j)=-rowg*(Bsoc0(i,j)-sealevel_2d(i,j)) 113 113 endif 114 114 end do -
trunk/SOURCES/no_lakes.f90
r4 r237 31 31 subroutine lake_to_slv 32 32 33 s lv(:,:)=sealevel33 sealevel_2d(:,:)=sealevel 34 34 35 35 end subroutine lake_to_slv -
trunk/SOURCES/prescribe-H-i2s_mod.f90
r124 r237 198 198 199 199 where (MK_gl0(:,:) .eq. 1) 200 hp(:,:) = (sealevel -Bsoc(:,:))*row/ro +1.200 hp(:,:) = (sealevel_2d(:,:)-Bsoc(:,:))*row/ro +1. 201 201 end where 202 202 -
trunk/SOURCES/steps_time_loop.f90
r182 r237 68 68 where(.not.flot(:,:)) 69 69 S(:,:)=Bsoc(:,:)+H(:,:) 70 S(:,:)=max(S(:,:),sealevel )70 S(:,:)=max(S(:,:),sealevel_2d(:,:)) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 71 71 B(:,:)=Bsoc(:,:) 72 72 end where -
trunk/SOURCES/steps_time_loop_avec_iterbeta.f90
r180 r237 69 69 where(.not.flot(:,:)) 70 70 S(:,:)=Bsoc(:,:)+H(:,:) 71 S(:,:)=max(S(:,:),sealevel )71 S(:,:)=max(S(:,:),sealevel_2d(:,:)) !afq -- WARNING: est-ce qu'on veut vraiment mettre S a la valeur locale du niveau marin? 72 72 B(:,:)=Bsoc(:,:) 73 73 end where -
trunk/SOURCES/taubed-0.3.f90
r181 r237 60 60 do J=1,NY 61 61 do I=1,NX 62 if (RO*H(I,J).ge.-ROW*(BSOC(I,J)-SEALEVEL )) then62 if (RO*H(I,J).ge.-ROW*(BSOC(I,J)-SEALEVEL_2D(i,j))) then 63 63 ! glace ou terre 64 64 CHARGE(I,J)=ROG*H(I,J) 65 65 else 66 66 ! ocean 67 CHARGE(I,J)=-ROWG*(BSOC(I,J)-SEALEVEL )67 CHARGE(I,J)=-ROWG*(BSOC(I,J)-SEALEVEL_2D(i,j)) 68 68 endif 69 69 end do … … 102 102 do J=1,NY 103 103 do I=1,NX 104 if (RO*H(I,J).ge.-ROW*(BSOC(I,J)-SEALEVEL )) then104 if (RO*H(I,J).ge.-ROW*(BSOC(I,J)-SEALEVEL_2D(i,j))) then 105 105 ! glace ou terre 106 106 W1(I,J)=RO/ROM*H(I,J) 107 107 else 108 108 ! ocean 109 W1(I,J)=-ROW/ROM*(BSOC(I,J)-SEALEVEL )109 W1(I,J)=-ROW/ROM*(BSOC(I,J)-SEALEVEL_2D(i,j)) 110 110 endif 111 111 end do
Note: See TracChangeset
for help on using the changeset viewer.