Changeset 14644 for NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/OCE/stpmlf.F90
- Timestamp:
- 2021-03-26T15:33:49+01:00 (3 years ago)
- Location:
- NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/stpmlf.F90 (modified) (19 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final
- Property svn:externals
-
old new 9 9 10 10 # SETTE 11 ^/utils/CI/sette _wave@13990sette11 ^/utils/CI/sette@14244 sette
-
- Property svn:externals
-
NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/OCE/stpmlf.F90
r14286 r14644 35 35 !! 4.x ! 2020-08 (S. Techene, G. Madec) quasi eulerian coordinate time stepping 36 36 !!---------------------------------------------------------------------- 37 38 37 #if defined key_qco || defined key_linssh 39 38 !!---------------------------------------------------------------------- … … 42 41 !! 'key_linssh Fixed in time vertical coordinate 43 42 !!---------------------------------------------------------------------- 44 43 !! 45 44 !!---------------------------------------------------------------------- 46 45 !! stp_MLF : NEMO modified Leap Frog time-stepping with qco or linssh … … 49 48 ! 50 49 USE domqco ! quasi-eulerian coordinate 51 USE traatf_qco ! time filtering (tra_atf_qco routine) 52 USE dynatf_qco ! time filtering (dyn_atf_qco routine) 53 54 USE bdydyn ! ocean open boundary conditions (define bdy_dyn) 55 56 #if defined key_agrif 57 USE agrif_oce_interp 58 #endif 50 USE traatf_qco ! time filtering (tra_atf_qco routine) 51 USE dynatf_qco ! time filtering (dyn_atf_qco routine) 59 52 60 53 IMPLICIT NONE 61 54 PRIVATE 62 55 63 56 PUBLIC stp_MLF ! called by nemogcm.F90 64 57 … … 67 60 68 61 !! * Substitutions 62 # include "do_loop_substitute.h90" 69 63 # include "domzgr_substitute.h90" 70 64 # include "single_precision_substitute.h90" 65 # include "do_loop_substitute.h90" 71 66 !!---------------------------------------------------------------------- 72 67 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 86 81 !! *** ROUTINE stp_MLF *** 87 82 !! 88 !! ** Purpose : - Time stepping of O PA(momentum and active tracer eqs.)83 !! ** Purpose : - Time stepping of OCE (momentum and active tracer eqs.) 89 84 !! - Time stepping of SI3 (dynamic and thermodynamic eqs.) 90 85 !! - Time stepping of TRC (passive tracer eqs.) … … 99 94 !! -8- Outputs and diagnostics 100 95 !!---------------------------------------------------------------------- 101 INTEGER :: ji, jj, jk ! dummy loop indice 102 INTEGER :: indic ! error indicator if < 0 103 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgdept 104 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zssh_f 96 INTEGER :: ji, jj, jk, jtile ! dummy loop indice 97 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zgdept 105 98 !! --------------------------------------------------------------------- 106 99 #if defined key_agrif 100 IF( nstop > 0 ) RETURN ! avoid to go further if an error was detected during previous time step (child grid) 107 101 kstp = nit000 + Agrif_Nb_Step() 108 102 Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs ! agrif_oce module copies of time level indices … … 112 106 ENDIF 113 107 IF( kstp == nit000 + 1 ) lk_agrif_fstep = .FALSE. 114 # if defined key_ iomput108 # if defined key_xios 115 109 IF( Agrif_Nbstepint() == 0 ) CALL iom_swap( cxios_context ) 116 110 # endif … … 131 125 ! update I/O and calendar 132 126 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 133 indic = 0 ! reset to no error condition 134 127 ! 135 128 IF( kstp == nit000 ) THEN ! initialize IOM context (must be done after nemo_init for AGRIF+XIOS+OASIS) 136 CALL iom_init( cxios_context, ld_closedef=.FALSE. ) ! for model grid (including p assible AGRIF zoom)129 CALL iom_init( cxios_context, ld_closedef=.FALSE. ) ! for model grid (including possible AGRIF zoom) 137 130 IF( lk_diamlr ) CALL dia_mlr_iom_init ! with additional setup for multiple-linear-regression analysis 138 131 CALL iom_init_closedef 139 132 IF( ln_crs ) CALL iom_init( TRIM(cxios_context)//"_crs" ) ! for coarse grid 133 ENDIF 134 IF( kstp == nitrst .AND. lwxios ) THEN 135 CALL iom_swap( cw_ocerst_cxt ) 136 CALL iom_init_closedef( cw_ocerst_cxt ) 137 CALL iom_setkt( kstp - nit000 + 1, cw_ocerst_cxt ) 138 #if defined key_top 139 CALL iom_swap( cw_toprst_cxt ) 140 CALL iom_init_closedef( cw_toprst_cxt ) 141 CALL iom_setkt( kstp - nit000 + 1, cw_toprst_cxt ) 142 #endif 143 ENDIF 144 IF( kstp + nn_fsbc - 1 == nitrst .AND. lwxios ) THEN 145 #if defined key_si3 146 CALL iom_swap( cw_icerst_cxt ) 147 CALL iom_init_closedef( cw_icerst_cxt ) 148 CALL iom_setkt( kstp - nit000 + 1, cw_icerst_cxt ) 149 #endif 150 IF( ln_abl ) THEN 151 CALL iom_swap( cw_ablrst_cxt ) 152 CALL iom_init_closedef( cw_ablrst_cxt ) 153 CALL iom_setkt( kstp - nit000 + 1, cw_ablrst_cxt ) 154 ENDIF 140 155 ENDIF 141 156 IF( kstp /= nit000 ) CALL day( kstp ) ! Calendar (day was already called at nit000 in day_init) … … 155 170 ! Update stochastic parameters and random T/S fluctuations 156 171 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 157 IF( ln_sto_eos ) CALL sto_par( kstp )! Stochastic parameters158 IF( ln_sto_eos ) CALL sto_pts( ts(:,:,:,:,Nnn) )! Random T/S fluctuations172 IF( ln_sto_eos ) CALL sto_par( kstp ) ! Stochastic parameters 173 IF( ln_sto_eos ) CALL sto_pts( ts(:,:,:,:,Nnn) ) ! Random T/S fluctuations 159 174 160 175 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> … … 195 210 ! Ocean dynamics : hdiv, ssh, e3, u, v, w 196 211 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 197 DO jk = 1, jpk 198 zgdept(:,:,jk) = gdept(:,:,jk,Nnn) 199 END DO 200 CALL ssh_nxt ( kstp, Nbb, Nnn, ssh, Naa ) ! after ssh (includes call to div_hor) 201 IF( .NOT.lk_linssh ) THEN 202 CALL dom_qco_r3c( CASTWP(ssh(:,:,Naa)), r3t(:,:,Naa), r3u(:,:,Naa), r3v(:,:,Naa) ) ! "after" ssh/h_0 ratio at t,u,v pts 203 IF( ln_dynspg_exp ) CALL dom_qco_r3c( CASTWP(ssh(:,:,Nnn)), r3t(:,:,Nnn), r3u(:,:,Nnn), r3v(:,:,Nnn), r3f(:,:) ) ! spg_exp : needed only for "now" ssh/h_0 ratio at f point 204 ENDIF 205 CALL wzv ( kstp, Nbb, Nnn, Naa, ww ) ! Nnn cross-level velocity 206 IF( ln_zad_Aimp ) CALL wAimp ( kstp, Nnn ) ! Adaptive-implicit vertical advection partitioning 207 CALL eos ( CASTWP(ts(:,:,:,:,Nnn)), rhd, rhop, zgdept ) ! now in situ density for hpg computation 208 212 213 CALL ssh_nxt ( kstp, Nbb, Nnn, ssh, Naa ) ! after ssh (includes call to div_hor) 214 IF( .NOT.lk_linssh ) THEN 215 CALL dom_qco_r3c( CASTWP(ssh(:,:,Naa)), r3t(:,:,Naa), r3u(:,:,Naa), r3v(:,:,Naa) ) ! "after" ssh/h_0 ratio at t,u,v pts 216 IF( ln_dynspg_exp ) & 217 & CALL dom_qco_r3c( CASTWP(ssh(:,:,Nnn)), r3t(:,:,Nnn), r3u(:,:,Nnn), r3v(:,:,Nnn), r3f(:,:) ) ! spg_exp : needed only for "now" ssh/h_0 ratio at f point 218 ENDIF 219 CALL wzv ( kstp, Nbb, Nnn, Naa, ww ) ! Nnn cross-level velocity 220 IF( ln_zad_Aimp ) CALL wAimp ( kstp, Nnn ) ! Adaptive-implicit vertical advection partitioning 221 ALLOCATE( zgdept(jpi,jpj,jpk) ) 222 DO jk = 1, jpk 223 zgdept(:,:,jk) = gdept(:,:,jk,Nnn) 224 END DO 225 CALL eos ( CASTWP(ts(:,:,:,:,Nnn)), rhd, rhop, zgdept ) ! now in situ density for hpg computation 226 DEALLOCATE( zgdept ) 209 227 210 228 uu(:,:,:,Nrhs) = 0._wp ! set dynamics trends to zero … … 224 242 CALL dyn_hpg( kstp, Nnn , uu, vv, Nrhs ) ! horizontal gradient of Hydrostatic pressure 225 243 CALL dyn_spg( kstp, Nbb, Nnn, Nrhs, uu, vv, ssh, uu_b, vv_b, Naa ) ! surface pressure gradient 226 ! With split-explicit free surface, since now transports have been updated and ssh(:,:,Nrhs) as well227 228 IF( ln_dynspg_ts ) THEN !vertical scale factors and vertical velocity need to be updated244 245 IF( ln_dynspg_ts ) THEN ! With split-explicit free surface, since now transports have been updated and ssh(:,:,Nrhs) 246 ! as well as vertical scale factors and vertical velocity need to be updated 229 247 CALL div_hor ( kstp, Nbb, Nnn ) ! Horizontal divergence (2nd call in time-split case) 230 IF(.NOT.lk_linssh) CALL dom_qco_r3c 248 IF(.NOT.lk_linssh) CALL dom_qco_r3c( CASTWP(ssh(:,:,Naa)), r3t(:,:,Naa), r3u(:,:,Naa), r3v(:,:,Naa), r3f(:,:) ) ! update ssh/h_0 ratio at t,u,v,f pts 231 249 ENDIF 232 250 CALL dyn_zdf ( kstp, Nbb, Nnn, Nrhs, uu, vv, Naa ) ! vertical diffusion … … 271 289 ! Active tracers 272 290 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 273 ts(:,:,:,:,Nrhs) = 0._wp ! set tracer trends to zero 274 275 IF( lk_asminc .AND. ln_asmiau .AND. & 276 & ln_trainc ) CALL tra_asm_inc( kstp, Nbb, Nnn, ts, Nrhs ) ! apply tracer assimilation increment 277 CALL tra_sbc ( kstp, Nnn, ts, Nrhs ) ! surface boundary condition 278 IF( ln_traqsr ) CALL tra_qsr ( kstp, Nnn, ts, Nrhs ) ! penetrative solar radiation qsr 279 IF( ln_isf ) CALL tra_isf ( kstp, Nnn, ts, Nrhs ) ! ice shelf heat flux 280 IF( ln_trabbc ) CALL tra_bbc ( kstp, Nnn, ts, Nrhs ) ! bottom heat flux 281 IF( ln_trabbl ) CALL tra_bbl ( kstp, Nbb, Nnn, ts, Nrhs ) ! advective (and/or diffusive) bottom boundary layer scheme 282 IF( ln_tradmp ) CALL tra_dmp ( kstp, Nbb, Nnn, ts, Nrhs ) ! internal damping trends 283 IF( ln_bdy ) CALL bdy_tra_dmp( kstp, Nbb, ts, Nrhs ) ! bdy damping trends 284 #if defined key_agrif 285 IF(.NOT. Agrif_Root()) & 286 & CALL Agrif_Sponge_tra ! tracers sponge 287 #endif 288 CALL tra_adv ( kstp, Nbb, Nnn, ts, Nrhs ) ! hor. + vert. advection ==> RHS 289 IF( ln_zdfosm ) CALL tra_osm ( kstp, Nnn, ts, Nrhs ) ! OSMOSIS non-local tracer fluxes ==> RHS 290 IF( lrst_oce .AND. ln_zdfosm ) & 291 & CALL osm_rst ( kstp, Nnn, 'WRITE' ) ! write OSMOSIS outputs + ww (so must do here) to restarts 292 CALL tra_ldf ( kstp, Nbb, Nnn, ts, Nrhs ) ! lateral mixing 293 294 CALL tra_zdf ( kstp, Nbb, Nnn, Nrhs, ts, Naa ) ! vertical mixing and after tracer fields 295 IF( ln_zdfnpc ) CALL tra_npc ( kstp, Nnn, Nrhs, ts, Naa ) ! update after fields by non-penetrative convection 296 291 ! Loop over tile domains 292 DO jtile = 1, nijtile 293 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = jtile ) 294 295 DO_3D( 0, 0, 0, 0, 1, jpk ) 296 ts(ji,jj,jk,:,Nrhs) = 0._wp ! set tracer trends to zero 297 END_3D 298 299 IF( lk_asminc .AND. ln_asmiau .AND. & 300 & ln_trainc ) CALL tra_asm_inc( kstp, Nbb, Nnn, ts, Nrhs ) ! apply tracer assimilation increment 301 CALL tra_sbc ( kstp, Nnn, ts, Nrhs ) ! surface boundary condition 302 IF( ln_traqsr ) CALL tra_qsr ( kstp, Nnn, ts, Nrhs ) ! penetrative solar radiation qsr 303 IF( ln_isf ) CALL tra_isf ( kstp, Nnn, ts, Nrhs ) ! ice shelf heat flux 304 IF( ln_trabbc ) CALL tra_bbc ( kstp, Nnn, ts, Nrhs ) ! bottom heat flux 305 IF( ln_trabbl ) CALL tra_bbl ( kstp, Nbb, Nnn, ts, Nrhs ) ! advective (and/or diffusive) bottom boundary layer scheme 306 IF( ln_tradmp ) CALL tra_dmp ( kstp, Nbb, Nnn, ts, Nrhs ) ! internal damping trends 307 IF( ln_bdy ) CALL bdy_tra_dmp( kstp, Nbb, ts, Nrhs ) ! bdy damping trends 308 END DO 309 310 #if defined key_agrif 311 IF(.NOT. Agrif_Root() ) THEN 312 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) 313 CALL Agrif_Sponge_tra ! tracers sponge 314 ENDIF 315 #endif 316 317 ! TEMP: [tiling] Separate loop over tile domains (due to tra_adv workarounds for tiling) 318 DO jtile = 1, nijtile 319 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = jtile ) 320 321 CALL tra_adv ( kstp, Nbb, Nnn, ts, Nrhs ) ! hor. + vert. advection ==> RHS 322 IF( ln_zdfmfc ) CALL tra_mfc ( kstp, Nbb, ts, Nrhs ) ! Mass Flux Convection 323 IF( ln_zdfosm ) THEN 324 CALL tra_osm ( kstp, Nnn, ts, Nrhs ) ! OSMOSIS non-local tracer fluxes ==> RHS 325 IF( lrst_oce ) CALL osm_rst ( kstp, Nnn, 'WRITE' ) ! write OSMOSIS outputs + ww (so must do here) to restarts 326 ENDIF 327 CALL tra_ldf ( kstp, Nbb, Nnn, ts, Nrhs ) ! lateral mixing 328 329 CALL tra_zdf ( kstp, Nbb, Nnn, Nrhs, ts, Naa ) ! vertical mixing and after tracer fields 330 IF( ln_zdfnpc ) CALL tra_npc ( kstp, Nnn, Nrhs, ts, Naa ) ! update after fields by non-penetrative convection 331 END DO 332 333 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Revert to tile over full domain 297 334 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 298 335 ! Set boundary conditions, time filter and swap time levels … … 320 357 r3v(:,:,Nnn) = r3v_f(:,:) 321 358 ENDIF 322 323 359 ! 324 360 ! Swap time levels … … 341 377 #if defined key_agrif 342 378 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 343 ! AGRIF 379 ! AGRIF recursive integration 344 380 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 345 381 Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs ! agrif_oce module copies of time level indices 346 382 CALL Agrif_Integrate_ChildGrids( stp_MLF ) ! allows to finish all the Child Grids before updating 347 383 348 IF( Agrif_NbStepint() == 0 ) THEN 349 CALL Agrif_update_all( ) ! Update all components 350 ENDIF 351 #endif 352 IF( ln_diaobs ) CALL dia_obs ( kstp, Nnn ) ! obs-minus-model (assimilation) diagnostics (call after dynamics update) 353 384 #endif 354 385 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 355 386 ! Control … … 357 388 CALL stp_ctl ( kstp, Nnn ) 358 389 390 #if defined key_agrif 391 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 392 ! AGRIF update 393 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 394 IF( Agrif_NbStepint() == 0 .AND. nstop == 0 ) & 395 & CALL Agrif_update_all( ) ! Update all components 396 397 #endif 398 IF( ln_diaobs .AND. nstop == 0 ) & 399 & CALL dia_obs( kstp, Nnn ) ! obs-minus-model (assimilation) diags (after dynamics update) 400 401 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 402 ! File manipulation at the end of the first time step 403 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 359 404 IF( kstp == nit000 ) THEN ! 1st time step only 360 405 CALL iom_close( numror ) ! close input ocean restart file 406 IF( lrxios ) CALL iom_context_finalize( cr_ocerst_cxt ) 361 407 IF(lwm) CALL FLUSH ( numond ) ! flush output namelist oce 362 408 IF(lwm .AND. numoni /= -1 ) CALL FLUSH ( numoni ) ! flush output namelist ice (if exist) … … 366 412 ! Coupled mode 367 413 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 368 !!gm why lk_oasis and not lk_cpl ???? 369 IF( lk_oasis ) CALL sbc_cpl_snd( kstp, Nbb, Nnn ) ! coupled mode : field exchanges 370 ! 371 #if defined key_iomput 372 IF( kstp == nitend .OR. indic < 0 ) THEN 414 IF( lk_oasis .AND. nstop == 0 ) CALL sbc_cpl_snd( kstp, Nbb, Nnn ) ! coupled mode : field exchanges 415 ! 416 #if defined key_xios 417 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 418 ! Finalize contextes if end of simulation or error detected 419 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 420 IF( kstp == nitend .OR. nstop > 0 ) THEN 373 421 CALL iom_context_finalize( cxios_context ) ! needed for XIOS+AGRIF 374 422 IF( ln_crs ) CALL iom_context_finalize( trim(cxios_context)//"_crs" ) ! … … 377 425 ! 378 426 IF( l_1st_euler ) THEN ! recover Leap-frog timestep 379 rDt = 2._wp * rn_Dt 427 rDt = 2._wp * rn_Dt 380 428 r1_Dt = 1._wp / rDt 381 l_1st_euler = .FALSE. 429 l_1st_euler = .FALSE. 382 430 ENDIF 383 431 ! … … 448 496 !! ** Action : puu(Kaa),pvv(Kaa) after horizontal velocity and tracers 449 497 !!---------------------------------------------------------------------- 498 #if defined key_agrif 499 USE agrif_oce_interp 500 #endif 501 USE bdydyn ! ocean open boundary conditions (define bdy_dyn) 502 !! 450 503 INTEGER , INTENT(in ) :: kt ! ocean time-step index 451 504 INTEGER , INTENT(in ) :: Kbb, Kaa ! before and after time level indices … … 461 514 # endif 462 515 ! ! local domain boundaries (T-point, unchanged sign) 463 CALL lbc_lnk _multi( 'finalize_lbc', puu(:,:,:, Kaa), 'U', -1._wp, pvv(:,:,: ,Kaa), 'V', -1._wp &464 & 516 CALL lbc_lnk( 'finalize_lbc', puu(:,:,:, Kaa), 'U', -1._wp, pvv(:,:,: ,Kaa), 'V', -1._wp & 517 & , pts(:,:,:,jp_tem,Kaa), 'T', 1._wp, pts(:,:,:,jp_sal,Kaa), 'T', 1._wp ) 465 518 ! 466 519 ! !* BDY open boundaries
Note: See TracChangeset
for help on using the changeset viewer.
