Changeset 2648
- Timestamp:
- 2011-03-03T17:13:18+01:00 (14 years ago)
- Location:
- branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC
- Files:
-
- 1 added
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/dommsk.F90
r2528 r2648 12 12 USE oce ! ocean dynamics and tracers 13 13 USE dom_oce ! ocean space and time domain 14 USE in_out_manager ! I/O manager14 USE lib_mpp 15 15 16 16 IMPLICIT NONE … … 19 19 PUBLIC dom_msk ! routine called by inidom.F90 20 20 21 #if defined key_degrad 22 !! ------------------------------------------------ 23 !! Degradation method 24 !! -------------------------------------------------- 25 REAL(wp), PUBLIC, DIMENSION (jpi,jpj,jpk) :: facvol !! volume for degraded regions 26 #endif 21 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: facvol !! volume for degraded regions 27 22 28 23 !! * Substitutions … … 34 29 !!---------------------------------------------------------------------- 35 30 CONTAINS 31 36 32 37 33 SUBROUTINE dom_msk … … 49 45 !! tpol : ??? 50 46 !!---------------------------------------------------------------------- 47 USE wrk_nemo, ONLY: iwrk_in_use, iwrk_not_released 48 USE wrk_nemo, ONLY: imsk => iwrk_2d_1 51 49 INTEGER :: ji, jk ! dummy loop indices 52 50 INTEGER :: iif, iil, ijf, ijl ! local integers 53 INTEGER, DIMENSION(jpi,jpj) :: imsk ! 2D workspace54 51 !!--------------------------------------------------------------------- 55 52 ! 53 IF( iwrk_in_use(2, 1) ) THEN 54 CALL ctl_stop('dom_msk: ERROR: requested workspace arrays unavailable.') ; RETURN 55 END IF 56 ! 57 CALL dom_msk_alloc 58 56 59 ! Interior domain mask (used for global sum) 57 60 ! -------------------- … … 95 98 ENDIF 96 99 ! 100 IF( iwrk_not_released(2, 1) ) CALL ctl_stop('dom_msk: failed to release workspace arrays.') 101 ! 97 102 END SUBROUTINE dom_msk 103 104 SUBROUTINE dom_msk_alloc() 105 !!--------------------------------------------------------------------- 106 !! *** ROUTINE dom_msk_alloc *** 107 !!--------------------------------------------------------------------- 108 #if defined key_degrad 109 INTEGER :: ierr 110 111 ALLOCATE( facvol(jpi,jpj,jpk), STAT=ierr ) 112 IF( ierr /= 0 ) & 113 & CALL ctl_stop('STOP', 'dom_msk : unable to allocate facvol array') 114 #endif 115 116 END SUBROUTINE dom_msk_alloc 98 117 99 118 !!====================================================================== -
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/domrea.F90
r2528 r2648 16 16 USE dommsk ! domain: masks 17 17 USE lbclnk ! lateral boundary condition - MPP exchanges 18 USE in_out_manager ! I/O manager 18 USE lib_mpp 19 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 19 20 20 21 IMPLICIT NONE … … 53 54 !!---------------------------------------------------------------------- 54 55 USE iom 56 USE wrk_nemo, ONLY: zmbk => wrk_2d_1, zprt => wrk_2d_2, zprw => wrk_2d_3 55 57 !! 56 58 INTEGER :: ji, jj, jk ! dummy loop indices 57 59 INTEGER :: ik, inum0 , inum1 , inum2 , inum3 , inum4 ! local integers 58 60 REAL(wp) :: zrefdep ! local real 59 REAL(wp), DIMENSION(jpi,jpj) :: zmbk, zprt, zprw ! 2D workspace60 61 !!---------------------------------------------------------------------- 61 62 … … 63 64 IF(lwp) WRITE(numout,*) 'dom_rea : read NetCDF mesh and mask information file(s)' 64 65 IF(lwp) WRITE(numout,*) '~~~~~~~' 66 67 IF( wrk_in_use(2, 1,2,3) ) THEN 68 CALL ctl_stop('dom_rea: ERROR: requested workspace arrays unavailable.') ; RETURN 69 END IF 65 70 66 71 zmbk(:,:) = 0._wp … … 141 146 CALL iom_get( inum3, jpdom_data, 'e2u', e2u ) 142 147 CALL iom_get( inum3, jpdom_data, 'e2v', e2v ) 148 149 e1e2t(:,:) = e1t(:,:) * e2t(:,:) ! surface at T-points 143 150 144 151 CALL iom_get( inum3, jpdom_data, 'ff', ff ) … … 314 321 END SELECT 315 322 ! 323 IF( wrk_not_released(2, 1,2,3) ) CALL ctl_stop('dom_rea:failed to release workspace arrays.') 324 ! 316 325 END SUBROUTINE dom_rea 317 326 … … 327 336 !! ** Action : - update mbathy: level bathymetry (in level index) 328 337 !!---------------------------------------------------------------------- 338 USE wrk_nemo, ONLY: zmbk => wrk_2d_4 339 ! 329 340 INTEGER :: ji, jj ! dummy loop indices 330 REAL(wp), DIMENSION(jpi,jpj) :: zmbk ! 2D workspace331 !!---------------------------------------------------------------------- 341 !!---------------------------------------------------------------------- 342 332 343 ! 333 344 IF(lwp) WRITE(numout,*) 334 345 IF(lwp) WRITE(numout,*) ' zgr_bot_level : ocean bottom k-index of T-, U-, V- and W-levels ' 335 346 IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~' 347 ! 348 IF( wrk_in_use(2, 4) ) THEN 349 CALL ctl_stop('dom_rea: ERROR: requested workspace arrays unavailable.') ; RETURN 350 END IF 336 351 ! 337 352 mbkt(:,:) = MAX( mbathy(:,:) , 1 ) ! bottom k-index of T-level (=1 over land) … … 347 362 zmbk(:,:) = REAL( mbkv(:,:), wp ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 348 363 ! 364 IF( wrk_not_released(2, 4) ) CALL ctl_stop('dom_rea:failed to release workspace arrays.') 365 ! 349 366 END SUBROUTINE zgr_bot_level 350 367 -
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/dtadyn.F90
r2559 r2648 63 63 INTEGER :: numfl_t, numfl_u, numfl_v, numfl_w 64 64 65 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: tdta ! temperature at two consecutive times66 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: sdta ! salinity at two consecutive times67 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: udta ! zonal velocity at two consecutive times68 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: vdta ! meridional velocity at two consecutive times69 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wdta ! vertical velocity at two consecutive times70 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: avtdta ! vertical diffusivity coefficient71 72 REAL(wp), DIMENSION(jpi,jpj ,2) :: hmlddta ! mixed layer depth at two consecutive times73 REAL(wp), DIMENSION(jpi,jpj ,2) :: wspddta ! wind speed at two consecutive times74 REAL(wp), DIMENSION(jpi,jpj ,2) :: frlddta ! sea-ice fraction at two consecutive times75 REAL(wp), DIMENSION(jpi,jpj ,2) :: empdta ! E-P at two consecutive times76 REAL(wp), DIMENSION(jpi,jpj ,2) :: qsrdta ! short wave heat flux at two consecutive times77 REAL(wp), DIMENSION(jpi,jpj ,2) :: bblxdta ! frequency of bbl in the x direction at 2 consecutive times78 REAL(wp), DIMENSION(jpi,jpj ,2) :: bblydta ! frequency of bbl in the y direction at 2 consecutive times65 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: tdta ! temperature at two consecutive times 66 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: sdta ! salinity at two consecutive times 67 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: udta ! zonal velocity at two consecutive times 68 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: vdta ! meridional velocity at two consecutive times 69 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wdta ! vertical velocity at two consecutive times 70 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: avtdta ! vertical diffusivity coefficient 71 72 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hmlddta ! mixed layer depth at two consecutive times 73 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: wspddta ! wind speed at two consecutive times 74 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: frlddta ! sea-ice fraction at two consecutive times 75 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: empdta ! E-P at two consecutive times 76 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: qsrdta ! short wave heat flux at two consecutive times 77 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: bblxdta ! frequency of bbl in the x direction at 2 consecutive times 78 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: bblydta ! frequency of bbl in the y direction at 2 consecutive times 79 79 LOGICAL :: l_offbbl 80 80 #if defined key_ldfslp 81 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: uslpdta ! zonal isopycnal slopes82 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: vslpdta ! meridional isopycnal slopes83 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wslpidta ! zonal diapycnal slopes84 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wslpjdta ! meridional diapycnal slopes81 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: uslpdta ! zonal isopycnal slopes 82 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: vslpdta ! meridional isopycnal slopes 83 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wslpidta ! zonal diapycnal slopes 84 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wslpjdta ! meridional diapycnal slopes 85 85 #endif 86 86 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 87 REAL(wp), DIMENSION(jpi,jpj ,2) :: aeiwdta ! G&M coefficient87 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: aeiwdta ! G&M coefficient 88 88 #endif 89 89 #if defined key_degrad 90 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: ahtudta, ahtvdta, ahtwdta ! Lateral diffusivity90 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: ahtudta, ahtvdta, ahtwdta ! Lateral diffusivity 91 91 # if defined key_traldf_eiv 92 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: aeiudta, aeivdta, aeiwdta ! G&M coefficient92 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: aeiudta, aeivdta, aeiwdta ! G&M coefficient 93 93 # endif 94 94 #endif … … 296 296 END SUBROUTINE dta_dyn 297 297 298 INTEGER FUNCTION dta_dyn_alloc() 299 !!--------------------------------------------------------------------- 300 !! *** ROUTINE dta_dyn_alloc *** 301 !!--------------------------------------------------------------------- 302 303 ALLOCATE( tdta (jpi,jpj,jpk,2), sdta (jpi,jpj,jpk,2), & 304 & udta (jpi,jpj,jpk,2), vdta (jpi,jpj,jpk,2), & 305 & wdta (jpi,jpj,jpk,2), avtdta (jpi,jpj,jpk,2), & 306 #if defined key_ldfslp 307 & uslpdta (jpi,jpj,jpk,2), vslpdta (jpi,jpj,jpk,2), & 308 & wslpidta(jpi,jpj,jpk,2), wslpjdta(jpi,jpj,jpk,2), & 309 #endif 310 #if defined key_degrad 311 & ahtudta (jpi,jpj,jpk,2), ahtvdta (jpi,jpj,jpk,2), & 312 & ahtwdta (jpi,jpj,jpk,2), & 313 # if defined key_traldf_eiv 314 & aeiudta (jpi,jpj,jpk,2), aeivdta (jpi,jpj,jpk,2), & 315 & aeiwdta (jpi,jpj,jpk,2), & 316 # endif 317 #endif 318 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 319 & aeiwdta (jpi,jpj, 2), & 320 #endif 321 322 & hmlddta (jpi,jpj, 2), wspddta (jpi,jpj, 2), & 323 & frlddta (jpi,jpj, 2), qsrdta (jpi,jpj, 2), & 324 & empdta (jpi,jpj, 2), STAT=dta_dyn_alloc ) 325 326 IF( dta_dyn_alloc /= 0 ) CALL ctl_warn('dta_dyn_alloc: failed to allocate facvol array.') 327 328 END FUNCTION dta_dyn_alloc 298 329 299 330 SUBROUTINE dynrea( kt, kenr ) … … 305 336 !! ** Method : READ the kenr records of DATA and store in udta(...,2), .... 306 337 !!---------------------------------------------------------------------- 338 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 339 USE wrk_nemo, ONLY: zu => wrk_3d_1 , zv => wrk_3d_2 , zw => wrk_3d_3 340 USE wrk_nemo, ONLY: zt => wrk_3d_4 , zs => wrk_3d_5 341 USE wrk_nemo, ONLY: zavt => wrk_3d_6 , zhdiv => wrk_3d_7 342 USE wrk_nemo, ONLY: zahtu => wrk_3d_8 , zahtv => wrk_3d_9 , zahtw => wrk_3d_10 343 USE wrk_nemo, ONLY: zaeiu => wrk_3d_11, zaeiv => wrk_3d_12, zaeiw => wrk_3d_13 344 ! 345 USE wrk_nemo, ONLY: zemp => wrk_2d_1 , zqsr => wrk_2d_2 , zmld => wrk_2d_3 346 USE wrk_nemo, ONLY: zice => wrk_2d_4 , zwspd => wrk_2d_5 347 USE wrk_nemo, ONLY: ztaux => wrk_2d_6 , ztauy => wrk_2d_7 348 USE wrk_nemo, ONLY: zbblx => wrk_2d_8 , zbbly => wrk_2d_9 349 USE wrk_nemo, ONLY: zaeiw2d => wrk_2d_10 350 ! 307 351 INTEGER, INTENT(in) :: kt, kenr ! time index 308 352 !! 309 353 INTEGER :: jkenr 310 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zu, zv, zw, zt, zs, zavt , zhdiv ! 3D workspace 311 REAL(wp), DIMENSION(jpi,jpj) :: zemp, zqsr, zmld, zice, zwspd, ztaux, ztauy ! 2D workspace 312 REAL(wp), DIMENSION(jpi,jpj) :: zbblx, zbbly 313 314 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 315 REAL(wp), DIMENSION(jpi,jpj) :: zaeiw 316 #endif 317 #if defined key_degrad 318 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zahtu, zahtv, zahtw ! Lateral diffusivity 319 # if defined key_traldf_eiv 320 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zaeiu, zaeiv, zaeiw ! G&M coefficient 321 # endif 322 #endif 323 !!---------------------------------------------------------------------- 324 325 ! 0. Initialization 354 !!---------------------------------------------------------------------- 355 356 ! 0. Memory allocation 357 IF( wrk_in_use(3, 1,2,3,4,5,6,7,8,9,10,11,12,13) .OR. & 358 wrk_in_use(2, 1,2,3,4,5,6,7,8,9,10) ) THEN 359 CALL ctl_stop('domrea/dta_dyn: requested workspace arrays unavailable') ; RETURN 360 END IF 326 361 327 362 ! cas d'un fichier non periodique : on utilise deux fois le premier et … … 390 425 391 426 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 392 CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw 427 CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw2d(:,: ), jkenr ) 393 428 #endif 394 429 … … 413 448 414 449 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 415 aeiwdta(:,:,2) = zaeiw (:,:) * tmask(:,:,1)450 aeiwdta(:,:,2) = zaeiw2d(:,:) * tmask(:,:,1) 416 451 #endif 417 452 … … 451 486 ENDIF 452 487 ! 488 IF( wrk_not_released(3, 1,2,3,4,5,6,7,8,9,10,11,12,13) .OR. & 489 wrk_not_released(2, 1,2,3,4,5,6,7,8,9,10) ) THEN 490 CALL ctl_stop('domrea/dta_dyn: failed to release workspace arrays.') 491 END IF 492 ! 453 493 END SUBROUTINE dynrea 454 494 … … 462 502 !! ** Method : 463 503 !!---------------------------------------------------------------------- 464 REAL(wp) :: znspyr !: number of time step per year 504 REAL(wp) :: znspyr !: number of time step per year 505 INTEGER :: ierr 465 506 !! 466 507 NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, lperdyn, & 467 508 & cfile_grid_T, cfile_grid_U, cfile_grid_V, cfile_grid_W 468 509 !!---------------------------------------------------------------------- 510 511 ierr = dta_dyn_alloc() 512 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'dta_dyn_alloc : unable to allocate standard ocean arrays' ) 469 513 470 514 ! Define the dynamical input parameters -
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90
r2574 r2648 5 5 !!====================================================================== 6 6 !! History : 3.3 ! 2010-05 (C. Ethe) Full reorganization of the off-line: phasing with the on-line 7 !! 4.0 ! 2011-01 (C. Ethe, A. R. Porter, STFC Daresbury) dynamical allocation 7 8 !!---------------------------------------------------------------------- 8 9 … … 26 27 USE trabbl ! bottom boundary layer (tra_bbl_init routine) 27 28 USE zdfini ! vertical physics: initialization 29 USE sbcmod ! surface boundary condition (sbc_init routine) 28 30 USE phycst ! physical constant (par_cst routine) 29 31 USE dtadyn ! Lecture and Interpolation of the dynamical fields … … 137 139 ! !--------------------------------------------! 138 140 #if defined key_iomput 139 CALL init_ioclient( ilocal_comm ) ! nemo local communicator (used or not) given bythe io_server140 narea = mynode( cltxt, ilocal_comm )! Nodes selection141 CALL init_ioclient( ilocal_comm ) ! exchange io_server nemo local communicator with the io_server 142 narea = mynode( cltxt, numnam, nstop, ilocal_comm ) ! Nodes selection 141 143 #else 142 narea = mynode( cltxt ) ! Nodes selection (control print return in cltxt) 144 ilocal_comm = 0 145 narea = mynode( cltxt, numnam, nstop ) ! Nodes selection (control print return in cltxt) 143 146 #endif 147 144 148 narea = narea + 1 ! mynode return the rank of proc (0 --> jpnij -1 ) 145 149 146 150 lwp = (narea == 1) .OR. ln_ctl ! control of all listing output print 151 152 ! Decide on size of grid now that we have our communicator size 153 ! If we're not using dynamic memory then mpp_partition does nothing. 154 155 #if defined key_mpp_mpi || defined key_mpp_shmem 156 CALL nemo_partition(mppsize) 157 #else 158 jpni = 1 159 jpnj = 1 160 jpnij = jpni*jpnj 161 #endif 162 ! Calculate domain dimensions given calculated jpni and jpnj 163 ! This used to be done in par_oce.F90 when they were parameters rather 164 ! than variables 165 jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci !: first dim. 166 jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj !: second dim. 167 jpim1 = jpi-1 !: inner domain indices 168 jpjm1 = jpj-1 !: " " 169 jpkm1 = jpk-1 !: " " 170 jpij = jpi*jpj !: jpi x j 171 172 ! Now we know the dimensions of the grid, allocate arrays 173 CALL nemo_alloc() 147 174 148 175 IF(lwp) THEN ! open listing units … … 182 209 183 210 ! ! Ocean physics 211 CALL sbc_init ! Forcings : surface module 184 212 #if ! defined key_degrad 185 213 CALL ldf_tra_init ! Lateral ocean tracer physics … … 307 335 END SUBROUTINE nemo_closefile 308 336 337 SUBROUTINE nemo_alloc 338 !!---------------------------------------------------------------------- 339 !! *** ROUTINE nemo_alloc *** 340 !! 341 !! ** Purpose : Allocate all the dynamic arrays of the OPA modules 342 !! 343 !! ** Method : 344 !!---------------------------------------------------------------------- 345 USE diawri, ONLY: dia_wri_alloc 346 USE dom_oce, ONLY: dom_oce_alloc 347 USE zdf_oce, ONLY: zdf_oce_alloc 348 USE zdfmxl, ONLY: zdf_mxl_alloc 349 USE ldftra_oce, ONLY: ldftra_oce_alloc 350 USE trc_oce, ONLY: trc_oce_alloc 351 352 USE wrk_nemo, ONLY: wrk_alloc 353 354 INTEGER :: ierr 355 !!---------------------------------------------------------------------- 356 357 ierr = oce_alloc () ! ocean 358 ierr = ierr + dia_wri_alloc () 359 ierr = ierr + dom_oce_alloc () ! ocean domain 360 ierr = ierr + ldftra_oce_alloc() ! ocean lateral physics : tracers 361 ierr = ierr + zdf_oce_alloc () ! ocean vertical physics 362 ierr = ierr + zdf_mxl_alloc () ! ocean vertical physics 363 ! 364 ierr = ierr + lib_mpp_alloc (numout) ! mpp exchanges 365 ierr = ierr + trc_oce_alloc () ! shared TRC / TRA arrays 366 ierr = ierr + wrk_alloc(numout, lwp) 367 368 IF( lk_mpp ) CALL mpp_sum( ierr ) 369 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'nemo_alloc : unable to allocate standard ocean arrays' ) 370 ! 371 END SUBROUTINE nemo_alloc 372 373 SUBROUTINE nemo_partition( num_pes ) 374 !!---------------------------------------------------------------------- 375 !! *** ROUTINE nemo_partition *** 376 !! 377 !! ** Purpose : 378 !! 379 !! ** Method : 380 !!---------------------------------------------------------------------- 381 USE par_oce 382 INTEGER, INTENT(in) :: num_pes ! The number of MPI processes we have 383 ! Local variables 384 INTEGER, PARAMETER :: nfactmax = 20 385 INTEGER :: nfact ! The no. of factors returned 386 INTEGER :: ierr ! Error flag 387 INTEGER :: i 388 INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are closest in value 389 INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors 390 !!---------------------------------------------------------------------- 391 392 ierr = 0 393 394 CALL factorise( ifact, nfactmax, nfact, num_pes, ierr ) 395 396 IF( nfact <= 1 ) THEN 397 WRITE (numout, *) 'WARNING: factorisation of number of PEs failed' 398 WRITE (numout, *) ' : using grid of ',num_pes,' x 1' 399 jpnj = 1 400 jpni = num_pes 401 ELSE 402 ! Search through factors for the pair that are closest in value 403 mindiff = 1000000 404 imin = 1 405 DO i=1,nfact-1,2 406 idiff = ABS(ifact(i) - ifact(i+1)) 407 IF(idiff < mindiff)THEN 408 mindiff = idiff 409 imin = i 410 END IF 411 END DO 412 jpnj = ifact(imin) 413 jpni = ifact(imin + 1) 414 ENDIF 415 jpnij = jpni*jpnj 416 417 WRITE(*,*) 'ARPDBG: jpni = ',jpni,'jpnj = ',jpnj,'jpnij = ',jpnij 418 ! 419 END SUBROUTINE nemo_partition 420 421 422 SUBROUTINE factorise( ifax, maxfax, nfax, n, ierr ) 423 !!---------------------------------------------------------------------- 424 !! *** ROUTINE factorise *** 425 !! 426 !! ** Purpose : return the prime factors of n. 427 !! nfax factors are returned in array ifax which is of 428 !! maximum dimension maxfax. 429 !! ** Method : 430 !!---------------------------------------------------------------------- 431 INTEGER, INTENT(in) :: n, maxfax 432 INTEGER, INTENT(Out) :: ierr, nfax 433 INTEGER, INTENT(out) :: ifax(maxfax) 434 ! Local variables. 435 INTEGER :: i, ifac, l, nu 436 INTEGER, PARAMETER :: ntest = 14 437 INTEGER :: lfax(ntest) 438 439 ! lfax contains the set of allowed factors. 440 data (lfax(i),i=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256, & 441 & 128, 64, 32, 16, 8, 4, 2 / 442 !!---------------------------------------------------------------------- 443 444 ! Clear the error flag and initialise output vars 445 ierr = 0 446 ifax = 1 447 nfax = 0 448 449 ! Find the factors of n. 450 IF( n == 1 ) GOTO 20 451 452 ! nu holds the unfactorised part of the number. 453 ! nfax holds the number of factors found. 454 ! l points to the allowed factor list. 455 ! ifac holds the current factor. 456 457 nu = n 458 nfax = 0 459 460 DO l = ntest, 1, -1 461 ! 462 ifac = lfax(l) 463 IF(ifac > nu)CYCLE 464 465 ! Test whether the factor will divide. 466 467 IF( MOD(nu,ifac) == 0 ) THEN 468 ! 469 nfax = nfax+1 ! Add the factor to the list 470 IF( nfax > maxfax ) THEN 471 ierr = 6 472 write (*,*) 'FACTOR: insufficient space in factor array ',nfax 473 return 474 ENDIF 475 ifax(nfax) = ifac 476 ! Store the other factor that goes with this one 477 nfax = nfax + 1 478 ifax(nfax) = nu / ifac 479 !WRITE (*,*) 'ARPDBG, factors ',nfax-1,' & ',nfax,' are ', & 480 ! ifax(nfax-1),' and ',ifax(nfax) 481 ENDIF 482 ! 483 END DO 484 485 20 CONTINUE ! Label 20 is the exit point from the factor search loop. 486 ! 487 RETURN 488 ! 489 END SUBROUTINE factorise 309 490 !!====================================================================== 310 491 END MODULE nemogcm
Note: See TracChangeset
for help on using the changeset viewer.