Changeset 789 for trunk/NEMO/OPA_SRC/LDF
- Timestamp:
- 2008-01-11T19:04:56+01:00 (16 years ago)
- Location:
- trunk/NEMO/OPA_SRC/LDF
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/LDF/ldfeiv.F90
r719 r789 10 10 !!---------------------------------------------------------------------- 11 11 !! ldf_eiv : compute the eddy induced velocity coefficients 12 !! Same results but not same routine if 'key_mpp_omp'13 !! is defined or not14 12 !!---------------------------------------------------------------------- 15 13 !! * Modules used … … 40 38 41 39 CONTAINS 42 43 # if defined key_mpp_omp44 !!----------------------------------------------------------------------45 !! 'key_mpp_omp' : OpenMP / NEC autotasking (j-slab)46 !!----------------------------------------------------------------------47 48 SUBROUTINE ldf_eiv( kt )49 !!----------------------------------------------------------------------50 !! *** ROUTINE ldf_eiv ***51 !!52 !! ** Purpose : Compute the eddy induced velocity coefficient from the53 !! growth rate of baroclinic instability.54 !!55 !! ** Method :56 !!57 !! ** Action : uslp(), : i- and j-slopes of neutral surfaces58 !! vslp() at u- and v-points, resp.59 !! wslpi(), : i- and j-slopes of neutral surfaces60 !! wslpj() at w-points.61 !!62 !! History :63 !! 8.1 ! 99-03 (G. Madec, A. Jouzeau) Original code64 !! 8.5 ! 02-06 (G. Madec) Free form, F9065 !!----------------------------------------------------------------------66 !! * Arguments67 INTEGER, INTENT( in ) :: kt ! ocean time-step inedx68 69 !! * Local declarations70 INTEGER :: ji, jj, jk ! dummy loop indices71 REAL(wp) :: &72 zfw, ze3w, zn2, zf20, & ! temporary scalars73 zaht, zaht_min74 REAL(wp), DIMENSION(jpi,jpj) :: &75 zn, zah, zhw, zross ! workspace76 !!----------------------------------------------------------------------77 78 IF( kt == nit000 ) THEN79 IF(lwp) WRITE(numout,*)80 IF(lwp) WRITE(numout,*) 'ldf_eiv : eddy induced velocity coefficients'81 IF(lwp) WRITE(numout,*) '~~~~~~~ NEC autotasking / OpenMP : j-slab'82 ENDIF83 84 ! ! ===============85 DO jj = 2, jpjm1 ! Vertical slab86 ! ! ===============87 88 ! 0. Local initialization89 ! -----------------------90 zn (:,jj) = 0.e091 zhw (:,jj) = 5.e092 zah (:,jj) = 0.e093 zross(:,jj) = 0.e094 95 ! 1. Compute lateral diffusive coefficient96 ! ----------------------------------------97 98 !CDIR NOVERRCHK99 DO jk = 1, jpk100 !CDIR NOVERRCHK101 DO ji = 2, jpim1102 ! Take the max of N^2 and zero then take the vertical sum103 ! of the square root of the resulting N^2 ( required to compute104 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f105 zn2 = MAX( rn2(ji,jj,jk), 0.e0 )106 ze3w = fse3w(ji,jj,jk) * tmask(ji,jj,jk)107 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * fse3w(ji,jj,jk)108 ! Compute elements required for the inverse time scale of baroclinic109 ! eddies using the isopycnal slopes calculated in ldfslp.F :110 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w))111 zah(ji,jj) = zah(ji,jj) + zn2 &112 * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk) &113 + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) &114 * ze3w115 zhw(ji,jj) = zhw(ji,jj) + ze3w116 END DO117 END DO118 119 !CDIR NOVERRCHK120 DO ji = 2, jpim1121 zfw = MAX( ABS( 2. * omega * SIN( rad * gphit(ji,jj) ) ) , 1.e-10 )122 ! Rossby radius at w-point taken < 40km and > 2km123 zross(ji,jj) = MAX( MIN( .4 * zn(ji,jj) / zfw, 40.e3 ), 2.e3 )124 ! Compute aeiw by multiplying Ro^2 and T^-1125 aeiw(ji,jj) = zross(ji,jj) * zross(ji,jj) * SQRT( zah(ji,jj) / zhw(ji,jj) ) * tmask(ji,jj,1)126 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA R02127 ! Take the minimum between aeiw and aeiv0 for depth levels128 ! lower than 20 (21 in w- point)129 IF( mbathy(ji,jj) <= 21. ) aeiw(ji,jj) = MIN( aeiw(ji,jj), 1000. )130 ENDIF131 END DO132 133 ! Decrease the coefficient in the tropics (20N-20S)134 zf20 = 2. * omega * sin( rad * 20. )135 DO ji = 2, jpim1136 aeiw(ji,jj) = MIN( 1., ABS( ff(ji,jj) / zf20 ) ) * aeiw(ji,jj)137 END DO138 139 ! ORCA R05: Take the minimum between aeiw and aeiv0140 IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN ! ORCA R05141 DO ji = 2, jpim1142 aeiw(ji,jj) = MIN( aeiw(ji,jj), aeiv0 )143 END DO144 ENDIF145 ! ! ===============146 END DO ! End of slab147 ! ! ===============148 149 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,150 151 ! lateral boundary condition on aeiw152 CALL lbc_lnk( aeiw, 'W', 1. )153 154 ! Average the diffusive coefficient at u- v- points155 DO jj = 2, jpjm1156 DO ji = fs_2, fs_jpim1 ! vector opt.157 aeiu(ji,jj) = .5 * (aeiw(ji,jj) + aeiw(ji+1,jj ))158 aeiv(ji,jj) = .5 * (aeiw(ji,jj) + aeiw(ji ,jj+1))159 END DO160 END DO161 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,162 163 ! lateral boundary condition on aeiu, aeiv164 CALL lbc_lnk( aeiu, 'U', 1. )165 CALL lbc_lnk( aeiv, 'V', 1. )166 167 IF(ln_ctl) THEN168 CALL prt_ctl(tab2d_1=aeiu, clinfo1=' eiv - u: ', ovlap=1)169 CALL prt_ctl(tab2d_1=aeiv, clinfo1=' eiv - v: ', ovlap=1)170 ENDIF171 172 ! ORCA R05: add a space variation on aht (=aeiv except at the equator and river mouth)173 IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN174 zf20 = 2. * omega * SIN( rad * 20. )175 zaht_min = 100. ! minimum value for aht176 DO jj = 1, jpj177 DO ji = 1, jpi178 zaht = ( 1. - MIN( 1., ABS( ff(ji,jj) / zf20 ) ) ) * ( aht0 - zaht_min ) &179 & + aht0 * upsrnfh(ji,jj) ! enhanced near river mouths180 ahtu(ji,jj) = MAX( MAX( zaht_min, aeiu(ji,jj) ) + zaht, aht0 )181 ahtv(ji,jj) = MAX( MAX( zaht_min, aeiv(ji,jj) ) + zaht, aht0 )182 ahtw(ji,jj) = MAX( MAX( zaht_min, aeiw(ji,jj) ) + zaht, aht0 )183 END DO184 END DO185 IF(ln_ctl) THEN186 CALL prt_ctl(tab2d_1=ahtu, clinfo1=' aht - u: ', ovlap=1)187 CALL prt_ctl(tab2d_1=ahtv, clinfo1=' aht - v: ', ovlap=1)188 CALL prt_ctl(tab2d_1=ahtw, clinfo1=' aht - w: ', ovlap=1)189 ENDIF190 ENDIF191 192 IF( aeiv0 == 0.e0 ) THEN193 aeiu(:,:) = 0.e0194 aeiv(:,:) = 0.e0195 aeiw(:,:) = 0.e0196 ENDIF197 198 END SUBROUTINE ldf_eiv199 200 # else201 !!----------------------------------------------------------------------202 !! Default key k-j-i loops203 !!----------------------------------------------------------------------204 40 205 41 SUBROUTINE ldf_eiv( kt ) … … 251 87 252 88 DO jk = 1, jpk 253 # if defined key_vectopt_loop && ! defined key_mpp_omp89 # if defined key_vectopt_loop 254 90 !CDIR NOVERRCHK 255 91 DO ji = 1, jpij ! vector opt. … … 373 209 END SUBROUTINE ldf_eiv 374 210 375 # endif376 377 211 #else 378 212 !!---------------------------------------------------------------------- -
trunk/NEMO/OPA_SRC/LDF/ldfslp.F90
r719 r789 139 139 140 140 IF( ln_zps ) THEN ! partial steps correction at the bottom ocean level (zps_hde routine) 141 # if defined key_vectopt_loop && ! defined key_mpp_omp141 # if defined key_vectopt_loop 142 142 jj = 1 143 143 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) … … 151 151 zgru(ji,jj,iku) = gru(ji,jj) 152 152 zgrv(ji,jj,ikv) = grv(ji,jj) 153 # if ! defined key_vectopt_loop || defined key_mpp_omp153 # if ! defined key_vectopt_loop 154 154 END DO 155 155 # endif … … 490 490 ! mask for mixed layer 491 491 DO jk = 1, jpk 492 # if defined key_vectopt_loop && ! defined key_mpp_omp492 # if defined key_vectopt_loop 493 493 jj = 1 494 494 DO ji = 1, jpij ! vector opt. (forced unrolling) … … 504 504 omlmask(ji,jj,jk) = 0.e0 505 505 ENDIF 506 # if ! defined key_vectopt_loop || defined key_mpp_omp506 # if ! defined key_vectopt_loop 507 507 END DO 508 508 # endif … … 522 522 zwy(:,jpj) = 0.e0 523 523 zwy(jpi,:) = 0.e0 524 # if defined key_vectopt_loop && ! defined key_mpp_omp524 # if defined key_vectopt_loop 525 525 jj = 1 526 526 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) … … 535 535 & * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) ) & 536 536 & / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 537 # if ! defined key_vectopt_loop || defined key_mpp_omp537 # if ! defined key_vectopt_loop 538 538 END DO 539 539 # endif … … 543 543 544 544 ! Slope at u points 545 # if defined key_vectopt_loop && ! defined key_mpp_omp545 # if defined key_vectopt_loop 546 546 jj = 1 547 547 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) … … 560 560 ! uslpml 561 561 uslpml (ji,jj) = zau / ( zbu - zeps ) * umask (ji,jj,ik) 562 # if ! defined key_vectopt_loop || defined key_mpp_omp562 # if ! defined key_vectopt_loop 563 563 END DO 564 564 # endif … … 572 572 zwy ( :, jpj) = 0.e0 573 573 zwy ( jpi, :) = 0.e0 574 # if defined key_vectopt_loop && ! defined key_mpp_omp574 # if defined key_vectopt_loop 575 575 jj = 1 576 576 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) … … 584 584 & * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) ) & 585 585 & / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 586 # if ! defined key_vectopt_loop || defined key_mpp_omp586 # if ! defined key_vectopt_loop 587 587 END DO 588 588 # endif … … 593 593 594 594 ! Slope at v points 595 # if defined key_vectopt_loop && ! defined key_mpp_omp595 # if defined key_vectopt_loop 596 596 jj = 1 597 597 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) … … 610 610 ! vslpml 611 611 vslpml (ji,jj) = zav / ( zbv - zeps ) * vmask (ji,jj,ik) 612 # if ! defined key_vectopt_loop || defined key_mpp_omp612 # if ! defined key_vectopt_loop 613 613 END DO 614 614 # endif … … 624 624 ! Local vertical density gradient evaluated from N^2 625 625 ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 626 # if defined key_vectopt_loop && ! defined key_mpp_omp626 # if defined key_vectopt_loop 627 627 jj = 1 628 628 DO ji = 1, jpij ! vector opt. (forced unrolling) … … 636 636 zwy (ji,jj) = zm05g * pn2 (ji,jj,ik) * & 637 637 & ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 638 # if ! defined key_vectopt_loop || defined key_mpp_omp638 # if ! defined key_vectopt_loop 639 639 END DO 640 640 # endif … … 642 642 643 643 ! Slope at w point 644 # if defined key_vectopt_loop && ! defined key_mpp_omp644 # if defined key_vectopt_loop 645 645 jj = 1 646 646 DO ji = jpi+2, jpij-jpi-1 ! vector opt. (forced unrolling) … … 672 672 wslpiml (ji,jj) = zai / ( zbi - zeps) * tmask (ji,jj,ik) 673 673 wslpjml (ji,jj) = zaj / ( zbj - zeps) * tmask (ji,jj,ik) 674 # if ! defined key_vectopt_loop || defined key_mpp_omp674 # if ! defined key_vectopt_loop 675 675 END DO 676 676 # endif
Note: See TracChangeset
for help on using the changeset viewer.