Changeset 2917 for branches/2011/dev_r2739_LOCEAN8_ZTC
- Timestamp:
- 2011-10-13T18:28:09+02:00 (13 years ago)
- Location:
- branches/2011/dev_r2739_LOCEAN8_ZTC/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2739_LOCEAN8_ZTC/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
r2905 r2917 90 90 END IF 91 91 ! 92 ! - ML - only used in domvvl but could be usefull in many routines92 ! - ML - Used in domvvl and traldf_(lab/bilap/iso) but could be usefull in many other modules 93 93 e12t (:,:) = e1t(:,:) * e2t(:,:) 94 94 e12u (:,:) = e1u(:,:) * e2u(:,:) … … 96 96 e12v_1(:,:) = 0.5 / ( e1v(:,:) * e2v(:,:) ) 97 97 e12f_1(:,:) = 0.5 / ( e1f(:,:) * e2f(:,:) ) 98 ! - ML - used in domvvl and traldf_(lab/bilap/iso)99 98 e1ur (:,:) = e2u(:,:) / e1u(:,:) 100 99 e2vr (:,:) = e1v(:,:) / e2v(:,:) -
branches/2011/dev_r2739_LOCEAN8_ZTC/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r2905 r2917 72 72 !!---------------------------------------------------------------------- 73 73 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 74 ALLOCATE( & 75 & un_td (jpi,jpj,jpk) , vn_td (jpi,jpj,jpk) , hdiv_lf(jpi,jpj,jpk) , & 76 & e3t_t_b(jpi,jpj,jpk) , e3t_t_n(jpi,jpj,jpk) , e3t_t_a(jpi,jpj,jpk) , & 77 & frq_restore_e3t(jpi,jpj) , frq_restore_hdv(jpi,jpj) , STAT= dom_vvl_alloc ) 74 ALLOCATE( e3t_t_n(jpi,jpj,jpk) , e3t_t_a(jpi,jpj,jpk) , e3t_t_b(jpi,jpj,jpk) , & 75 & un_td (jpi,jpj,jpk) , vn_td (jpi,jpj,jpk) , STAT = dom_vvl_alloc ) 76 IF( lk_mpp ) CALL mpp_sum ( dom_vvl_alloc ) 77 IF( dom_vvl_alloc /= 0 ) CALL ctl_warn('dom_vvl_alloc: failed to allocate arrays') 78 ENDIF 79 IF( ln_vvl_ztilde ) THEN 80 ALLOCATE( frq_restore_e3t(jpi,jpj) , frq_restore_hdv(jpi,jpj) , hdiv_lf(jpi,jpj,jpk) , STAT= dom_vvl_alloc ) 78 81 IF( lk_mpp ) CALL mpp_sum ( dom_vvl_alloc ) 79 82 IF( dom_vvl_alloc /= 0 ) CALL ctl_warn('dom_vvl_alloc: failed to allocate arrays') … … 166 169 END DO 167 170 171 ! Restoring frequencies for z_tilde coordinate 172 ! ============================================ 173 IF( ln_vvl_ztilde ) THEN 174 ! - ML - In the future, this should be tunable by the user 175 ! DO jj = 1, jpj 176 ! DO ji = 1, jpi 177 ! frq_restore_hdv(ji,jj) = 2.e0 * rpi / 86400.e0 / 5.e0 & 178 ! & * MAX( SIN( ABS( gphit(ji,jj) ) / .5e0, 1.e0 / 6.e0) ) 179 ! END DO 180 ! END DO 181 ! frq_restore_e3t(:,:) = ( 1.e0 / 6.e0 ) * frq_restore_hdv(:,:) 182 frq_restore_e3t(:,:) = 2.e0 * rpi / ( 30.e0 * 86400.e0 ) 183 frq_restore_hdv(:,:) = 2.e0 * rpi / ( 5.e0 * 86400.e0 ) 184 ! frq_restore_hdv(:,:) = 2.e0 * rpi / ( 2.e0 * 86400.e0 ) 185 ENDIF 186 168 187 END SUBROUTINE dom_vvl_init 169 188 … … 220 239 ! After acale factors at t-points ! 221 240 ! ******************************* ! 222 ! !----------------------------! 223 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! ZTILDE or LAYER coordinate ! 224 ! !----------------------------! 225 226 ! I - Initialisations 227 ! =================== 228 IF( kt == nit000 ) THEN 229 ! - ML - In the future, this should be tunable by the user 230 IF( ln_vvl_ztilde ) THEN ! ZTILDE CASE 231 ! DO jj = 1, jpj 232 ! DO ji = 1, jpi 233 ! frq_restore_hdv(ji,jj) = 2.e0 * rpi / 86400.e0 / 5.e0 & 234 ! & * MAX( SIN( ABS( gphit(ji,jj) ) / .5e0, 1.e0 / 6.e0) ) 235 ! END DO 236 ! END DO 237 ! frq_restore_e3t(:,:) = ( 1.e0 / 6.e0 ) * frq_restore_hdv(:,:) 238 frq_restore_e3t(:,:) = 2.e0 * rpi / ( 30.e0 * 86400.e0 ) 239 frq_restore_hdv(:,:) = 2.e0 * rpi / ( 5.e0 * 86400.e0 ) 240 ! frq_restore_hdv(:,:) = 2.e0 * rpi / ( 2.e0 * 86400.e0 ) 241 ELSE ! LAYER CASE 242 frq_restore_e3t(:,:) = 0.e0 243 frq_restore_hdv(:,:) = 0.e0 244 ENDIF 245 246 ENDIF 247 248 ! II - Low frequency horizontal divergence 249 ! ======================================== 241 242 ! ! ----------------- ! 243 IF( ln_vvl_zstar ) THEN ! z_star coordinate ! 244 ! ! ----------------- ! 245 z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * tmask(:,:,1) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask(:,:,1) ) 246 DO jk = 1, jpkm1 247 fse3t_a(:,:,jk) = fse3t_b(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) 248 END DO 249 250 ! ! --------------------------- ! 251 ELSE ! z_tilde or layer coordinate ! 252 ! ! --------------------------- ! 253 254 ! I - Low frequency horizontal divergence 255 ! ======================================= 250 256 251 257 ! 1 - barotropic divergence … … 259 265 zhdiv(:,:) = zhdiv(:,:) / ( zht(:,:) + 1. - tmask(:,:,1) ) 260 266 261 ! 2 - restoring equation 267 ! 2 - restoring equation (z-tilde case only) 262 268 ! ---------------------- 263 IF( kt .GT. nit000 ) THEN 264 DO jk = 1, jpkm1 265 hdiv_lf(:,:,jk) = hdiv_lf(:,:,jk) - rdt * frq_restore_hdv(:,:) & 266 & * ( hdiv_lf(:,:,jk) - fse3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) ) 267 END DO 268 ENDIF 269 270 ! III - after z_tilde increments of vertical scale factors 271 ! ========================================================= 269 IF( ln_vvl_ztilde ) THEN 270 IF( kt .GT. nit000 ) THEN 271 DO jk = 1, jpkm1 272 hdiv_lf(:,:,jk) = hdiv_lf(:,:,jk) - rdt * frq_restore_hdv(:,:) & 273 & * ( hdiv_lf(:,:,jk) - fse3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) ) 274 END DO 275 ENDIF 276 END IF 277 278 ! II - after z_tilde increments of vertical scale factors 279 ! ======================================================= 272 280 e3t_t_a(:,:,:) = 0.e0 ! e3t_t_a used to store tendency terms 273 281 274 282 ! 1 - High frequency divergence term 275 283 ! ---------------------------------- 276 DO jk = 1, jpkm1 277 e3t_t_a(:,:,jk) = e3t_t_a(:,:,jk) - ( fse3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) - hdiv_lf(:,:,jk) ) 278 END DO 279 280 ! 2 - Restoring term 284 IF( ln_vvl_ztilde ) THEN ! z_tilde case 285 DO jk = 1, jpkm1 286 e3t_t_a(:,:,jk) = e3t_t_a(:,:,jk) - ( fse3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) - hdiv_lf(:,:,jk) ) 287 END DO ! layer case 288 ELSE 289 DO jk = 1, jpkm1 290 e3t_t_a(:,:,jk) = e3t_t_a(:,:,jk) - fse3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) 291 END DO 292 END IF 293 294 ! 2 - Restoring term (z-tilde case only) 281 295 ! ------------------ 282 DO jk = 1, jpk 283 e3t_t_a(:,:,jk) = e3t_t_a(:,:,jk) - frq_restore_e3t(:,:) * e3t_t_b(:,:,jk) 284 END DO 296 IF( ln_vvl_ztilde ) THEN 297 DO jk = 1, jpk 298 e3t_t_a(:,:,jk) = e3t_t_a(:,:,jk) - frq_restore_e3t(:,:) * e3t_t_b(:,:,jk) 299 END DO 300 END IF 285 301 286 302 ! 3 - Thickness diffusion term … … 288 304 zwu(:,:) = 0.e0 289 305 zwv(:,:) = 0.e0 290 291 306 ! a - first derivative: diffusive fluxes 292 307 DO jk = 1, jpkm1 … … 318 333 END DO 319 334 ! d - thickness diffusion equivalent transport: boundary conditions 320 ! (stored for tracer adv action and continuity equation)335 ! (stored for tracer advction and continuity equation) 321 336 CALL lbc_lnk( un_td , 'U' , -1.) 322 337 CALL lbc_lnk( vn_td , 'V' , -1.) … … 331 346 z2dt = 2.e0 * rdt 332 347 ENDIF 333 CALL lbc_lnk( e3t_t_a(:,:,:), 'T', 1. ) ! - ML - Do not use lbc_lnk_e3: e3t_t_a is a tendency term at this stage348 CALL lbc_lnk( e3t_t_a(:,:,:), 'T', 1. ) 334 349 e3t_t_a(:,:,:) = e3t_t_b(:,:,:) + z2dt * tmask(:,:,:) * e3t_t_a(:,:,:) 335 350 fse3t_a(:,:,:) = fse3t_0(:,:,:) + e3t_t_a(:,:,:) … … 360 375 e3t_t_a(:,:,:) = MAX( e3t_t_a(:,:,:), ( 1.e0 - z_def_max ) * fse3t_0(:,:,:) ) 361 376 362 ! Boundary conditions 363 ! ~~~~~~~~~~~~~~~~~~~ 364 ! - ML - think it's not necessary 365 ! CALL lbc_lnk( e3t_t_a(:,:,:), 'T', 1. ) ! - ML - Do not use lbc_lnk_e3: e3t_t_a is a level thickness ANOMALY 366 367 ! IV - Barotropic repartition of the sea surface height over the baroclinic profile 368 ! ================================================================================= 377 ! III - Barotropic repartition of the sea surface height over the baroclinic profile 378 ! ================================================================================== 369 379 ! add e3t(n-1) "star" Asselin-filtered 370 380 DO jk = 1, jpkm1 371 ! - ML - : multiplication by tmask not necessary a priori. Just to be sure for the moment. 372 fse3t_a(:,:,jk) = fse3t_a(:,:,jk) + ( fse3t_b(:,:,jk) - fse3t_0(:,:,jk) - e3t_t_b(:,:,jk) ) * tmask(:,:,jk) 381 fse3t_a(:,:,jk) = fse3t_a(:,:,jk) + fse3t_b(:,:,jk) - fse3t_0(:,:,jk) - e3t_t_b(:,:,jk) 373 382 END DO 374 383 ! add ( ssh increment + "baroclinicity error" ) proportionnaly to e3t(n) … … 382 391 DO jk = 1, jpkm1 383 392 fse3t_a(:,:,jk) = fse3t_a(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 384 END DO385 386 ! !------------------!387 ELSEIF( ln_vvl_zstar ) THEN ! Zstar coordinate !388 ! !------------------!389 z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * tmask(:,:,1) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask(:,:,1) )390 DO jk = 1, jpkm1391 fse3t_a(:,:,jk) = fse3t_b(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:)392 393 END DO 393 394 … … 659 660 id4 = iom_varid( numror, 'e3t_t_n', ldstop = .FALSE. ) 660 661 id5 = iom_varid( numror, 'hdif_lf', ldstop = .FALSE. ) 661 662 ! ! ----------- ! 663 IF( ln_vvl_zstar ) THEN ! z_star case ! 664 ! ! ----------- ! 665 IF( MIN( id1, id2 ) > 0 ) THEN ! all required arrays exist 666 CALL iom_get( numror, jpdom_autoglo, 'fse3t_b', fse3t_b(:,:,:) ) 667 CALL iom_get( numror, jpdom_autoglo, 'fse3t_n', fse3t_n(:,:,:) ) 668 IF( neuler == 0 ) THEN 669 fse3t_b(:,:,:) = fse3t_n(:,:,:) 670 ENDIF 671 ELSE ! one at least array is missing 672 CALL ctl_stop( 'vvl cannot restart from a non vvl run' ) 662 ! ! --------- ! 663 ! ! all cases ! 664 ! ! --------- ! 665 IF( MIN( id1, id2 ) > 0 ) THEN ! all required arrays exist 666 CALL iom_get( numror, jpdom_autoglo, 'fse3t_b', fse3t_b(:,:,:) ) 667 CALL iom_get( numror, jpdom_autoglo, 'fse3t_n', fse3t_n(:,:,:) ) 668 IF( neuler == 0 ) THEN 669 fse3t_b(:,:,:) = fse3t_n(:,:,:) 673 670 ENDIF 674 IF( MIN( id3, id4, id5 ) > 0 ) CALL ctl_stop( 'z_star coordinate cannot restart from a z_tilde run' ) 675 ! ! ------------ ! 676 ELSE ! z_tilde case ! 677 ! ! ------------ ! 678 IF( MIN( id1, id2 ) > 0 ) THEN ! all required arrays exist 679 CALL iom_get( numror, jpdom_autoglo, 'fse3t_b', fse3t_b(:,:,:) ) 680 CALL iom_get( numror, jpdom_autoglo, 'fse3t_n', fse3t_n(:,:,:) ) 681 IF( neuler == 0 ) THEN 682 fse3t_b(:,:,:) = fse3t_n(:,:,:) 683 ENDIF 684 ELSE ! one at least array is missing 685 CALL ctl_stop( 'vvl cannot restart from a non vvl run' ) 671 ELSE ! one at least array is missing 672 CALL ctl_stop( 'dom_vvl_rst: vvl cannot restart from a non vvl run' ) 673 ENDIF 674 ! ! ----------- ! 675 IF( ln_vvl_zstar ) THEN ! z_star case ! 676 ! ! ----------- ! 677 IF( MIN( id3, id4 ) > 0 ) THEN 678 CALL ctl_stop( 'dom_vvl_rst: z_star cannot restart from a z_tilde or layer run' ) 686 679 ENDIF 687 IF( MIN( id3, id4, id5 ) > 0 ) THEN ! all required arrays exist 680 ! ! ----------------------- ! 681 ELSE ! z_tilde and layer cases ! 682 ! ! ----------------------- ! 683 IF( MIN( id3, id4 ) > 0 ) THEN ! all required arrays exist 688 684 CALL iom_get( numror, jpdom_autoglo, 'e3t_t_b', e3t_t_b(:,:,:) ) 689 685 CALL iom_get( numror, jpdom_autoglo, 'e3t_t_n', e3t_t_n(:,:,:) ) 690 CALL iom_get( numror, jpdom_autoglo, 'hdiv_lf', hdiv_lf(:,:,:) ) 691 ELSE ! one at least array is missing 686 ELSE ! one at least array is missing 692 687 e3t_t_b(:,:,:) = 0.e0 693 688 e3t_t_n(:,:,:) = 0.e0 694 hdiv_lf(:,:,:) = 0.e0695 689 ENDIF 696 690 ! ! ------------ ! 691 IF( ln_vvl_ztilde ) THEN ! z_tilde case ! 692 ! ! ------------ ! 693 IF( id5 > 0 ) THEN ! required array exists 694 CALL iom_get( numror, jpdom_autoglo, 'hdiv_lf', hdiv_lf(:,:,:) ) 695 ELSE ! array is missing 696 hdiv_lf(:,:,:) = 0.e0 697 ENDIF 698 ENDIF 697 699 ENDIF 698 700 ! 699 701 ELSE !* Initialize at "rest" 700 702 fse3t_b(:,:,:) = fse3t_0(:,:,:) … … 708 710 ! ! =================== 709 711 IF(lwp) WRITE(numout,*) '---- dom_vvl_rst ----' 710 ! ! ---------------------- !711 ! ! z_star & z_tildecases !712 ! ! ---------------------- !712 ! ! --------- ! 713 ! ! all cases ! 714 ! ! --------- ! 713 715 CALL iom_rstput( kt, nitrst, numrow, 'fse3t_b', fse3t_b(:,:,:) ) 714 716 CALL iom_rstput( kt, nitrst, numrow, 'fse3t_n', fse3t_n(:,:,:) ) 715 ! ! ----------------- !716 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! z_tilde case only!717 ! ! ----------------- !717 ! ! ----------------------- ! 718 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! z_tilde and layer cases ! 719 ! ! ----------------------- ! 718 720 CALL iom_rstput( kt, nitrst, numrow, 'e3t_t_b', e3t_t_b(:,:,:) ) 719 721 CALL iom_rstput( kt, nitrst, numrow, 'e3t_t_n', e3t_t_n(:,:,:) ) 722 END IF 723 ! ! -------------! 724 IF( ln_vvl_ztilde ) THEN ! z_tilde case ! 725 ! ! ------------ ! 720 726 CALL iom_rstput( kt, nitrst, numrow, 'hdiv_lf', hdiv_lf(:,:,:) ) 721 727 ENDIF -
branches/2011/dev_r2739_LOCEAN8_ZTC/NEMOGCM/NEMO/OPA_SRC/step.F90
r2905 r2917 107 107 IF( lk_vvl ) CALL dom_vvl_sf_nxt( kstp ) ! after vertical scale factors 108 108 CALL wzv ( kstp ) ! now cross-level velocity 109 109 110 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 110 111 ! Ocean physics update (ua, va, ta, sa used as workspace)
Note: See TracChangeset
for help on using the changeset viewer.