Changeset 14179
- Timestamp:
- 2020-12-15T23:17:09+01:00 (4 years ago)
- Location:
- NEMO/trunk/src
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/DIA/diawri.F90
r14143 r14179 215 215 ENDIF 216 216 217 CALL iom_put( "rhop", rhop(:,:,:) ) ! 3D potential density (sigma0)217 IF( .NOT.lk_SWE ) CALL iom_put( "rhop", rhop(:,:,:) ) ! 3D potential density (sigma0) 218 218 219 219 IF ( iom_use("taubot") ) THEN ! bottom stress -
NEMO/trunk/src/OCE/DOM/domqco.F90
r14053 r14179 116 116 ! !== Set of all other vertical scale factors ==! (now and before) 117 117 ! ! Horizontal interpolation of e3t 118 #if defined key_RK3 119 CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb), r3f(:,:) ) 120 CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm) ) 121 #else 118 122 CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) ) 119 123 CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm), r3f(:,:) ) 124 #endif 120 125 ! 121 126 END SUBROUTINE dom_qco_zgr -
NEMO/trunk/src/OCE/IOM/restart.F90
r14145 r14179 47 47 !! * Substitutions 48 48 # include "do_loop_substitute.h90" 49 # include "domzgr_substitute.h90" 49 50 !!---------------------------------------------------------------------- 50 51 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 169 170 ENDIF 170 171 ! 172 #if ! defined key_RK3 171 173 CALL iom_rstput( kt, nitrst, numrow, 'sshn' , ssh(:,: ,Kmm) ) ! now fields 172 174 CALL iom_rstput( kt, nitrst, numrow, 'un' , uu(:,:,: ,Kmm) ) … … 177 179 CALL iom_rstput( kt, nitrst, numrow, 'rhop', rhop ) 178 180 ENDIF 181 #endif 179 182 ENDIF 180 183 … … 264 267 INTEGER :: jk 265 268 REAL(wp), DIMENSION(jpi, jpj, jpk) :: w3d 269 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zgdept ! 3D workspace for QCO 266 270 !!---------------------------------------------------------------------- 267 271 ! … … 279 283 ENDIF 280 284 ! 281 ! !* Read Kmm fields 285 #if defined key_RK3 286 ! !* Read Kbb fields (NB: in RK3 Kmm = Kbb = Nbb) 287 IF(lwp) WRITE(numout,*) ' Kbb u, v and T-S fields read in the restart file' 288 CALL iom_get( numror, jpdom_auto, 'ub' , uu(:,:,: ,Kbb), cd_type = 'U', psgn = -1._wp ) 289 CALL iom_get( numror, jpdom_auto, 'vb' , vv(:,:,: ,Kbb), cd_type = 'V', psgn = -1._wp ) 290 IF( .NOT.lk_SWE ) THEN 291 CALL iom_get( numror, jpdom_auto, 'tb', ts(:,:,:,jp_tem,Kbb) ) 292 CALL iom_get( numror, jpdom_auto, 'sb', ts(:,:,:,jp_sal,Kbb) ) 293 ENDIF 294 #else 295 ! !* Read Kmm fields (MLF only) 282 296 IF(lwp) WRITE(numout,*) ' Kmm u, v and T-S fields read in the restart file' 283 297 CALL iom_get( numror, jpdom_auto, 'un' , uu(:,:,: ,Kmm), cd_type = 'U', psgn = -1._wp ) … … 288 302 ENDIF 289 303 ! 290 IF( l_1st_euler ) THEN !* Euler restart 304 IF( l_1st_euler ) THEN !* Euler restart (MLF only) 291 305 IF(lwp) WRITE(numout,*) ' Kbb u, v and T-S fields set to Kmm values' 292 306 uu(:,:,: ,Kbb) = uu(:,:,: ,Kmm) ! all before fields set to now values … … 294 308 IF( .NOT.lk_SWE ) ts(:,:,:,:,Kbb) = ts(:,:,:,:,Kmm) 295 309 ! 296 ELSE !* Leap frog restart 310 ELSE !* Leap frog restart (MLF only) 297 311 IF(lwp) WRITE(numout,*) ' Kbb u, v and T-S fields read in the restart file' 298 312 CALL iom_get( numror, jpdom_auto, 'ub' , uu(:,:,: ,Kbb), cd_type = 'U', psgn = -1._wp ) … … 303 317 ENDIF 304 318 ENDIF 319 #endif 305 320 ! 306 321 IF( .NOT.lk_SWE ) THEN … … 308 323 CALL iom_get( numror, jpdom_auto, 'rhop' , rhop ) ! now potential density 309 324 ELSE 325 #if defined key_qco 326 ALLOCATE( zgdept(jpi,jpj,jpk) ) 327 DO jk = 1, jpk 328 zgdept(:,:,jk) = gdept(:,:,jk,Kmm) 329 END DO 330 CALL eos( ts(:,:,:,:,Kmm), rhd, rhop, zgdept ) 331 DEALLOCATE( zgdept ) 332 #else 310 333 CALL eos( ts(:,:,:,:,Kmm), rhd, rhop, gdept(:,:,:,Kmm) ) 334 #endif 311 335 ENDIF 312 336 ENDIF … … 347 371 ! !=============================! 348 372 ! 349 ! !* Read ssh at Kmm 373 #if defined key_RK3 374 ! !* RK3: Read ssh at Kbb 375 IF(lwp) WRITE(numout,*) 376 IF(lwp) WRITE(numout,*) ' Kbb sea surface height read in the restart file' 377 CALL iom_get( numror, jpdom_auto, 'sshb' , ssh(:,:,Kbb) ) 378 ! 379 ! !* RK3: Set ssh at Kmm for AGRIF 380 ssh(:,:,Kmm) = 0._wp 381 #else 382 ! !* MLF: Read ssh at Kmm 350 383 IF(lwp) WRITE(numout,*) 351 384 IF(lwp) WRITE(numout,*) ' Kmm sea surface height read in the restart file' 352 385 CALL iom_get( numror, jpdom_auto, 'sshn' , ssh(:,:,Kmm) ) 353 386 ! 354 IF( l_1st_euler ) THEN !* Euler at first time-step387 IF( l_1st_euler ) THEN !* MLF: Euler at first time-step 355 388 IF(lwp) WRITE(numout,*) 356 389 IF(lwp) WRITE(numout,*) ' Euler first time step : ssh(Kbb) = ssh(Kmm)' 357 390 ssh(:,:,Kbb) = ssh(:,:,Kmm) 358 391 ! 359 ELSE !* read ssh at Kbb392 ELSE !* MLF: read ssh at Kbb 360 393 IF(lwp) WRITE(numout,*) 361 394 IF(lwp) WRITE(numout,*) ' Kbb sea surface height read in the restart file' 362 395 CALL iom_get( numror, jpdom_auto, 'sshb', ssh(:,:,Kbb) ) 363 396 ENDIF 397 #endif 364 398 ! !============================! 365 399 ELSE !== Initialize at "rest" ==! … … 390 424 ENDIF 391 425 ! 392 ssh(:,:,Kmm) = ssh(:,:,Kbb) !* set now values from to before ones 426 #if defined key_RK3 427 ssh(:,:,Kmm) = 0._wp !* RK3: set Kmm to 0 for AGRIF 428 #else 429 ssh(:,:,Kmm) = ssh(:,:,Kbb) !* MLF: set now values from to before ones 430 #endif 393 431 ENDIF 394 432 ! 395 433 ! !==========================! 396 ssh(:,:,Kaa) = 0._wp !== Set to 0 for agrif==!434 ssh(:,:,Kaa) = 0._wp !== Set to 0 for AGRIF ==! 397 435 ! !==========================! 398 436 ! -
NEMO/trunk/src/SWE/stprk3.F90
r14143 r14179 129 129 CALL dyn_ldf( kstp, Nbb, Nbb, uu, vv, Nrhs ) ! lateral mixing 130 130 #endif 131 !!st !132 !!st DO_3D( 0,0, 0,0, 1,jpkm1 )133 !!st ! ! horizontal pressure gradient134 !!st uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) - grav * ( ssh(ji+1,jj,Nbb) - ssh(ji,jj,Nbb) ) * r1_e1u(ji,jj)135 !!st vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) - grav * ( ssh(ji,jj+1,Nbb) - ssh(ji,jj,Nbb) ) * r1_e2v(ji,jj)136 !!st END_3D137 !!st !138 !!st #if defined key_RK3all139 !!st ! ! wind stress and layer friction140 !!st z5_6 = 5._wp/6._wp141 !!st DO_3D( 0, 0, 0, 0,1,jpkm1)142 !!st uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) + r1_rho0 * ( z5_6*utau_b(ji,jj) + (1._wp - z5_6)*utau(ji,jj) ) / e3u(ji,jj,jk,Nbb) &143 !!st & - rn_rfr * uu(ji,jj,jk,Nbb)144 !!st vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) + r1_rho0 * ( z5_6*vtau_b(ji,jj) + (1._wp - z5_6)*vtau(ji,jj) ) / e3v(ji,jj,jk,Nbb) &145 !!st & - rn_rfr * vv(ji,jj,jk,Nbb)146 !!st END_3D147 !!st #endif148 !!st why not ?149 131 z5_6 = 5._wp/6._wp 150 132 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) … … 163 145 vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) + zrhs_v 164 146 END_3D 165 !!st end166 147 ! 167 148 ! !== Time stepping of ssh Eq. ==! (and update r3_Naa) … … 229 210 vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) + zrhs_v 230 211 END_3D 231 !!st !232 !!st DO_3D( 0, 0, 0, 0, 1, jpkm1 )233 !!st ! ! horizontal pressure gradient234 !!st uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) - grav * ( ssh(ji+1,jj,Nnn) - ssh(ji,jj,Nnn) ) * r1_e1u(ji,jj)235 !!st vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) - grav * ( ssh(ji,jj+1,Nnn) - ssh(ji,jj,Nnn) ) * r1_e2v(ji,jj)236 !!st END_3D237 !!st !238 !!st #if defined key_RK3all239 !!st ! ! wind stress and layer friction240 !!st z3_4 = 3._wp/4._wp241 !!st DO_3D( 0, 0, 0, 0,1,jpkm1)242 !!st uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) + r1_rho0 * ( z3_4*utau_b(ji,jj) + (1._wp - z3_4)*utau(ji,jj) ) / e3u(ji,jj,jk,Nbb) &243 !!st & - rn_rfr * uu(ji,jj,jk,Nbb)244 !!st vv(ji,jj,jk,Nrhs) = vv(ji,jj,jk,Nrhs) + r1_rho0 * ( z3_4*vtau_b(ji,jj) + (1._wp - z3_4)*vtau(ji,jj) ) / e3v(ji,jj,jk,Nbb) &245 !!st & - rn_rfr * vv(ji,jj,jk,Nbb)246 !!st END_3D247 !!st #endif248 212 ! 249 213 ! !== Time stepping of ssh Eq. ==! (and update r3_Naa) … … 344 308 345 309 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 346 ! diagnostics and outputs 310 ! diagnostics and outputs at Nbb (i.e. the just computed time step) 347 311 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 348 312 349 IF( ln_diacfl ) CALL dia_cfl ( kstp, N nn) ! Courant number diagnostics350 CALL dia_wri ( kstp, N nn) ! ocean model: outputs351 ! 352 IF( lrst_oce ) CALL rst_write ( kstp, Nbb, N nn) ! write output ocean restart file313 IF( ln_diacfl ) CALL dia_cfl ( kstp, Nbb ) ! Courant number diagnostics 314 CALL dia_wri ( kstp, Nbb ) ! ocean model: outputs 315 ! 316 IF( lrst_oce ) CALL rst_write ( kstp, Nbb, Nbb ) ! write output ocean restart file 353 317 354 318 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 355 319 ! Control 356 320 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 357 CALL stp_ctl_SWE ( kstp , N nn)321 CALL stp_ctl_SWE ( kstp , Nbb ) 358 322 359 323 IF( kstp == nit000 ) THEN ! 1st time step only
Note: See TracChangeset
for help on using the changeset viewer.