Changeset 3325 for branches/2012/dev_r3309_LOCEAN12_Ediag
- Timestamp:
- 2012-03-12T15:44:43+01:00 (12 years ago)
- Location:
- branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD
- Files:
-
- 5 edited
- 1 copied
Legend:
- Unmodified
- Added
- Removed
-
branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trddyn.F90
r3318 r3325 5 5 !!===================================================================== 6 6 !! History : 3.5 ! 2012-02 (G. Madec) creation from trdmod: split DYN and TRA trends 7 !! and manage 3D trends output for U, V, and KE7 !! and manage 3D trends output for U, V, and KE 8 8 !!---------------------------------------------------------------------- 9 9 … … 21 21 USE sbc_oce ! surface boundary condition: ocean 22 22 USE phycst ! physical constants 23 USE trdvor ! ocean vorticity trends 24 USE trdglo ! trends:global domain averaged 25 USE trdmld ! ocean active mixed layer tracers trends 23 USE trdken ! trends: Kinetic ENergy 24 USE trdglo ! trends: global domain averaged 25 USE trdvor ! trends: vertical averaged vorticity 26 USE trdmld ! trends: mixed layer averaged 26 27 USE in_out_manager ! I/O manager 27 28 USE iom ! I/O manager library … … 31 32 IMPLICIT NONE 32 33 PRIVATE 33 34 REAL(wp) :: r2dt ! time-step, = 2 rdttra except at nit000 (=rdttra) if neuler=035 34 36 35 PUBLIC trd_dyn ! called by all dynXX modules … … 57 56 INTEGER , INTENT(in ) :: ktrd ! trend index 58 57 INTEGER , INTENT(in ) :: kt ! time step 59 !!60 INTEGER :: ji, jj ! dummy loop indices61 REAL(wp), POINTER, DIMENSION(:,:) :: ztswu, ztswv ! 2D workspace62 58 !!---------------------------------------------------------------------- 59 ! 60 putrd(:,:,:) = putrd(:,:,:) * umask(:,:,:) ! mask the trends 61 pvtrd(:,:,:) = pvtrd(:,:,:) * vmask(:,:,:) 62 ! 63 63 64 CALL wrk_alloc( jpi, jpj, ztswu, ztswv ) 65 66 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt = rdt ! = rdtra (restart with Euler time stepping) 67 ELSEIF( kt <= nit000 + 1) THEN ; r2dt = 2. * rdt ! = 2 rdttra (leapfrog) 68 ENDIF 64 !!gm NB : here a lbc_lnk should probably be added 69 65 70 66 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< … … 81 77 ! Kinetic Energy trends 82 78 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 83 !!gm IF( ln_KE_trd ) CALL trd_KE( putrd, pvtrd, ktrd, kt )79 IF( ln_KE_trd ) CALL trd_ken( putrd, pvtrd, ktrd, kt ) 84 80 85 81 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< … … 92 88 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 93 89 !!gm IF( ln_dyn_mld ) CALL trd_mld_dyn 94 !95 CALL wrk_dealloc( jpi, jpj, ztswu, ztswv )96 90 ! 97 91 END SUBROUTINE trd_dyn … … 113 107 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3dx, z3dy ! 3D workspace 114 108 !!---------------------------------------------------------------------- 115 !116 putrd(:,:,:) = putrd(:,:,:) * umask(:,:,:) ! mask the trends117 pvtrd(:,:,:) = pvtrd(:,:,:) * vmask(:,:,:)118 119 !!gm NB : here a lbc_lnk should probably be added120 121 109 ! 122 110 SELECT CASE( ktrd ) … … 151 139 CALL iom_put( "vtrd_zdf", pvtrd ) 152 140 ! ! wind stress trends 153 z2dx(:,:) = ( utau_b( ji,jj) + utau(ji,jj) ) / ( fse3u(:,:,1) * rau0 )154 z2dy(:,:) = ( vtau_b( ji,jj) + vtau(ji,jj) ) / ( fse3v(:,:,1) * rau0 )141 z2dx(:,:) = ( utau_b(:,:) + utau(:,:) ) / ( fse3u(:,:,1) * rau0 ) 142 z2dy(:,:) = ( vtau_b(:,:) + vtau(:,:) ) / ( fse3v(:,:,1) * rau0 ) 155 143 CALL iom_put( "utrd_tau", z2dx ) 156 144 CALL iom_put( "vtrd_tau", z2dy ) -
branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdglo.F90
r3318 r3325 184 184 185 185 186 SUBROUTINE trd_glo_init187 !!---------------------------------------------------------------------188 !! *** ROUTINE trd_glo_init ***189 !!190 !! ** Purpose : Read the namtrd namelist191 !!----------------------------------------------------------------------192 INTEGER :: ji, jj, jk ! dummy loop indices193 !!----------------------------------------------------------------------194 195 IF(lwp) THEN196 WRITE(numout,*)197 WRITE(numout,*) 'trd_glo_init : integral constraints properties trends'198 WRITE(numout,*) '~~~~~~~~~~~~~'199 ENDIF200 201 ! Total volume at t-points:202 tvolt = 0._wp203 DO jk = 1, jpkm1204 tvolt = SUM( e1e2t(:,:) * fse3t(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:) )205 END DO206 IF( lk_mpp ) CALL mpp_sum( tvolt ) ! sum over the global domain207 208 IF(lwp) WRITE(numout,*) ' total ocean volume at T-point tvolt = ',tvolt209 210 ! Initialization of potential to kinetic energy conversion211 rpktrd = 0._wp212 213 ! Total volume at u-, v- points:214 tvolu = 0._wp215 tvolv = 0._wp216 217 DO jk = 1, jpk218 DO jj = 2, jpjm1219 DO ji = fs_2, fs_jpim1 ! vector opt.220 tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) * tmask_i(ji+1,jj ) * tmask_i(ji,jj) * umask(ji,jj,jk)221 tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) * tmask_i(ji ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)222 END DO223 END DO224 END DO225 IF( lk_mpp ) CALL mpp_sum( tvolu ) ! sums over the global domain226 IF( lk_mpp ) CALL mpp_sum( tvolv )227 228 IF(lwp) THEN229 WRITE(numout,*) ' total ocean volume at U-point tvolu = ',tvolu230 WRITE(numout,*) ' total ocean volume at V-point tvolv = ',tvolv231 ENDIF232 !233 END SUBROUTINE trd_glo_init234 235 236 186 SUBROUTINE glo_dyn_wri( kt ) 237 187 !!--------------------------------------------------------------------- … … 556 506 END SUBROUTINE glo_tra_wri 557 507 508 509 SUBROUTINE trd_glo_init 510 !!--------------------------------------------------------------------- 511 !! *** ROUTINE trd_glo_init *** 512 !! 513 !! ** Purpose : Read the namtrd namelist 514 !!---------------------------------------------------------------------- 515 INTEGER :: ji, jj, jk ! dummy loop indices 516 !!---------------------------------------------------------------------- 517 518 IF(lwp) THEN 519 WRITE(numout,*) 520 WRITE(numout,*) 'trd_glo_init : integral constraints properties trends' 521 WRITE(numout,*) '~~~~~~~~~~~~~' 522 ENDIF 523 524 ! Total volume at t-points: 525 tvolt = 0._wp 526 DO jk = 1, jpkm1 527 tvolt = SUM( e1e2t(:,:) * fse3t(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:) ) 528 END DO 529 IF( lk_mpp ) CALL mpp_sum( tvolt ) ! sum over the global domain 530 531 IF(lwp) WRITE(numout,*) ' total ocean volume at T-point tvolt = ',tvolt 532 533 ! Initialization of potential to kinetic energy conversion 534 rpktrd = 0._wp 535 536 ! Total volume at u-, v- points: 537 !!gm : bug? je suis quasi sur que le produit des tmask_i ne correspond pas exactement au umask_i et vmask_i ! 538 tvolu = 0._wp 539 tvolv = 0._wp 540 541 DO jk = 1, jpk 542 DO jj = 2, jpjm1 543 DO ji = fs_2, fs_jpim1 ! vector opt. 544 tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) * tmask_i(ji+1,jj ) * tmask_i(ji,jj) * umask(ji,jj,jk) 545 tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) * tmask_i(ji ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk) 546 END DO 547 END DO 548 END DO 549 IF( lk_mpp ) CALL mpp_sum( tvolu ) ! sums over the global domain 550 IF( lk_mpp ) CALL mpp_sum( tvolv ) 551 552 IF(lwp) THEN 553 WRITE(numout,*) ' total ocean volume at U-point tvolu = ',tvolu 554 WRITE(numout,*) ' total ocean volume at V-point tvolv = ',tvolv 555 ENDIF 556 ! 557 END SUBROUTINE trd_glo_init 558 558 559 !!====================================================================== 559 560 END MODULE trdglo -
branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdini.F90
r3318 r3325 11 11 !!---------------------------------------------------------------------- 12 12 USE trd_oce ! trends: ocean variables 13 ! USE ldftra_oce ! ocean active tracers lateral physics 14 USE trdglo ! ocean bassin integral constraints properties13 USE trdken ! trends: 3D kinetic energy 14 USE trdglo ! trends: global domain integral constraints properties 15 15 USE trdmld ! ocean active mixed layer tracers trends 16 16 USE trdvor ! ocean vorticity trends … … 79 79 ! 80 80 81 IF( ln_ KE_trd .OR. ln_PE_trd .OR. ln_dyn_mld ) &82 CALL ctl_stop( ' KE, PE, aur ML on momentum are not yet coded we stop' )81 IF( ln_PE_trd .OR. ln_dyn_mld ) & 82 CALL ctl_stop( 'PE, or ML on momentum are not yet coded we stop' ) 83 83 84 84 ! 85 IF( ln_glo_trd ) CALL trd_glo_init ! integral constraints trends 86 IF( ln_tra_mld ) CALL trd_mld_init ! mixed-layer trends (active tracers) 87 IF( ln_vor_trd ) CALL trd_vor_init ! vorticity trends 85 IF( ln_glo_trd ) CALL trd_glo_init ! integral constraints trends 86 IF( ln_tra_mld ) CALL trd_mld_init ! mixed-layer trends (active tracers) 87 IF( ln_vor_trd ) CALL trd_vor_init ! vorticity trends 88 IF( ln_KE_trd ) CALL trd_ken_init ! 3D Kinetic energy diagnostics 89 88 90 ! 89 91 !!gm : Potential BUG : 3D output only for vector invariant form! add a ctl_stop or code the flux form case -
branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90
r3318 r3325 1 MODULE trd dyn1 MODULE trdken 2 2 !!====================================================================== 3 !! *** MODULE trd dyn ***4 !! Ocean diagnostics: ocean dynamictrends3 !! *** MODULE trdken *** 4 !! Ocean diagnostics: compute and output 3D kinetic energy trends 5 5 !!===================================================================== 6 !! History : 3.5 ! 2012-02 (G. Madec) creation from trdmod: split DYN and TRA trends 7 !! and manage 3D trends output for U, V, and KE 8 !!---------------------------------------------------------------------- 9 10 !!---------------------------------------------------------------------- 11 !! trd_dyn : manage the type of momentum trend diagnostics (3D I/O, domain averaged, KE) 12 !! trd_dyn_iom : output 3D momentum and/or tracer trends using IOM 13 !! trd_dyn_init : initialization step 6 !! History : 3.5 ! 2012-02 (G. Madec) original code 7 !!---------------------------------------------------------------------- 8 9 !!---------------------------------------------------------------------- 10 !! trd_ken : compute and output 3D Kinetic energy trends using IOM 11 !! trd_ken_init : initialisation 14 12 !!---------------------------------------------------------------------- 15 13 USE oce ! ocean dynamics and tracers variables … … 32 30 PRIVATE 33 31 34 REAL(wp) :: r2dt ! time-step, = 2 rdttra except at nit000 (=rdttra) if neuler=0 35 36 PUBLIC trd_dyn ! called by all dynXX modules 32 PUBLIC trd_ken ! called by trddyn module 33 PUBLIC trd_ken_init ! called by trdini module 34 35 INTEGER :: nkstp ! current time step 36 REAL(wp):: r1_2_rau0 ! = 0.5 / rau0 37 38 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: bu, bv ! volume of u- and v-boxes 39 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: r1_bt ! inverse of t-box volume 37 40 38 41 !! * Substitutions … … 46 49 CONTAINS 47 50 48 SUBROUTINE trd_dyn( putrd, pvtrd, ktrd, kt ) 49 !!--------------------------------------------------------------------- 50 !! *** ROUTINE trd_mod *** 51 INTEGER FUNCTION trd_ken_alloc() 52 !!--------------------------------------------------------------------- 53 !! *** FUNCTION trd_ken_alloc *** 54 !!--------------------------------------------------------------------- 55 ALLOCATE( bu(jpi,jpj,jpk) , bv(jpi,jpj,jpk) , r1_bt(jpi,jpj,jpk) , STAT= trd_ken_alloc ) 56 ! 57 IF( lk_mpp ) CALL mpp_sum ( trd_ken_alloc ) 58 IF( trd_ken_alloc /= 0 ) CALL ctl_warn('trd_ken_alloc: failed to allocate arrays') 59 END FUNCTION trd_ken_alloc 60 61 62 SUBROUTINE trd_ken( putrd, pvtrd, ktrd, kt ) 63 !!--------------------------------------------------------------------- 64 !! *** ROUTINE trd_ken *** 51 65 !! 52 !! ** Purpose : Dispatch momentum trend computation, e.g. 3D output, 53 !! integral constraints, barotropic vorticity, kinetic enrgy, 54 !! and/or mixed layer budget. 55 !!---------------------------------------------------------------------- 56 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: putrd, pvtrd ! U and V trends 66 !! ** Purpose : output 3D Kinetic Energy trends using IOM 67 !! 68 !! ** Method : - apply lbc to the input masked velocity trends 69 !! - compute the associated KE trend: 70 !! zke = 0.5 * ( mi-1[ un * putrd * bu ] + mj-1[ vn * pvtrd * bv] ) / bt 71 !! where bu, bv, bt are the volume of u-, v- and t-boxes. 72 !! - vertical diffusion case (jpdyn_zdf): 73 !! diagnose separately the KE trend associated with wind stress 74 !! - bottom friction case (jpdyn_bfr): 75 !! explicit case (ln_bfrimp=F): bottom trend put in the 1st level 76 !! of putrd, pvtrd 77 ! 78 ! 79 !!---------------------------------------------------------------------- 80 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: putrd, pvtrd ! U and V masked trends 57 81 INTEGER , INTENT(in ) :: ktrd ! trend index 58 82 INTEGER , INTENT(in ) :: kt ! time step 59 !! 60 INTEGER :: ji, jj ! dummy loop indices 61 REAL(wp), POINTER, DIMENSION(:,:) :: ztswu, ztswv ! 2D workspace 62 !!---------------------------------------------------------------------- 63 64 CALL wrk_alloc( jpi, jpj, ztswu, ztswv ) 65 66 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt = rdt ! = rdtra (restart with Euler time stepping) 67 ELSEIF( kt <= nit000 + 1) THEN ; r2dt = 2. * rdt ! = 2 rdttra (leapfrog) 83 ! 84 INTEGER :: ji, jj, jk ! dummy loop indices 85 INTEGER :: ikbu , ikbv ! local integers 86 INTEGER :: ikbum1, ikbvm1 ! - - 87 REAL(wp), POINTER, DIMENSION(:,:) :: z2dx, z2dy, zke2d ! 2D workspace 88 REAL(wp), POINTER, DIMENSION(:,:,:) :: zke ! 3D workspace 89 !!---------------------------------------------------------------------- 90 ! 91 CALL lbc_lnk( putrd, 'U', -1. ) ; CALL lbc_lnk( pvtrd, 'V', -1. ) ! lateral boundary conditions 92 ! 93 IF ( lk_vvl .AND. kt /= nkstp ) THEN ! Variable volume: set box volume at the 1st call of kt time step 94 nkstp = kt 95 DO jk = 1, jpkm1 96 bu (:,:,jk) = e1u(:,:) * e2u(:,:) * fse3u_n(:,:,jk) 97 bv (:,:,jk) = e1v(:,:) * e2v(:,:) * fse3v_n(:,:,jk) 98 r1_bt(:,:,jk) = 1._wp / ( e1e2t(:,:) * fse3t_n(:,:,jk) ) * tmask(:,:,jk) 99 END DO 68 100 ENDIF 69 70 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 71 ! 3D output of momentum and/or tracers trends using IOM interface 72 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 73 IF( ln_dyn_trd ) CALL trd_dyn_iom( putrd, pvtrd, ktrd, kt ) 74 75 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 76 ! Integral Constraints Properties for momentum and/or tracers trends 77 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 78 IF( ln_glo_trd ) CALL trd_glo( putrd, pvtrd, ktrd, 'DYN', kt ) 79 80 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 81 ! Kinetic Energy trends 82 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 83 !!gm IF( ln_KE_trd ) CALL trd_KE( putrd, pvtrd, ktrd, kt ) 84 85 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 86 ! Vorticity trends 87 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 88 IF( ln_vor_trd ) CALL trd_vor( putrd, pvtrd, ktrd, kt ) 89 90 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 91 ! Mixed layer trends for active tracers 92 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 93 !!gm IF( ln_dyn_mld ) CALL trd_mld_dyn 94 ! 95 CALL wrk_dealloc( jpi, jpj, ztswu, ztswv ) 96 ! 97 END SUBROUTINE trd_dyn 98 99 100 SUBROUTINE trd_dyn_iom( putrd, pvtrd, ktrd, kt ) 101 !!--------------------------------------------------------------------- 102 !! *** ROUTINE trd_dyn_iom *** 103 !! 104 !! ** Purpose : output 3D trends using IOM 105 !!---------------------------------------------------------------------- 106 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: putrd, pvtrd ! U and V trends 107 INTEGER , INTENT(in ) :: ktrd ! trend index 108 INTEGER , INTENT(in ) :: kt ! time step 109 ! 110 INTEGER :: ji, jj, jk ! dummy loop indices 111 INTEGER :: ikbu, ikbv ! local integers 112 REAL(wp), POINTER, DIMENSION(:,:) :: z2dx, z2dy, ztswu, ztswv ! 2D workspace 113 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3dx, z3dy ! 3D workspace 114 !!---------------------------------------------------------------------- 115 ! 116 putrd(:,:,:) = putrd(:,:,:) * umask(:,:,:) ! mask the trends 117 pvtrd(:,:,:) = pvtrd(:,:,:) * vmask(:,:,:) 118 119 !!gm NB : here a lbc_lnk should probably be added 120 121 ! 122 SELECT CASE( ktrd ) 123 CASE( jpdyn_hpg ) ; CALL iom_put( "utrd_hpg", putrd ) ! hydrostatic pressure gradient 124 CALL iom_put( "vtrd_hpg", pvtrd ) 125 CASE( jpdyn_spg ) ; CALL iom_put( "utrd_spg", putrd ) ! surface pressure gradient 126 CALL iom_put( "vtrd_spg", pvtrd ) 127 CASE( jpdyn_pvo ) ; CALL iom_put( "utrd_pvo", putrd ) ! planetary vorticity 128 CALL iom_put( "vtrd_pvo", pvtrd ) 129 CASE( jpdyn_rvo ) ; CALL iom_put( "utrd_rvo", putrd ) ! relative vorticity (or metric term) 130 CALL iom_put( "vtrd_rvo", pvtrd ) 131 CASE( jpdyn_keg ) ; CALL iom_put( "utrd_keg", putrd ) ! Kinetic Energy gradient (or had) 132 CALL iom_put( "vtrd_keg", pvtrd ) 133 z3dx(:,:,:) = 0._wp ! U.dxU & V.dyV (approximation) 134 z3dy(:,:,:) = 0._wp 135 DO jk = 1, jpkm1 ! no mask as un,vn are masked 136 DO jj = 2, jpjm1 137 DO ji = 2, jpim1 138 z3dx(ji,jj,jk) = un(ji,jj,jk) * ( un(ji+1,jj,jk) - un(ji-1,jj,jk) ) / ( 2._wp * e1u(ji,jj) ) 139 z3dy(ji,jj,jk) = vn(ji,jj,jk) * ( vn(ji,jj+1,jk) - vn(ji,jj-1,jk) ) / ( 2._wp * e2v(ji,jj) ) 140 END DO 101 ! 102 zke(:,:,jpk) = 0._wp 103 zke(1,:, : ) = 0._wp 104 zke(:,1, : ) = 0._wp 105 DO jk = 1, jpkm1 106 DO jj = 2, jpj 107 DO ji = 2, jpj 108 zke(ji,jj,jk) = 0.5_wp * ( un(ji ,jj,jk) * putrd(ji ,jj,jk) * bu(ji ,jj,jk) & 109 & + un(ji-1,jj,jk) * putrd(ji-1,jj,jk) * bu(ji-1,jj,jk) & 110 & + vn(ji,jj ,jk) * pvtrd(ji,jj ,jk) * bv(ji,jj ,jk) & 111 & + vn(ji,jj-1,jk) * pvtrd(ji,jj-1,jk) * bv(ji,jj-1,jk) ) * r1_bt(ji,jj,jk) 141 112 END DO 142 113 END DO 143 CALL lbc_lnk( z3dx, 'U', -1. ) ; CALL lbc_lnk( z3dy, 'V', -1. ) 144 CALL iom_put( "utrd_udx", z3dx ) 145 CALL iom_put( "vtrd_vdy", z3dy ) 146 CASE( jpdyn_zad ) ; CALL iom_put( "utrd_zad", putrd ) ! vertical advection 147 CALL iom_put( "vtrd_zad", pvtrd ) 148 CASE( jpdyn_ldf ) ; CALL iom_put( "utrd_ldf", putrd ) ! lateral diffusion 149 CALL iom_put( "vtrd_ldf", pvtrd ) 150 CASE( jpdyn_zdf ) ; CALL iom_put( "utrd_zdf", putrd ) ! vertical diffusion 151 CALL iom_put( "vtrd_zdf", pvtrd ) 152 ! ! wind stress trends 153 z2dx(:,:) = ( utau_b(ji,jj) + utau(ji,jj) ) / ( fse3u(:,:,1) * rau0 ) 154 z2dy(:,:) = ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( fse3v(:,:,1) * rau0 ) 155 CALL iom_put( "utrd_tau", z2dx ) 156 CALL iom_put( "vtrd_tau", z2dy ) 157 CASE( jpdyn_bfr ) 158 IF( .NOT.ln_bfrimp ) CALL iom_put( "utrd_bfr", putrd ) ! bottom friction (explicit case) 159 IF( .NOT.ln_bfrimp ) CALL iom_put( "vtrd_bfr", pvtrd ) 114 END DO 115 ! 116 SELECT CASE( ktrd ) 117 CASE( jpdyn_hpg ) ; CALL iom_put( "ketrd_hpg", zke ) ! hydrostatic pressure gradient 118 CASE( jpdyn_spg ) ; CALL iom_put( "ketrd_spg", zke ) ! surface pressure gradient 119 CASE( jpdyn_pvo ) ; CALL iom_put( "ketrd_pvo", zke ) ! planetary vorticity 120 CASE( jpdyn_rvo ) ; CALL iom_put( "ketrd_rvo", zke ) ! relative vorticity (or metric term) 121 CASE( jpdyn_keg ) ; CALL iom_put( "ketrd_keg", zke ) ! Kinetic Energy gradient (or had) 122 CASE( jpdyn_zad ) ; CALL iom_put( "ketrd_zad", zke ) ! vertical advection 123 CASE( jpdyn_ldf ) ; CALL iom_put( "ketrd_ldf", zke ) ! lateral diffusion 124 CASE( jpdyn_zdf ) ; CALL iom_put( "ketrd_zdf", zke ) ! vertical diffusion 125 ! ! wind stress trends 126 z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1u(:,:) * e2u(:,:) * umask(:,:,1) 127 z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1v(:,:) * e2v(:,:) * vmask(:,:,1) 128 zke2d(1,:) = 0._wp ; zke2d(:,1) = 0._wp 129 DO jj = 2, jpj 130 DO ji = 2, jpj 131 zke2d(ji,jj) = r1_2_rau0 * ( z2dx(ji,jj) + z2dx(ji-1,jj) & 132 & + z2dy(ji,jj) + z2dy(ji,jj-1) ) * r1_bt(ji,jj,1) 133 END DO ! 134 END DO ! 135 CALL iom_put( "ketrd_tau", zke2d ) 136 CASE( jpdyn_bfr ) ; CALL iom_put( "ketrd_bfr", zke ) ! bottom friction (explicit) 137 !!gm TO BE DONE properly 138 ! IF( .NOT.ln_bfrimp ) THEN 139 ! DO jj = 1, jpj ! 140 ! DO ji = 1, jpi 141 ! ikbu = mbku(ji,jj) ! deepest ocean u- & v-levels 142 ! ikbv = mbkv(ji,jj) 143 ! z2dx(ji,jj) = un(ji,jj,ikbu) * bfrua(ji,jj) * un(ji,jj,ikbu) 144 ! z2dy(ji,jj) = vn(ji,jj,ikbu) * bfrva(ji,jj) * vn(ji,jj,ikbv) 145 ! END DO 146 ! END DO 147 ! zke2d(1,:) = 0._wp ; zke2d(:,1) = 0._wp 148 ! DO jj = 2, jpj 149 ! DO ji = 2, jpj 150 ! zke2d(ji,jj) = 0.5_wp * ( z2dx(ji,jj) + z2dx(ji-1,jj) & 151 ! & + z2dy(ji,jj) + z2dy(ji,jj-1) ) * r1_bt(ji,jj, BEURK!!! 152 ! END DO 153 ! END DO 154 ! CALL iom_put( "ketrd_bfr", zke2d ) ! bottom friction (explicit case) 160 155 !!gm only valid if ln_bfrimp=T otherwise the bottom stress as to be recomputed at the end of the compuation.... 161 162 CASE( jpdyn_atf ) ; CALL iom_put( "utrd_atf", putrd ) ! asselin filter trends 163 CALL iom_put( "vtrd_atf", pvtrd ) 164 IF( ln_bfrimp ) THEN ! bottom friction (implicit case) 165 z3dx(:,:,:) = 0._wp ; z3dy(:,:,:) = 0._wp ! after velocity known (now filed at this stage) 166 DO jk = 1, jpkm1 167 DO jj = 2, jpjm1 168 DO ji = 2, jpim1 169 ikbu = mbku(ji,jj) ! deepest ocean u- & v-levels 170 ikbv = mbkv(ji,jj) 171 z3dx(ji,jj,jk) = bfrua(ji,jj) * un(ji,jj,ikbu) / fse3u(ji,jj,ikbu) 172 z3dy(ji,jj,jk) = bfrva(ji,jj) * vn(ji,jj,ikbv) / fse3v(ji,jj,ikbv) 173 END DO 174 END DO 175 END DO 176 CALL iom_put( "utrd_bfr", z3dx ) ! bottom friction (implicit) 177 CALL iom_put( "vtrd_bfr", z3dy ) 178 ENDIF 156 ! ENDIF 157 !!gm end 158 CASE( jpdyn_atf ) ; CALL iom_put( "ketrd_atf", zke ) ! asselin filter trends 159 160 161 !! a faire !!!! idee changer dynnxt pour avoir un appel a jpdynbfr avant le swap !!! 162 !! reflechir à une possible sauvegarde du "vrai" un,vn pour le calcul de atf.... 163 ! 164 ! IF( ln_bfrimp ) THEN ! bottom friction (implicit case) 165 ! DO jj = 1, jpj ! after velocity known (now filed at this stage) 166 ! DO ji = 1, jpi 167 ! ikbu = mbku(ji,jj) ! deepest ocean u- & v-levels 168 ! ikbv = mbkv(ji,jj) 169 ! z2dx(ji,jj) = un(ji,jj,ikbu) * bfrua(ji,jj) * un(ji,jj,ikbu) / fse3u(ji,jj,ikbu) 170 ! z2dy(ji,jj) = un(ji,jj,ikbu) * bfrva(ji,jj) * vn(ji,jj,ikbv) / fse3v(ji,jj,ikbv) 171 ! END DO 172 ! END DO 173 ! zke2d(1,:) = 0._wp ; zke2d(:,1) = 0._wp 174 ! DO jj = 2, jpj 175 ! DO ji = 2, jpj 176 ! zke2d(ji,jj) = 0.5_wp * ( z2dx(ji,jj) + z2dx(ji-1,jj) & 177 ! & + z2dy(ji,jj) + z2dy(ji,jj-1) ) 178 ! END DO 179 ! END DO 180 ! CALL iom_put( "ketrd_bfr", zke2d ) 181 ! ENDIF 179 182 ! 180 183 END SELECT 181 184 ! 182 CALL wrk_dealloc( jpi, jpj , z2dx, z2dy, ztswu, ztswv ) 183 CALL wrk_dealloc( jpi, jpj, jpk, z3dx, z3dy ) 184 ! 185 END SUBROUTINE trd_dyn_iom 185 CALL wrk_dealloc( jpi, jpj , z2dx, z2dy, zke2d ) 186 CALL wrk_dealloc( jpi, jpj, jpk, zke ) 187 ! 188 END SUBROUTINE trd_ken 189 190 191 SUBROUTINE trd_ken_init 192 !!--------------------------------------------------------------------- 193 !! *** ROUTINE trd_ken_init *** 194 !! 195 !! ** Purpose : initialisation of 3D Kinetic Energy trend diagnostic 196 !!---------------------------------------------------------------------- 197 INTEGER :: ji, jj, jk ! dummy loop indices 198 !!---------------------------------------------------------------------- 199 ! 200 IF(lwp) THEN 201 WRITE(numout,*) 202 WRITE(numout,*) 'trd_ken_init : 3D Kinetic Energy trends' 203 WRITE(numout,*) '~~~~~~~~~~~~~' 204 ENDIF 205 ! ! allocate box volume arrays 206 IF ( trd_ken_alloc() /= 0 ) CALL ctl_stop('trd_ken_alloc: failed to allocate arrays') 207 ! 208 IF ( .NOT.lk_vvl ) THEN ! constant volume: bu, bv, 1/bt computed one for all 209 DO jk = 1, jpkm1 210 bu (:,:,jk) = e1u(:,:) * e2u(:,:) * fse3u_n(:,:,jk) 211 bv (:,:,jk) = e1v(:,:) * e2v(:,:) * fse3v_n(:,:,jk) 212 r1_bt(:,:,jk) = 1._wp / ( e1e2t(:,:) * fse3t_n(:,:,jk) ) 213 END DO 214 ENDIF 215 ! 216 nkstp = 0 217 r1_2_rau0 = 0.5_wp / rau0 218 ! 219 END SUBROUTINE trd_ken_init 186 220 187 221 !!====================================================================== 188 END MODULE trd dyn222 END MODULE trdken -
branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdvor.F90
r3318 r3325 116 116 CALL trd_vor_zint( putrd, pvtrd, jpvor_zdf ) ! zdf trend including surf./bot. stresses 117 117 CALL trd_vor_zint( ztswu, ztswv, jpvor_swf ) ! surface wind stress 118 CASE 118 CASE( jpdyn_bfr ) 119 119 CALL trd_vor_zint( putrd, pvtrd, jpvor_bfr ) ! Bottom stress 120 120 ! … … 176 176 ! ===================================== 177 177 178 SELECT CASE (ktrd)179 ! 180 CASE (jpvor_bfr) ! bottom friction178 SELECT CASE( ktrd ) 179 ! 180 CASE( jpvor_bfr ) ! bottom friction 181 181 DO jj = 2, jpjm1 182 182 DO ji = fs_2, fs_jpim1 … … 188 188 END DO 189 189 ! 190 CASE (jpvor_swf) ! wind stress190 CASE( jpvor_swf ) ! wind stress 191 191 zudpvor(:,:) = putrdvor(:,:) * fse3u(:,:,1) * e1u(:,:) * umask(:,:,1) 192 192 zvdpvor(:,:) = pvtrdvor(:,:) * fse3v(:,:,1) * e2v(:,:) * vmask(:,:,1) … … 199 199 200 200 ! Curl 201 DO ji =1,jpim1202 DO jj =1,jpjm1201 DO ji = 1, jpim1 202 DO jj = 1, jpjm1 203 203 vortrd(ji,jj,ktrd) = ( zvdpvor(ji+1,jj) - zvdpvor(ji,jj) & 204 204 & - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) -
branches/2012/dev_r3309_LOCEAN12_Ediag/NEMOGCM/NEMO/OPA_SRC/TRD/trdvor_oce.F90
r3318 r3325 4 4 !! Ocean trends : set vorticity trend variables 5 5 !!====================================================================== 6 !! History : 9.0 ! 04-2006 (L. Brunier, A-M. Treguier) Original code6 !! History : 1.0 ! 04-2006 (L. Brunier, A-M. Treguier) Original code 7 7 !!---------------------------------------------------------------------- 8 9 !!---------------------------------------------------------------------- 8 10 9 USE par_oce ! ocean parameters 11 10
Note: See TracChangeset
for help on using the changeset viewer.