- Timestamp:
- 2015-07-09T18:07:16+02:00 (9 years ago)
- Location:
- branches/2015/dev_r5546_CNRS19_HPC_scalability/NEMOGCM/NEMO/LIM_SRC_3
- Files:
-
- 2 edited
-
limhdf.F90 (modified) (4 diffs)
-
limtrp.F90 (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5546_CNRS19_HPC_scalability/NEMOGCM/NEMO/LIM_SRC_3/limhdf.F90
r5429 r5579 28 28 29 29 PUBLIC lim_hdf ! called by lim_trp 30 PUBLIC lim_hdf_multiple ! called by lim_trp 30 31 PUBLIC lim_hdf_init ! called by sbc_lim_init 31 32 … … 124 125 CALL lbc_lnk( zrlx, 'T', 1. ) ! lateral boundary condition 125 126 ! 126 IF ( MOD( iter , nn_convfrq ) == 0 ) THEN ! convergence test every nn_convfrq iterations (perf. optimization)127 IF ( MOD( iter - 1 , nn_convfrq ) == 0 ) THEN ! convergence test every nn_convfrq iterations (perf. optimization) 127 128 zconv = 0._wp 128 129 DO jj = 2, jpjm1 … … 166 167 ! 167 168 END SUBROUTINE lim_hdf 169 170 171 SUBROUTINE lim_hdf_multiple( ptab , ihdf_vars , jpl , nlay_i ) 172 !!------------------------------------------------------------------- 173 !! *** ROUTINE lim_hdf *** 174 !! 175 !! ** purpose : Compute and add the diffusive trend on sea-ice variables 176 !! 177 !! ** method : Second order diffusive operator evaluated using a 178 !! Cranck-Nicholson time Scheme. 179 !! 180 !! ** Action : update ptab with the diffusive contribution 181 !!------------------------------------------------------------------- 182 INTEGER :: jpl, nlay_i, isize, ihdf_vars 183 REAL(wp), DIMENSION(:,:,:), INTENT( inout ),TARGET :: ptab ! Field on which the diffusion is applied 184 REAL(wp), POINTER, DIMENSION(:,:,:) :: pahu3D , pahv3D 185 ! 186 INTEGER :: ji, jj, jk, jl , jm ! dummy loop indices 187 INTEGER :: iter, ierr ! local integers 188 REAL(wp) :: zrlxint ! local scalars 189 REAL(wp), POINTER , DIMENSION ( : ) :: zconv ! local scalars 190 REAL(wp), POINTER , DIMENSION(:,:,:) :: zrlx,zdiv0, ztab0 191 REAL(wp), POINTER , DIMENSION(:,:) :: zflu, zflv, zdiv 192 CHARACTER(lc) :: charout ! local character 193 REAL(wp), PARAMETER :: zrelax = 0.5_wp ! relaxation constant for iterative procedure 194 REAL(wp), PARAMETER :: zalfa = 0.5_wp ! =1.0/0.5/0.0 = implicit/Cranck-Nicholson/explicit 195 INTEGER , PARAMETER :: its = 100 ! Maximum number of iteration 196 !!------------------------------------------------------------------- 197 TYPE(arrayptr) , ALLOCATABLE, DIMENSION(:) :: pt2d_array, zrlx_array 198 CHARACTER(len=1) , ALLOCATABLE, DIMENSION(:) :: type_array ! define the nature of ptab array grid-points 199 ! ! = T , U , V , F , W and I points 200 REAL(wp) , ALLOCATABLE, DIMENSION(:) :: psgn_array ! =-1 the sign change across the north fold boundary 201 202 !!--------------------------------------------------------------------- 203 204 ! !== Initialisation ==! 205 isize = jpl*(ihdf_vars+nlay_i) 206 ALLOCATE( zconv (isize) ) 207 ALLOCATE( pt2d_array(isize) , zrlx_array(isize) ) 208 ALLOCATE( type_array(isize) ) 209 ALLOCATE( psgn_array(isize) ) 210 211 CALL wrk_alloc( jpi, jpj, isize, zrlx, zdiv0, ztab0 ) 212 CALL wrk_alloc( jpi, jpj, zflu, zflv, zdiv ) 213 CALL wrk_alloc( jpi, jpj, jpl, pahu3D , pahv3D ) 214 215 216 DO jl = 1 , jpl 217 jm = (jl-1)*(ihdf_vars+nlay_i)+1 218 DO jj = 1, jpjm1 ! NB: has not to be defined on jpj line and jpi row 219 DO ji = 1 , fs_jpim1 ! vector opt. 220 pahu3D(ji,jj,jl) = ( 1._wp - MAX( 0._wp, SIGN( 1._wp, -ptab(ji ,jj,jm) ) ) ) & 221 & * ( 1._wp - MAX( 0._wp, SIGN( 1._wp, -ptab(ji+1, jj, jm ) ) ) ) * ahiu(ji,jj) 222 pahv3D(ji,jj,jl) = ( 1._wp - MAX( 0._wp, SIGN( 1._wp, -ptab(ji, jj, jm ) ) ) ) & 223 & * ( 1._wp - MAX( 0._wp, SIGN( 1._wp,- ptab(ji, jj+1, jm ) ) ) ) * ahiv(ji,jj) 224 END DO 225 END DO 226 END DO 227 228 DO jk= 1 , isize 229 pt2d_array(jk)%pt2d=>ptab(:,:,jk) 230 zrlx_array(jk)%pt2d=>zrlx(:,:,jk) 231 type_array(jk)='T' 232 psgn_array(jk)=1. 233 END DO 234 235 ! 236 IF( linit ) THEN ! Metric coefficient (compute at the first call and saved in efact) 237 ALLOCATE( efact(jpi,jpj) , STAT=ierr ) 238 IF( lk_mpp ) CALL mpp_sum( ierr ) 239 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'lim_hdf : unable to allocate arrays' ) 240 DO jj = 2, jpjm1 241 DO ji = fs_2 , fs_jpim1 ! vector opt. 242 efact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj) + e1v(ji,jj) + e1v(ji,jj-1) ) * r1_e12t(ji,jj) 243 END DO 244 END DO 245 linit = .FALSE. 246 ENDIF 247 ! ! Time integration parameters 248 ! 249 zflu (jpi,: ) = 0._wp 250 zflv (jpi,: ) = 0._wp 251 252 DO jk=1 , isize 253 ztab0(:, : , jk ) = ptab(:,:,jk) ! Arrays initialization 254 zdiv0(:, 1 , jk ) = 0._wp 255 zdiv0(:,jpj, jk ) = 0._wp 256 zdiv0(1, :, jk ) = 0._wp 257 zdiv0(jpi,:, jk ) = 0._wp 258 END DO 259 260 zconv = 1._wp !== horizontal diffusion using a Crant-Nicholson scheme ==! 261 iter = 0 262 ! 263 DO WHILE( MAXVAL(zconv(:)) > ( 2._wp * 1.e-04 ) .AND. iter <= its ) ! Sub-time step loop 264 ! 265 iter = iter + 1 ! incrementation of the sub-time step number 266 ! 267 268 DO jk = 1 , isize 269 jl = (jk-1) /( ihdf_vars+nlay_i)+1 270 IF (zconv(jk) > ( 2._wp * 1.e-04 )) THEN 271 DO jj = 1, jpjm1 ! diffusive fluxes in U- and V- direction 272 DO ji = 1 , fs_jpim1 ! vector opt. 273 zflu(ji,jj) = pahu3D(ji,jj,jl) * e2u(ji,jj) * r1_e1u(ji,jj) * ( ptab(ji+1,jj,jk) - ptab(ji,jj,jk) ) 274 zflv(ji,jj) = pahv3D(ji,jj,jl) * e1v(ji,jj) * r1_e2v(ji,jj) * ( ptab(ji,jj+1,jk) - ptab(ji,jj,jk) ) 275 END DO 276 END DO 277 ! 278 DO jj= 2, jpjm1 ! diffusive trend : divergence of the fluxes 279 DO ji = fs_2 , fs_jpim1 ! vector opt. 280 zdiv(ji,jj) = ( zflu(ji,jj) - zflu(ji-1,jj) + zflv(ji,jj) - zflv(ji,jj-1) ) * r1_e12t(ji,jj) 281 END DO 282 END DO 283 ! 284 IF( iter == 1 ) zdiv0(:,:,jk) = zdiv(:,:) ! save the 1st evaluation of the diffusive trend in zdiv0 285 ! 286 DO jj = 2, jpjm1 ! iterative evaluation 287 DO ji = fs_2 , fs_jpim1 ! vector opt. 288 zrlxint = ( ztab0(ji,jj,jk) & 289 & + rdt_ice * ( zalfa * ( zdiv(ji,jj) + efact(ji,jj) * ptab(ji,jj,jk) ) & 290 & + ( 1.0 - zalfa ) * zdiv0(ji,jj,jk) ) & 291 & ) / ( 1.0 + zalfa * rdt_ice * efact(ji,jj) ) 292 zrlx(ji,jj,jk) = ptab(ji,jj,jk) + zrelax * ( zrlxint - ptab(ji,jj,jk) ) 293 END DO 294 END DO 295 END IF 296 297 END DO 298 299 CALL lbc_lnk_multi( zrlx_array, type_array , psgn_array , isize ) ! Multiple interchange of all the variables 300 ! 301 302 IF ( MOD( iter-1 , nn_convfrq ) == 0 ) THEN !Convergence test every nn_convfrq iterations (perf. optimization ) 303 DO jk=1,isize 304 zconv(jk) = 0._wp ! convergence test 305 DO jj = 2, jpjm1 306 DO ji = fs_2, fs_jpim1 307 zconv(jk) = MAX( zconv(jk), ABS( zrlx(ji,jj,jk) - ptab(ji,jj,jk) ) ) 308 END DO 309 END DO 310 END DO 311 IF( lk_mpp ) CALL mpp_max_multiple( zconv , isize ) ! max over the global domain for all the variables 312 ENDIF 313 ! 314 DO jk=1,isize 315 ptab(:,:,jk) = zrlx(:,:,jk) 316 END DO 317 ! 318 END DO ! end of sub-time step loop 319 320 ! ----------------------- 321 !!! final step (clem) !!! 322 DO jk = 1, isize 323 jl = (jk-1) /( ihdf_vars+nlay_i)+1 324 DO jj = 1, jpjm1 ! diffusive fluxes in U- and V- direction 325 DO ji = 1 , fs_jpim1 ! vector opt. 326 zflu(ji,jj) = pahu3D(ji,jj,jl) * e2u(ji,jj) * r1_e1u(ji,jj) * ( ptab(ji+1,jj,jk) - ptab(ji,jj,jk) ) 327 zflv(ji,jj) = pahv3D(ji,jj,jl) * e1v(ji,jj) * r1_e2v(ji,jj) * ( ptab(ji,jj+1,jk) - ptab(ji,jj,jk) ) 328 END DO 329 END DO 330 ! 331 DO jj= 2, jpjm1 ! diffusive trend : divergence of the fluxes 332 DO ji = fs_2 , fs_jpim1 ! vector opt. 333 zdiv(ji,jj) = ( zflu(ji,jj) - zflu(ji-1,jj) + zflv(ji,jj) - zflv(ji,jj-1) ) * r1_e12t(ji,jj) 334 ptab(ji,jj,jk) = ztab0(ji,jj,jk) + 0.5 * ( zdiv(ji,jj) + zdiv0(ji,jj,jk) ) 335 END DO 336 END DO 337 END DO 338 339 CALL lbc_lnk_multi( pt2d_array, type_array , psgn_array , isize ) ! Multiple interchange of all the variables 340 341 !!! final step (clem) !!! 342 ! ----------------------- 343 344 IF(ln_ctl) THEN 345 DO jk = 1 , isize 346 zrlx(:,:,jk) = ptab(:,:,jk) - ztab0(:,:,jk) 347 WRITE(charout,FMT="(' lim_hdf : zconv =',D23.16, ' iter =',I4,2X)") zconv, iter 348 CALL prt_ctl( tab2d_1=zrlx(:,:,jk), clinfo1=charout ) 349 END DO 350 ENDIF 351 ! 352 CALL wrk_dealloc( jpi, jpj, isize, zrlx, zdiv0, ztab0 ) 353 CALL wrk_dealloc( jpi, jpj, zflu, zflv, zdiv ) 354 CALL wrk_dealloc( jpi, jpj, jpl, pahu3D , pahv3D ) 355 356 DEALLOCATE( zconv ) 357 DEALLOCATE( pt2d_array , zrlx_array ) 358 DEALLOCATE( type_array ) 359 DEALLOCATE( psgn_array ) 360 ! 361 END SUBROUTINE lim_hdf_multiple 362 168 363 169 364 … … 179 374 !!------------------------------------------------------------------- 180 375 INTEGER :: ios ! Local integer output status for namelist read 181 NAMELIST/namicehdf/ nn_convfrq 376 NAMELIST/namicehdf/ nn_convfrq 182 377 !!------------------------------------------------------------------- 183 378 ! -
branches/2015/dev_r5546_CNRS19_HPC_scalability/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90
r5202 r5579 63 63 INTEGER, INTENT(in) :: kt ! number of iteration 64 64 ! 65 INTEGER :: ji, jj, jk, j l, jt ! dummy loop indices65 INTEGER :: ji, jj, jk, jm , jl, jt ! dummy loop indices 66 66 INTEGER :: initad ! number of sub-timestep for the advection 67 67 REAL(wp) :: zcfl , zusnit ! - - … … 75 75 REAL(wp), POINTER, DIMENSION(:,:,:) :: zhimax ! old ice thickness 76 76 REAL(wp), POINTER, DIMENSION(:,:) :: zatold, zeiold, zesold ! old concentration, enthalpies 77 REAL(wp), POINTER, DIMENSION(:,:,:) :: zhdfptab 77 78 REAL(wp) :: zdv, zvi, zvs, zsmv, zes, zei 78 79 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 80 !!--------------------------------------------------------------------- 81 INTEGER :: ihdf_vars = 6 !!Number of variables in which we apply horizontal diffusion 82 !! inside limtrp for each ice category , not counting the 83 !! variables corresponding to ice_layers 79 84 !!--------------------------------------------------------------------- 80 85 IF( nn_timing == 1 ) CALL timing_start('limtrp') … … 85 90 CALL wrk_alloc( jpi,jpj,nlay_i,jpl, z0ei ) 86 91 CALL wrk_alloc( jpi,jpj,jpl, zhimax, zviold, zvsold, zsmvold ) 92 CALL wrk_alloc( jpi,jpj,jpl*(ihdf_vars + nlay_i),zhdfptab) 87 93 88 94 IF( numit == nstart .AND. lwp ) THEN … … 170 176 z0oi (:,:,jl) = oa_i (:,:,jl) * e12t(:,:) ! Age content 171 177 z0es (:,:,jl) = e_s (:,:,1,jl) * e12t(:,:) ! Snow heat content 172 178 DO jk = 1, nlay_i 173 179 z0ei (:,:,jk,jl) = e_i (:,:,jk,jl) * e12t(:,:) ! Ice heat content 174 180 END DO … … 305 311 !------------------------------------ 306 312 DO jl = 1, jpl 307 ! ! Masked eddy diffusivity coefficient at ocean U- and V-points 308 DO jj = 1, jpjm1 ! NB: has not to be defined on jpj line and jpi row 309 DO ji = 1 , fs_jpim1 ! vector opt. 310 pahu(ji,jj) = ( 1._wp - MAX( 0._wp, SIGN( 1._wp, -a_i(ji ,jj,jl) ) ) ) & 311 & * ( 1._wp - MAX( 0._wp, SIGN( 1._wp, -a_i(ji+1,jj,jl) ) ) ) * ahiu(ji,jj) 312 pahv(ji,jj) = ( 1._wp - MAX( 0._wp, SIGN( 1._wp, -a_i(ji,jj ,jl) ) ) ) & 313 & * ( 1._wp - MAX( 0._wp, SIGN( 1._wp,- a_i(ji,jj+1,jl) ) ) ) * ahiv(ji,jj) 314 END DO 315 END DO 316 317 CALL lim_hdf( v_i (:,:, jl) ) 318 CALL lim_hdf( v_s (:,:, jl) ) 319 CALL lim_hdf( smv_i(:,:, jl) ) 320 CALL lim_hdf( oa_i (:,:, jl) ) 321 CALL lim_hdf( a_i (:,:, jl) ) 322 CALL lim_hdf( e_s (:,:,1,jl) ) 313 jm=(jl-1)*(ihdf_vars+nlay_i)+1 314 zhdfptab(:,:,jm)= a_i (:,:, jl); jm = jm + 1 ! IMPORTANT a_i must go in the first position because 315 ! it is needed at this position inside lim_hdf_multiple. 316 zhdfptab(:,:,jm)= v_i (:,:, jl); jm = jm + 1 317 zhdfptab(:,:,jm)= v_s (:,:, jl); jm = jm + 1 318 zhdfptab(:,:,jm)= smv_i(:,:, jl); jm = jm + 1 319 zhdfptab(:,:,jm)= oa_i (:,:, jl); jm = jm + 1 320 zhdfptab(:,:,jm)= e_s (:,:,1,jl); jm = jm + 1 321 ! Sample of adding more variables to apply lim_hdf using lim_hdf_multiple optimization--- 322 ! zhdfptab(:,:,jm) = variable_1 (:,:,1,jl); jm = jm + 1 323 ! zhdfptab(:,:,jm) = variable_2 (:,:,1,jl); jm = jm + 1 324 ! 325 ! and in this example the parameter ihdf_vars musb be changed to 8 (necessary for allocation) 326 !---------------------------------------------------------------------------------------- 323 327 DO jk = 1, nlay_i 324 CALL lim_hdf( e_i(:,:,jk,jl) ) 325 END DO 326 END DO 328 zhdfptab(:,:,jm)=e_i(:,:,jk,jl); jm= jm+1 329 END DO 330 END DO 331 CALL lim_hdf_multiple( zhdfptab, ihdf_vars, jpl, nlay_i) 332 333 DO jl = 1, jpl 334 jm=(jl-1)*(ihdf_vars+nlay_i)+1 335 a_i (:,:, jl) = zhdfptab(:,:,jm); jm = jm + 1 336 v_i (:,:, jl) = zhdfptab(:,:,jm); jm = jm + 1 337 v_s (:,:, jl) = zhdfptab(:,:,jm); jm = jm + 1 338 smv_i(:,:, jl) = zhdfptab(:,:,jm); jm = jm + 1 339 oa_i (:,:, jl) = zhdfptab(:,:,jm); jm = jm + 1 340 e_s (:,:,1,jl) = zhdfptab(:,:,jm); jm = jm + 1 341 ! Sample of adding more variables to apply lim_hdf--------- 342 ! variable_1 (:,:,1,jl) = zhdfptab(:,:, jm ) ; jm + 1 343 ! variable_2 (:,:,1,jl) = zhdfptab(:,:, jm ) ; jm + 1 344 !----------------------------------------------------------- 345 DO jk = 1, nlay_i 346 e_i(:,:,jk,jl) = zhdfptab(:,:,jm);jm= jm + 1 347 END DO 348 END DO 349 327 350 328 351 !------------------------------------------------------------------------------! … … 464 487 CALL wrk_dealloc( jpi,jpj,nlay_i,jpl, z0ei ) 465 488 CALL wrk_dealloc( jpi,jpj,jpl, zviold, zvsold, zhimax, zsmvold ) 489 CALL wrk_dealloc( jpi,jpj,jpl*(ihdf_vars+nlay_i),zhdfptab) 466 490 ! 467 491 IF( nn_timing == 1 ) CALL timing_stop('limtrp')
Note: See TracChangeset
for help on using the changeset viewer.
