Changeset 4338 for branches/2013/dev_MERGE_2013/NEMOGCM
- Timestamp:
- 2013-12-17T17:33:09+01:00 (10 years ago)
- Location:
- branches/2013/dev_MERGE_2013/NEMOGCM
- Files:
-
- 1 added
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/LIM_SRC_2/limdia_2.F90
r4325 r4338 242 242 243 243 ! opening "ice_evolu" file 244 CALL ctl_opn( numevo_ice, 'ice_evolu', ' REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )244 CALL ctl_opn( numevo_ice, 'ice_evolu', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 245 245 246 246 !- ecriture de 2 lignes d''entete : -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r4296 r4338 61 61 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdiv_lf ! low frequency part of hz divergence 62 62 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_b, tilde_e3t_n ! baroclinic scale factors 63 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a 63 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a, dtilde_e3t_a ! baroclinic scale factors 64 64 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: frq_rst_e3t ! retoring period for scale factors 65 65 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: frq_rst_hdv ! retoring period for low freq. divergence … … 82 82 IF( ln_vvl_zstar ) dom_vvl_alloc = 0 83 83 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 84 ALLOCATE( tilde_e3t_b(jpi,jpj,jpk) , tilde_e3t_n(jpi,jpj,jpk) , tilde_e3t_a(jpi,jpj,jpk) , & 85 & un_td (jpi,jpj,jpk) , vn_td (jpi,jpj,jpk) , STAT = dom_vvl_alloc ) 84 ALLOCATE( tilde_e3t_b(jpi,jpj,jpk) , tilde_e3t_n(jpi,jpj,jpk) , tilde_e3t_a(jpi,jpj,jpk) , & 85 & dtilde_e3t_a(jpi,jpj,jpk) , un_td (jpi,jpj,jpk) , vn_td (jpi,jpj,jpk) , & 86 & STAT = dom_vvl_alloc ) 86 87 IF( lk_mpp ) CALL mpp_sum ( dom_vvl_alloc ) 87 88 IF( dom_vvl_alloc /= 0 ) CALL ctl_warn('dom_vvl_alloc: failed to allocate arrays') … … 215 216 END DO 216 217 END DO 217 IF( cp_cfg == "orca" .AND. jp_cfg == 2) THEN218 IF( cp_cfg == "orca" .AND. jp_cfg == 3 ) THEN 218 219 ii0 = 103 ; ii1 = 111 ! Suppress ztilde in the Foxe Basin for ORCA2 219 220 ij0 = 128 ; ij1 = 135 ; … … 229 230 230 231 231 SUBROUTINE dom_vvl_sf_nxt( kt )232 SUBROUTINE dom_vvl_sf_nxt( kt, kcall ) 232 233 !!---------------------------------------------------------------------- 233 234 !! *** ROUTINE dom_vvl_sf_nxt *** … … 255 256 !! * Arguments 256 257 INTEGER, INTENT( in ) :: kt ! time step 258 INTEGER, INTENT( in ), OPTIONAL :: kcall ! optional argument indicating call sequence 257 259 !! * Local declarations 258 260 INTEGER :: ji, jj, jk ! dummy loop indices … … 260 262 REAL(wp) :: z2dt ! temporary scalars 261 263 REAL(wp) :: z_tmin, z_tmax ! temporary scalars 264 LOGICAL :: ll_do_bclinic ! temporary logical 262 265 !!---------------------------------------------------------------------- 263 266 IF( nn_timing == 1 ) CALL timing_start('dom_vvl_sf_nxt') … … 271 274 ENDIF 272 275 276 ll_do_bclinic = .TRUE. 277 IF( PRESENT(kcall) ) THEN 278 IF ( kcall == 2 .AND. ln_vvl_ztilde ) ll_do_bclinic = .FALSE. 279 ENDIF 280 273 281 ! ******************************* ! 274 282 ! After acale factors at t-points ! 275 283 ! ******************************* ! 276 284 277 ! ! ----------------- !278 IF( ln_vvl_zstar ) THEN ! z_star coordinate!279 ! !----------------- !280 281 282 283 fse3t_a(:,:,jk) = fse3t_b(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk)284 END DO285 286 ! ! --------------------------- ! 287 ELSEIF( ln_vvl_ztilde .OR. ln_vvl_layer) THEN ! z_tilde or layer coordinate !288 ! ! --------------------------- !285 ! ! --------------------------------------------- ! 286 ! z_star coordinate and barotropic z-tilde part ! 287 ! ! --------------------------------------------- ! 288 289 z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * tmask(:,:,1) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask(:,:,1) ) 290 DO jk = 1, jpkm1 291 ! formally this is the same as fse3t_a = e3t_0*(1+ssha/ht_0) 292 fse3t_a(:,:,jk) = fse3t_b(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 293 END DO 294 295 IF( ln_vvl_ztilde .OR. ln_vvl_layer .AND. ll_do_bclinic ) THEN ! z_tilde or layer coordinate ! 296 ! ! ------baroclinic part------ ! 289 297 290 298 ! I - initialization … … 423 431 tilde_e3t_a(:,:,:) = MAX( tilde_e3t_a(:,:,:), - rn_zdef_max * e3t_0(:,:,:) ) 424 432 425 ! Add "tilda" part to the after scale factor 433 ! 434 ! "tilda" change in the after scale factor 426 435 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 427 fse3t_a(:,:,:) = e3t_0(:,:,:) + tilde_e3t_a(:,:,:) 428 436 DO jk = 1, jpkm1 437 dtilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) - tilde_e3t_b(:,:,jk) 438 END DO 429 439 ! III - Barotropic repartition of the sea surface height over the baroclinic profile 430 440 ! ================================================================================== 431 ! add e3t(n-1) "star" Asselin-filtered 432 DO jk = 1, jpkm1 433 fse3t_a(:,:,jk) = fse3t_a(:,:,jk) + fse3t_b(:,:,jk) - e3t_0(:,:,jk) - tilde_e3t_b(:,:,jk) 434 END DO 435 ! add ( ssh increment + "baroclinicity error" ) proportionnaly to e3t(n) 441 ! add ( ssh increment + "baroclinicity error" ) proportionly to e3t(n) 436 442 ! - ML - baroclinicity error should be better treated in the future 437 443 ! i.e. locally and not spread over the water column. … … 441 447 zht(:,:) = zht(:,:) + tilde_e3t_a(:,:,jk) * tmask(:,:,jk) 442 448 END DO 443 z_scale(:,:) = ( ssha(:,:) - sshb(:,:) - zht(:,:) ) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask(:,:,1) ) 444 DO jk = 1, jpkm1 445 fse3t_a(:,:,jk) = fse3t_a(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 446 END DO 447 448 ENDIF 449 450 IF( ln_vvl_dbg ) THEN ! - ML - test: control prints for debuging 449 z_scale(:,:) = - zht(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask(:,:,1) ) 450 DO jk = 1, jpkm1 451 dtilde_e3t_a(:,:,jk) = dtilde_e3t_a(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 452 END DO 453 454 ENDIF 455 456 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! z_tilde or layer coordinate ! 457 ! ! ---baroclinic part--------- ! 458 DO jk = 1, jpkm1 459 fse3t_a(:,:,jk) = fse3t_a(:,:,jk) + dtilde_e3t_a(:,:,jk) 460 END DO 461 ENDIF 462 463 IF( ln_vvl_dbg .AND. .NOT. ll_do_bclinic ) THEN ! - ML - test: control prints for debuging 451 464 ! 452 465 IF( lwp ) WRITE(numout, *) 'kt =', kt -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r4328 r4338 252 252 ! ---------------------------------------------------- 253 253 IF (lk_dynspg_ts.AND.ln_bt_fw) THEN 254 ! Removeasselin filtering on thicknesses if forward time splitting254 ! No asselin filtering on thicknesses if forward time splitting 255 255 fse3t_b(:,:,:) = fse3t_n(:,:,:) 256 256 ELSE -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90
r4328 r4338 110 110 END DO 111 111 ! ! Sea surface elevation time stepping 112 ! In forward Euler time stepping case, the same formulation as in the leap-frog case can be used 113 ! because emp_b field is initialized with the vlaues of emp field. Hence, 0.5 * ( emp + emp_b ) = emp 112 ! In time-split case we need a first guess of the ssh after (using the baroclinic timestep) in order to 113 ! compute the vertical velocity which can be used to compute the non-linear terms of the momentum equations. 114 ! 114 115 z1_rau0 = 0.5_wp * r1_rau0 115 116 ssha(:,:) = ( sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) ) ) * tmask(:,:,1) … … 197 198 DO jk = 1, jpkm1 198 199 ! horizontal divergence of thickness diffusion transport ( velocity multiplied by e3t) 199 ! - ML - note: computation al lready done in dom_vvl_sf_nxt. Could be optimized (not critical and clearer this way)200 ! - ML - note: computation already done in dom_vvl_sf_nxt. Could be optimized (not critical and clearer this way) 200 201 DO jj = 2, jpjm1 201 202 DO ji = fs_2, fs_jpim1 ! vector opt. -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/step.F90
r4328 r4338 72 72 INTEGER :: jk ! dummy loop indice 73 73 INTEGER :: indic ! error indicator if < 0 74 INTEGER :: kcall ! optional integer argument (dom_vvl_sf_nxt) 74 75 !! --------------------------------------------------------------------- 75 76 … … 106 107 CALL ssh_nxt ( kstp ) ! after ssh (includes call to div_cur) 107 108 IF( lk_dynspg_ts ) THEN 109 IF( lk_vvl ) CALL dom_vvl_sf_nxt( kstp ) ! after vertical scale factors 108 110 CALL wzv ( kstp ) ! now cross-level velocity 109 111 ! In case the time splitting case, update almost all momentum trends here: 110 112 ! Note that the computation of vertical velocity above, hence "after" sea level 111 113 ! is necessary to compute momentum advection for the rhs of barotropic loop: 112 CALL eos ( tsn, rhd, rhop, fsdept_n(:,:,:) ) 113 IF( ln_zps ) CALL zps_hde( kstp, jpts, tsn, gtsu, gtsv, & ! zps: now hor. derivative114 & rhd, gru , grv ) ! of t, s, rd at the last ocean level114 CALL eos ( tsn, rhd, rhop, fsdept_n(:,:,:) ) ! now in situ density for hpg computation 115 IF( ln_zps ) CALL zps_hde( kstp, jpts, tsn, gtsu, gtsv, & ! zps: now hor. derivative 116 & rhd, gru , grv ) ! of t, s, rd at the last ocean level 115 117 116 118 ua(:,:,:) = 0.e0 ! set dynamics trends to zero … … 135 137 CALL div_cur( kstp ) ! Horizontal divergence & Relative vorticity (2nd call in time-split case) 136 138 ENDIF 137 IF( lk_vvl ) CALL dom_vvl_sf_nxt( kstp ) ! after vertical scale factors139 IF( lk_vvl ) CALL dom_vvl_sf_nxt( kstp, kcall=2 ) ! after vertical scale factors 138 140 CALL wzv ( kstp ) ! now cross-level velocity 139 141
Note: See TracChangeset
for help on using the changeset viewer.