- Timestamp:
- 2015-07-10T13:28:53+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/OPA_SRC/TRA/tranpc.F90
r4313 r5581 2 2 !!============================================================================== 3 3 !! *** MODULE tranpc *** 4 !! Ocean active tracers: non penetrative convecti onscheme4 !! Ocean active tracers: non penetrative convective adjustment scheme 5 5 !!============================================================================== 6 6 !! History : 1.0 ! 1990-09 (G. Madec) Original code … … 9 9 !! 3.0 ! 2008-06 (G. Madec) applied on ta, sa and called before tranxt in step.F90 10 10 !! 3.3 ! 2010-05 (C. Ethe, G. Madec) merge TRC-TRA 11 !! 3.6 ! 2015-05 (L. Brodeau) new algorithm based on local Brunt-Vaisala freq. 11 12 !!---------------------------------------------------------------------- 12 13 … … 14 15 !! tra_npc : apply the non penetrative convection scheme 15 16 !!---------------------------------------------------------------------- 16 USE oce ! ocean dynamics and active tracers 17 USE oce ! ocean dynamics and active tracers 17 18 USE dom_oce ! ocean space and time domain 19 USE phycst ! physical constants 18 20 USE zdf_oce ! ocean vertical physics 19 USE trd mod_oce! ocean active tracer trends21 USE trd_oce ! ocean active tracer trends 20 22 USE trdtra ! ocean active tracer trends 21 USE eosbn2 ! equation of state (eos routine) 23 USE eosbn2 ! equation of state (eos routine) 24 ! 22 25 USE lbclnk ! lateral boundary conditions (or mpp link) 23 26 USE in_out_manager ! I/O manager … … 29 32 PRIVATE 30 33 31 PUBLIC tra_npc 34 PUBLIC tra_npc ! routine called by step.F90 32 35 33 36 !! * Substitutions 34 37 # include "domzgr_substitute.h90" 35 !!---------------------------------------------------------------------- 36 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 37 !! $Id$ 38 # include "vectopt_loop_substitute.h90" 39 !!---------------------------------------------------------------------- 40 !! NEMO/OPA 3.6 , NEMO Consortium (2014) 41 !! $Id$ 38 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 39 43 !!---------------------------------------------------------------------- … … 44 48 !! *** ROUTINE tranpc *** 45 49 !! 46 !! ** Purpose : Non penetrative convective adjustment scheme. solve50 !! ** Purpose : Non-penetrative convective adjustment scheme. solve 47 51 !! the static instability of the water column on after fields 48 52 !! while conserving heat and salt contents. 49 53 !! 50 !! ** Method : The algorithm used converges in a maximium of jpk 51 !! iterations. instabilities are treated when the vertical density 52 !! gradient is less than 1.e-5. 53 !! l_trdtra=T: the trend associated with this algorithm is saved. 54 !! ** Method : updated algorithm able to deal with non-linear equation of state 55 !! (i.e. static stability computed locally) 54 56 !! 55 57 !! ** Action : - (ta,sa) after the application od the npc scheme 56 !! - s ave the associated trends (ttrd,strd) ('key_trdtra')58 !! - send the associated trends for on-line diagnostics (l_trdtra=T) 57 59 !! 58 !! References : Madec, et al., 1991, JPO, 21, 9, 1349-1371.60 !! References : Madec, et al., 1991, JPO, 21, 9, 1349-1371. 59 61 !!---------------------------------------------------------------------- 60 !61 62 INTEGER, INTENT(in) :: kt ! ocean time-step index 62 63 ! 63 64 INTEGER :: ji, jj, jk ! dummy loop indices 64 65 INTEGER :: inpcc ! number of statically instable water column 65 INTEGER :: inpci ! number of iteration for npc scheme 66 INTEGER :: jiter, jkdown, jkp ! ??? 67 INTEGER :: ikbot, ik, ikup, ikdown ! ??? 68 REAL(wp) :: ze3tot, zta, zsa, zraua, ze3dwn 69 REAL(wp), POINTER, DIMENSION(:,: ) :: zwx, zwy, zwz 70 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdt, ztrds, zrhop 66 INTEGER :: jiter, ikbot, ikp, ikup, ikdown, ilayer, ik_low ! local integers 67 LOGICAL :: l_bottom_reached, l_column_treated 68 REAL(wp) :: zta, zalfa, zsum_temp, zsum_alfa, zaw, zdz, zsum_z 69 REAL(wp) :: zsa, zbeta, zsum_sali, zsum_beta, zbw, zrw, z1_r2dt 70 REAL(wp), PARAMETER :: zn2_zero = 1.e-14_wp ! acceptance criteria for neutrality (N2==0) 71 REAL(wp), POINTER, DIMENSION(:) :: zvn2 ! vertical profile of N2 at 1 given point... 72 REAL(wp), POINTER, DIMENSION(:,:) :: zvts ! vertical profile of T and S at 1 given point... 73 REAL(wp), POINTER, DIMENSION(:,:) :: zvab ! vertical profile of alpha and beta 74 REAL(wp), POINTER, DIMENSION(:,:,:) :: zn2 ! N^2 75 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: zab ! alpha and beta 76 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdt, ztrds ! 3D workspace 77 ! 78 LOGICAL, PARAMETER :: l_LB_debug = .FALSE. ! set to true if you want to follow what is 79 INTEGER :: ilc1, jlc1, klc1, nncpu ! actually happening in a water column at point "ilc1, jlc1" 80 LOGICAL :: lp_monitor_point = .FALSE. ! in CPU domain "nncpu" 71 81 !!---------------------------------------------------------------------- 72 82 ! 73 83 IF( nn_timing == 1 ) CALL timing_start('tra_npc') 74 84 ! 75 CALL wrk_alloc(jpi, jpj, jpk, zrhop )76 CALL wrk_alloc(jpi, jpk, zwx, zwy, zwz )77 !78 85 IF( MOD( kt, nn_npc ) == 0 ) THEN 79 80 inpcc = 081 inpci = 082 83 CALL eos( tsa, rhd, zrhop, fsdept_n(:,:,:) ) ! Potential density84 85 IF( l_trdtra ) THEN !* Save ta and sa trends86 ! 87 CALL wrk_alloc( jpi, jpj, jpk, zn2 ) ! N2 88 CALL wrk_alloc( jpi, jpj, jpk, 2, zab ) ! Alpha and Beta 89 CALL wrk_alloc( jpk, 2, zvts, zvab ) ! 1D column vector at point ji,jj 90 CALL wrk_alloc( jpk, zvn2 ) ! 1D column vector at point ji,jj 91 92 IF( l_trdtra ) THEN !* Save initial after fields 86 93 CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds ) 87 94 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) … … 89 96 ENDIF 90 97 91 ! ! =============== 92 DO jj = 1, jpj ! Vertical slab 93 ! ! =============== 94 ! Static instability pointer 95 ! ---------------------------- 96 DO jk = 1, jpkm1 97 DO ji = 1, jpi 98 zwx(ji,jk) = ( zrhop(ji,jj,jk) - zrhop(ji,jj,jk+1) ) * tmask(ji,jj,jk+1) 99 END DO 100 END DO 101 102 ! 1.1 do not consider the boundary points 103 104 ! even if east-west cyclic b. c. do not considere ji=1 or jpi 105 DO jk = 1, jpkm1 106 zwx( 1 ,jk) = 0.e0 107 zwx(jpi,jk) = 0.e0 108 END DO 109 ! even if south-symmetric b. c. used, do not considere jj=1 110 IF( jj == 1 ) zwx(:,:) = 0.e0 111 112 DO jk = 1, jpkm1 113 DO ji = 1, jpi 114 zwx(ji,jk) = 1. 115 IF( zwx(ji,jk) < 1.e-5 ) zwx(ji,jk) = 0.e0 116 END DO 117 END DO 118 119 zwy(:,1) = 0.e0 120 DO ji = 1, jpi 121 DO jk = 1, jpkm1 122 zwy(ji,1) = zwy(ji,1) + zwx(ji,jk) 123 END DO 124 END DO 125 126 zwz(1,1) = 0.e0 127 DO ji = 1, jpi 128 zwz(1,1) = zwz(1,1) + zwy(ji,1) 129 END DO 130 131 inpcc = inpcc + NINT( zwz(1,1) ) 132 133 134 ! 2. Vertical mixing for each instable portion of the density profil 135 ! ------------------------------------------------------------------ 136 137 IF( zwz(1,1) /= 0.e0 ) THEN ! -->> the density profil is statically instable : 138 DO ji = 1, jpi 139 IF( zwy(ji,1) /= 0.e0 ) THEN 98 IF( l_LB_debug ) THEN 99 ! Location of 1 known convection site to follow what's happening in the water column 100 ilc1 = 45 ; jlc1 = 3 ; ! ORCA2 4x4, Antarctic coast, more than 2 unstable portions in the water column... 101 nncpu = 1 ; ! the CPU domain contains the convection spot 102 klc1 = mbkt(ilc1,jlc1) ! bottom of the ocean for debug point... 103 ENDIF 104 105 CALL eos_rab( tsa, zab ) ! after alpha and beta (given on T-points) 106 CALL bn2 ( tsa, zab, zn2 ) ! after Brunt-Vaisala (given on W-points) 107 108 inpcc = 0 109 110 DO jj = 2, jpjm1 ! interior column only 111 DO ji = fs_2, fs_jpim1 112 ! 113 IF( tmask(ji,jj,2) == 1 ) THEN ! At least 2 ocean points 114 ! ! consider one ocean column 115 zvts(:,jp_tem) = tsa(ji,jj,:,jp_tem) ! temperature 116 zvts(:,jp_sal) = tsa(ji,jj,:,jp_sal) ! salinity 117 118 zvab(:,jp_tem) = zab(ji,jj,:,jp_tem) ! Alpha 119 zvab(:,jp_sal) = zab(ji,jj,:,jp_sal) ! Beta 120 zvn2(:) = zn2(ji,jj,:) ! N^2 121 122 IF( l_LB_debug ) THEN !LB debug: 123 lp_monitor_point = .FALSE. 124 IF( ( ji == ilc1 ).AND.( jj == jlc1 ) ) lp_monitor_point = .TRUE. 125 ! writing only if on CPU domain where conv region is: 126 lp_monitor_point = (narea == nncpu).AND.lp_monitor_point 127 ENDIF !LB debug end 128 129 ikbot = mbkt(ji,jj) ! ikbot: ocean bottom T-level 130 ikp = 1 ! because N2 is irrelevant at the surface level (will start at ikp=2) 131 ilayer = 0 132 jiter = 0 133 l_column_treated = .FALSE. 134 135 DO WHILE ( .NOT. l_column_treated ) 140 136 ! 141 ikbot = mbkt(ji,jj) ! ikbot: ocean bottom T-level 137 jiter = jiter + 1 138 139 IF( jiter >= 400 ) EXIT 140 141 l_bottom_reached = .FALSE. 142 143 DO WHILE ( .NOT. l_bottom_reached ) 144 145 ikp = ikp + 1 146 147 !! Testing level ikp for instability 148 !! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 149 IF( zvn2(ikp) < -zn2_zero ) THEN ! Instability found! 150 151 ilayer = ilayer + 1 ! yet another instable portion of the water column found.... 152 153 IF( lp_monitor_point ) THEN 154 WRITE(numout,*) 155 IF( ilayer == 1 .AND. jiter == 1 ) THEN ! first time a column is spoted with an instability 156 WRITE(numout,*) 157 WRITE(numout,*) 'Time step = ',kt,' !!!' 158 ENDIF 159 WRITE(numout,*) ' * Iteration #',jiter,': found instable portion #',ilayer, & 160 & ' in column! Starting at ikp =', ikp 161 WRITE(numout,*) ' *** N2 for point (i,j) = ',ji,' , ',jj 162 DO jk = 1, klc1 163 WRITE(numout,*) jk, zvn2(jk) 164 END DO 165 WRITE(numout,*) 166 ENDIF 167 168 169 IF( jiter == 1 ) inpcc = inpcc + 1 170 171 IF( lp_monitor_point ) WRITE(numout, *) 'Negative N2 at ikp =',ikp,' for layer #', ilayer 172 173 !! ikup is the uppermost point where mixing will start: 174 ikup = ikp - 1 ! ikup is always "at most at ikp-1", less if neutral levels overlying 175 176 !! If the points above ikp-1 have N2 == 0 they must also be mixed: 177 IF( ikp > 2 ) THEN 178 DO jk = ikp-1, 2, -1 179 IF( ABS(zvn2(jk)) < zn2_zero ) THEN 180 ikup = ikup - 1 ! 1 more upper level has N2=0 and must be added for the mixing 181 ELSE 182 EXIT 183 ENDIF 184 END DO 185 ENDIF 186 187 IF( ikup < 1 ) CALL ctl_stop( 'tra_npc : PROBLEM #1') 188 189 zsum_temp = 0._wp 190 zsum_sali = 0._wp 191 zsum_alfa = 0._wp 192 zsum_beta = 0._wp 193 zsum_z = 0._wp 194 195 DO jk = ikup, ikbot ! Inside the instable (and overlying neutral) portion of the column 196 ! 197 zdz = fse3t(ji,jj,jk) 198 zsum_temp = zsum_temp + zvts(jk,jp_tem)*zdz 199 zsum_sali = zsum_sali + zvts(jk,jp_sal)*zdz 200 zsum_alfa = zsum_alfa + zvab(jk,jp_tem)*zdz 201 zsum_beta = zsum_beta + zvab(jk,jp_sal)*zdz 202 zsum_z = zsum_z + zdz 203 ! 204 IF( jk == ikbot ) EXIT ! avoid array-index overshoot in case ikbot = jpk, cause we're calling jk+1 next line 205 !! EXIT when we have reached the last layer that is instable (N2<0) or neutral (N2=0): 206 IF( zvn2(jk+1) > zn2_zero ) EXIT 207 END DO 208 209 ikdown = jk ! for the current unstable layer, ikdown is the deepest point with a negative or neutral N2 210 IF( ikup == ikdown ) CALL ctl_stop( 'tra_npc : PROBLEM #2') 211 212 ! Mixing Temperature, salinity, alpha and beta from ikup to ikdown included: 213 zta = zsum_temp/zsum_z 214 zsa = zsum_sali/zsum_z 215 zalfa = zsum_alfa/zsum_z 216 zbeta = zsum_beta/zsum_z 217 218 IF( lp_monitor_point ) THEN 219 WRITE(numout,*) 'MIXED T, S, alfa and beta between ikup =',ikup, & 220 & ' and ikdown =',ikdown,', in layer #',ilayer 221 WRITE(numout,*) ' => Mean temp. in that portion =', zta 222 WRITE(numout,*) ' => Mean sali. in that portion =', zsa 223 WRITE(numout,*) ' => Mean Alfa in that portion =', zalfa 224 WRITE(numout,*) ' => Mean Beta in that portion =', zbeta 225 ENDIF 226 227 !! Homogenaizing the temperature, salinity, alpha and beta in this portion of the column 228 DO jk = ikup, ikdown 229 zvts(jk,jp_tem) = zta 230 zvts(jk,jp_sal) = zsa 231 zvab(jk,jp_tem) = zalfa 232 zvab(jk,jp_sal) = zbeta 233 END DO 234 235 236 !! Updating N2 in the relvant portion of the water column 237 !! Temperature, Salinity, Alpha and Beta have been homogenized in the unstable portion 238 !! => Need to re-compute N2! will use Alpha and Beta! 239 240 ikup = MAX(2,ikup) ! ikup can never be 1 ! 241 ik_low = MIN(ikdown+1,ikbot) ! we must go 1 point deeper than ikdown! 242 243 DO jk = ikup, ik_low ! we must go 1 point deeper than ikdown! 244 245 !! Interpolating alfa and beta at W point: 246 zrw = (fsdepw(ji,jj,jk ) - fsdept(ji,jj,jk)) & 247 & / (fsdept(ji,jj,jk-1) - fsdept(ji,jj,jk)) 248 zaw = zvab(jk,jp_tem) * (1._wp - zrw) + zvab(jk-1,jp_tem) * zrw 249 zbw = zvab(jk,jp_sal) * (1._wp - zrw) + zvab(jk-1,jp_sal) * zrw 250 251 !! N2 at W point, doing exactly as in eosbn2.F90: 252 zvn2(jk) = grav*( zaw * ( zvts(jk-1,jp_tem) - zvts(jk,jp_tem) ) & 253 & - zbw * ( zvts(jk-1,jp_sal) - zvts(jk,jp_sal) ) ) & 254 & / fse3w(ji,jj,jk) * tmask(ji,jj,jk) 255 256 !! OR, faster => just considering the vertical gradient of density 257 !! as only the signa maters... 258 !zvn2(jk) = ( zaw * ( zvts(jk-1,jp_tem) - zvts(jk,jp_tem) ) & 259 ! & - zbw * ( zvts(jk-1,jp_sal) - zvts(jk,jp_sal) ) ) 260 261 END DO 262 263 ikp = MIN(ikdown+1,ikbot) 264 265 266 ENDIF !IF( zvn2(ikp) < 0. ) 267 268 269 IF( ikp == ikbot ) l_bottom_reached = .TRUE. 270 ! 271 END DO ! DO WHILE ( .NOT. l_bottom_reached ) 272 273 IF( ikp /= ikbot ) CALL ctl_stop( 'tra_npc : PROBLEM #3') 274 275 ! ******* At this stage ikp == ikbot ! ******* 276 277 IF( ilayer > 0 ) THEN !! least an unstable layer has been found 278 ! 279 IF( lp_monitor_point ) THEN 280 WRITE(numout,*) 281 WRITE(numout,*) 'After ',jiter,' iteration(s), we neutralized ',ilayer,' instable layer(s)' 282 WRITE(numout,*) ' ==> N2 at i,j=',ji,',',jj,' now looks like this:' 283 DO jk = 1, klc1 284 WRITE(numout,*) jk, zvn2(jk) 285 END DO 286 WRITE(numout,*) 287 ENDIF 288 ! 289 ikp = 1 ! starting again at the surface for the next iteration 290 ilayer = 0 291 ENDIF 142 292 ! 143 DO jiter = 1, jpk ! vertical iteration 144 ! 145 ! search of ikup : the first static instability from the sea surface 146 ! 147 ik = 0 148 220 CONTINUE 149 ik = ik + 1 150 IF( ik >= ikbot ) GO TO 200 151 zwx(ji,ik) = zrhop(ji,jj,ik) - zrhop(ji,jj,ik+1) 152 IF( zwx(ji,ik) <= 0.e0 ) GO TO 220 153 ikup = ik 154 ! the density profil is instable below ikup 155 ! ikdown : bottom of the instable portion of the density profil 156 ! search of ikdown and vertical mixing from ikup to ikdown 157 ! 158 ze3tot= fse3t(ji,jj,ikup) 159 zta = tsa (ji,jj,ikup,jp_tem) 160 zsa = tsa (ji,jj,ikup,jp_sal) 161 zraua = zrhop(ji,jj,ikup) 162 ! 163 DO jkdown = ikup+1, ikbot-1 164 IF( zraua <= zrhop(ji,jj,jkdown) ) THEN 165 ikdown = jkdown 166 GO TO 240 167 ENDIF 168 ze3dwn = fse3t(ji,jj,jkdown) 169 ze3tot = ze3tot + ze3dwn 170 zta = ( zta*(ze3tot-ze3dwn) + tsa(ji,jj,jkdown,jp_tem)*ze3dwn )/ze3tot 171 zsa = ( zsa*(ze3tot-ze3dwn) + tsa(ji,jj,jkdown,jp_sal)*ze3dwn )/ze3tot 172 zraua = ( zraua*(ze3tot-ze3dwn) + zrhop(ji,jj,jkdown)*ze3dwn )/ze3tot 173 inpci = inpci+1 174 END DO 175 ikdown = ikbot-1 176 240 CONTINUE 177 ! 178 DO jkp = ikup, ikdown-1 179 tsa (ji,jj,jkp,jp_tem) = zta 180 tsa (ji,jj,jkp,jp_sal) = zsa 181 zrhop(ji,jj,jkp ) = zraua 182 END DO 183 IF (ikdown == ikbot-1 .AND. zraua >= zrhop(ji,jj,ikdown) ) THEN 184 tsa (ji,jj,jkp,jp_tem) = zta 185 tsa (ji,jj,jkp,jp_sal) = zsa 186 zrhop(ji,jj,ikdown ) = zraua 187 ENDIF 188 END DO 189 ENDIF 190 200 CONTINUE 191 END DO 192 ! <<-- no more static instability on slab jj 193 ENDIF 194 ! ! =============== 195 END DO ! End of slab 196 ! ! =============== 197 ! 198 IF( l_trdtra ) THEN ! save the Non penetrative mixing trends for diagnostic 199 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:) 200 ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:) 201 CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_npc, ztrdt ) 202 CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_npc, ztrds ) 293 IF( ikp >= ikbot ) l_column_treated = .TRUE. 294 ! 295 END DO ! DO WHILE ( .NOT. l_column_treated ) 296 297 !! Updating tsa: 298 tsa(ji,jj,:,jp_tem) = zvts(:,jp_tem) 299 tsa(ji,jj,:,jp_sal) = zvts(:,jp_sal) 300 301 !! LB: Potentially some other global variable beside theta and S can be treated here 302 !! like BGC tracers. 303 304 IF( lp_monitor_point ) WRITE(numout,*) 305 306 ENDIF ! IF( tmask(ji,jj,3) == 1 ) THEN 307 308 END DO ! ji 309 END DO ! jj 310 ! 311 IF( l_trdtra ) THEN ! send the Non penetrative mixing trends for diagnostic 312 z1_r2dt = 1._wp / (2._wp * rdt) 313 ztrdt(:,:,:) = ( tsa(:,:,:,jp_tem) - ztrdt(:,:,:) ) * z1_r2dt 314 ztrds(:,:,:) = ( tsa(:,:,:,jp_sal) - ztrds(:,:,:) ) * z1_r2dt 315 CALL trd_tra( kt, 'TRA', jp_tem, jptra_npc, ztrdt ) 316 CALL trd_tra( kt, 'TRA', jp_sal, jptra_npc, ztrds ) 203 317 CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds ) 204 318 ENDIF 205 206 ! Lateral boundary conditions on ( ta, sa ) ( Unchanged sign) 207 ! ------------------------------============ 319 ! 208 320 CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1. ) ; CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1. ) 209 210 211 ! 2. non penetrative convective scheme statistics 212 ! ----------------------------------------------- 213 IF( nn_npcp /= 0 .AND. MOD( kt, nn_npcp ) == 0 ) THEN 214 IF(lwp) WRITE(numout,*)' kt=',kt, ' number of statically instable', & 215 & ' water column : ',inpcc, ' number of iteration : ',inpci 321 ! 322 IF( lwp .AND. l_LB_debug ) THEN 323 WRITE(numout,*) 'Exiting tra_npc , kt = ',kt,', => numb. of statically instable water-columns: ', inpcc 324 WRITE(numout,*) 216 325 ENDIF 217 326 ! 218 ENDIF 219 ! 220 CALL wrk_dealloc(jpi, jpj, jpk, zrhop ) 221 CALL wrk_dealloc(jpi, jpk, zwx, zwy, zwz ) 327 CALL wrk_dealloc(jpi, jpj, jpk, zn2 ) 328 CALL wrk_dealloc(jpi, jpj, jpk, 2, zab ) 329 CALL wrk_dealloc(jpk, zvn2 ) 330 CALL wrk_dealloc(jpk, 2, zvts, zvab ) 331 ! 332 ENDIF ! IF( MOD( kt, nn_npc ) == 0 ) THEN 222 333 ! 223 334 IF( nn_timing == 1 ) CALL timing_stop('tra_npc')
Note: See TracChangeset
for help on using the changeset viewer.