- Timestamp:
- 2011-03-04T12:04:28+01:00 (13 years ago)
- Location:
- branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO
- Files:
-
- 4 edited
-
OPA_SRC/DOM/dommsk.F90 (modified) (7 diffs)
-
OPA_SRC/OBS/obs_readmdt.F90 (modified) (10 diffs)
-
OPA_SRC/nemogcm.F90 (modified) (4 diffs)
-
TOP_SRC/TRP/trcsbc.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90
r2633 r2651 34 34 PRIVATE 35 35 36 PUBLIC dom_msk ! routine called by inidom.F9037 PUBLIC dom_msk_alloc ! routine called by nemogcm.F9036 PUBLIC dom_msk ! routine called by inidom.F90 37 PUBLIC dom_msk_alloc ! routine called by nemogcm.F90 38 38 39 39 ! !!* Namelist namlbc : lateral boundary condition * … … 51 51 CONTAINS 52 52 53 FUNCTION dom_msk_alloc()53 INTEGER FUNCTION dom_msk_alloc() 54 54 !!--------------------------------------------------------------------- 55 !! *** ROUTINEdom_msk_alloc ***55 !! *** FUNCTION dom_msk_alloc *** 56 56 !!--------------------------------------------------------------------- 57 INTEGER :: dom_msk_alloc58 59 57 dom_msk_alloc = 0 60 61 58 #if defined key_noslip_accurate 62 ALLOCATE(icoord(jpi*jpj*jpk,3), S tat=dom_msk_alloc)59 ALLOCATE(icoord(jpi*jpj*jpk,3), STAT=dom_msk_alloc) 63 60 #endif 64 65 IF(dom_msk_alloc /= 0)THEN 66 CALL ctl_warn('dom_msk_alloc: failed to allocate icoord array.') 67 END IF 68 61 IF( dom_msk_alloc /= 0 ) CALL ctl_warn('dom_msk_alloc: failed to allocate icoord array') 62 ! 69 63 END FUNCTION dom_msk_alloc 70 64 … … 131 125 !! tmask_i : interior ocean mask 132 126 !!---------------------------------------------------------------------- 133 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released, iwrk_in_use, iwrk_not_released134 USE wrk_nemo, ONLY: zwf => wrk_2d_1135 USE wrk_nemo, ONLY: imsk => iwrk_2d_1127 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released, iwrk_in_use, iwrk_not_released 128 USE wrk_nemo, ONLY: zwf => wrk_2d_1 129 USE wrk_nemo, ONLY: imsk => iwrk_2d_1 136 130 !! 137 131 INTEGER :: ji, jj, jk ! dummy loop indices … … 142 136 !!--------------------------------------------------------------------- 143 137 144 IF( wrk_in_use(2,1) .OR. iwrk_in_use(2,1) )THEN 145 CALL ctl_stop('dom_msk: ERROR: requested workspace arrays unavailable.') 146 RETURN 147 END IF 138 IF( wrk_in_use(2, 1) .OR. iwrk_in_use(2, 1) )THEN 139 CALL ctl_stop('dom_msk: ERROR: requested workspace arrays unavailable') ; RETURN 140 ENDIF 148 141 149 142 REWIND( numnam ) ! Namelist namlbc : lateral momentum boundary condition … … 443 436 ENDIF 444 437 ! 445 IF( wrk_not_released(2,1) .OR. iwrk_not_released(2,1) )THEN 446 CALL ctl_stop('dom_msk: ERROR: failed to release workspace arrays.') 447 END IF 438 IF( wrk_not_released(2, 1) .OR. & 439 iwrk_not_released(2, 1) ) CALL ctl_stop('dom_msk: ERROR: failed to release workspace arrays') 448 440 ! 449 441 END SUBROUTINE dom_msk … … 464 456 !! ** Action : 465 457 !!---------------------------------------------------------------------- 466 INTEGER :: ji, jj, jk, jl ! dummy loop indices458 INTEGER :: ji, jj, jk, jl ! dummy loop indices 467 459 INTEGER :: ine, inw, ins, inn, itest, ierror, iind, ijnd 468 460 REAL(wp) :: zaa 469 461 !!--------------------------------------------------------------------- 470 471 462 472 463 IF(lwp)WRITE(numout,*) 473 464 IF(lwp)WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition' 474 465 IF(lwp)WRITE(numout,*) '~~~~~~~~~~~ using Schchepetkin and O Brian scheme' 475 IF( lk_mpp ) CALL ctl_stop( ' mpp version is not yet implemented' )466 IF( lk_mpp ) CALL ctl_stop( ' mpp version is not yet implemented' ) 476 467 477 468 ! mask for second order calculation of vorticity … … 628 619 CALL ctl_stop( 'We stop...' ) 629 620 ENDIF 630 621 ! 631 622 END SUBROUTINE dom_msk_nsa 632 623 -
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OPA_SRC/OBS/obs_readmdt.F90
r2638 r2651 4 4 !! Observation diagnostics: Read the MDT for SLA data (skeleton for now) 5 5 !!====================================================================== 6 7 !!---------------------------------------------------------------------- 8 !! obs_rea_mdt : Driver for reading MDT 9 !!---------------------------------------------------------------------- 10 11 !! * Modules used 12 USE par_kind, ONLY : & ! Precision variables 13 & wp, & 14 & dp, & 15 & sp 16 USE par_oce, ONLY : & ! Domain parameters 17 & jpi, & 18 & jpj, & 19 & jpim1 20 USE in_out_manager, ONLY : & ! I/O manager 21 & lwp, & 22 & numout 23 USE obs_surf_def ! Surface observation definitions 24 USE dom_oce, ONLY : & ! Domain variables 25 & tmask, & 26 & tmask_i, & 27 & e1t, & 28 & e2t, & 29 & gphit, & 30 & glamt 31 USE obs_const, ONLY : & 32 & obfillflt ! Fillvalue 33 USE oce, ONLY : & ! Model variables 34 & sshn 35 USE obs_inter_sup ! Interpolation support routines 36 USE obs_inter_h2d ! 2D interpolation 37 USE obs_utils ! Various observation tools 38 USE lib_mpp, only: & ! MPP routines 39 & lk_mpp, & 40 & mpp_sum 41 USE iom_nf90 42 USE netcdf ! NetCDF library 43 USE lib_mpp ! For ctl_warn/stop 6 !! History : ! 2007-03 (K. Mogensen) Initial skeleton version 7 !! ! 2007-04 (E. Remy) migration and improvement from OPAVAR 8 !!---------------------------------------------------------------------- 9 10 !!---------------------------------------------------------------------- 11 !! obs_rea_mdt : Driver for reading MDT 12 !! obs_offset_mdt : Remove the offset between the model MDT and the used one 13 !!---------------------------------------------------------------------- 14 USE par_kind ! Precision variables 15 USE par_oce ! Domain parameters 16 USE in_out_manager ! I/O manager 17 USE obs_surf_def ! Surface observation definitions 18 USE obs_inter_sup ! Interpolation support routines 19 USE obs_inter_h2d ! 2D interpolation 20 USE obs_utils ! Various observation tools 21 USE iom_nf90 ! IOM NetCDF 22 USE netcdf ! NetCDF library 23 USE lib_mpp ! MPP library 24 USE dom_oce, ONLY : & ! Domain variables 25 & tmask, tmask_i, e1t, e2t, gphit, glamt 26 USE obs_const, ONLY : obfillflt ! Fillvalue 27 USE oce , ONLY : sshn ! Model variables 44 28 45 29 IMPLICIT NONE 46 47 !! * Routine accessibility48 30 PRIVATE 49 50 INTEGER, PUBLIC :: nmsshc = 1 ! MDT correction scheme51 REAL(wp), PUBLIC :: mdtcorr = 1.61 ! User specified MDT correction52 REAL(wp), PUBLIC :: mdtcutoff = 65.0 ! MDT cutoff for computed correction 53 PUBLIC obs_rea_mdt ! Read the MDT54 PUBLIC obs_offset_mdt ! Remove the offset between the model MDT and the55 ! used one31 32 PUBLIC obs_rea_mdt ! called by ? 33 PUBLIC obs_offset_mdt ! called by ? 34 35 INTEGER , PUBLIC :: nmsshc = 1 ! MDT correction scheme 36 REAL(wp), PUBLIC :: mdtcorr = 1.61_wp ! User specified MDT correction 37 REAL(wp), PUBLIC :: mdtcutoff = 65.0_wp ! MDT cutoff for computed correction 56 38 57 39 !!---------------------------------------------------------------------- 58 40 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 59 41 !! $Id$ 60 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 61 !!---------------------------------------------------------------------- 62 42 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 43 !!---------------------------------------------------------------------- 63 44 CONTAINS 64 45 … … 73 54 !! 74 55 !! ** Action : 75 !! 76 !! References : 77 !! 78 !! History : 79 !! ! : 2007-03 (K. Mogensen) Initial skeleton version 80 !! ! : 2007-04 (E. Remy) migration and improvement from OPAVAR 81 !!---------------------------------------------------------------------- 82 !! * Modules used 56 !!---------------------------------------------------------------------- 83 57 USE iom 84 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 85 USE wrk_nemo, ONLY: z_mdt => wrk_2d_1, & ! Array to store the MDT values 86 mdtmask => wrk_2d_2 ! Array to store the mask for the MDT 87 !! 88 !! * Arguments 89 INTEGER, INTENT(IN) :: kslano ! Number of SLA Products 90 TYPE(obs_surf), DIMENSION(kslano), INTENT(INOUT) :: & 91 & sladata ! SLA data 92 INTEGER, INTENT(IN) :: k2dint 93 94 !! * Local declarations 95 96 CHARACTER(LEN=12), PARAMETER :: & 97 & cpname = 'obs_rea_mdt' 98 CHARACTER(LEN=20), PARAMETER :: & 99 & mdtname = 'slaReferenceLevel.nc' 100 101 INTEGER :: jslano ! Data set loop variable 102 INTEGER :: jobs ! Obs loop variable 103 INTEGER :: jpimdt ! Number of grid point in latitude for the MDT 104 INTEGER :: jpjmdt ! Number of grid point in longitude for the MDT 105 INTEGER :: iico ! Grid point indicies 106 INTEGER :: ijco 107 INTEGER :: i_nx_id ! Index to read the NetCDF file 108 INTEGER :: i_ny_id ! 109 INTEGER :: i_file_id ! 110 INTEGER :: i_var_id 111 INTEGER :: i_stat 112 113 REAL(wp), DIMENSION(1) :: & 114 & zext, & 115 & zobsmask 116 REAL(wp), DIMENSION(2,2,1) :: & 117 & zweig 118 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: & 119 & zmask, & 120 & zmdtl, & 121 & zglam, & 122 & zgphi 58 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 59 USE wrk_nemo, ONLY: z_mdt => wrk_2d_1 ! Array to store the MDT values 60 USE wrk_nemo, ONLY: mdtmask => wrk_2d_2 ! Array to store the mask for the MDT 61 ! 62 INTEGER , INTENT(IN) :: kslano ! Number of SLA Products 63 TYPE(obs_surf), DIMENSION(kslano), INTENT(inout) :: sladata ! SLA data 64 INTEGER , INTENT(in) :: k2dint ! ? 65 ! 66 CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_rea_mdt' 67 CHARACTER(LEN=20), PARAMETER :: mdtname = 'slaReferenceLevel.nc' 68 69 INTEGER :: jslano ! Data set loop variable 70 INTEGER :: jobs ! Obs loop variable 71 INTEGER :: jpimdt, jpjmdt ! Number of grid point in lat/lon for the MDT 72 INTEGER :: iico, ijco ! Grid point indicies 73 INTEGER :: i_nx_id, i_ny_id, i_file_id, i_var_id, i_stat 74 INTEGER :: nummdt 75 ! 76 REAL(wp), DIMENSION(1) :: zext, zobsmask 77 REAL(wp), DIMENSION(2,2,1) :: zweig 78 ! 79 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zmask, zmdtl, zglam, zgphi 80 INTEGER , DIMENSION(:,:,:), ALLOCATABLE :: igrdi, igrdj 123 81 124 REAL(wp) :: zlam 125 REAL(wp) :: zphi 126 REAL(wp) :: zfill 127 REAL(sp) :: zinfill 128 INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: & 129 & igrdi, & 130 & igrdj 131 INTEGER :: nummdt 132 !!---------------------------------------------------------------------- 133 134 IF(wrk_in_use(2, 1,2))THEN 135 CALL ctl_stop('obs_rea_mdt : requested workspace array unavailable.') 136 RETURN 137 END IF 82 REAL(wp) :: zlam, zphi, zfill, zinfill ! local scalar 83 !!---------------------------------------------------------------------- 84 85 IF( wrk_in_use(2, 1,2) ) THEN 86 CALL ctl_stop('obs_rea_mdt : requested workspace array unavailable') ; RETURN 87 ENDIF 138 88 139 89 IF(lwp)WRITE(numout,*) 140 IF(lwp)WRITE(numout,*) ' obs_rea_mdt : '90 IF(lwp)WRITE(numout,*) ' obs_rea_mdt : Read MDT for referencing altimeter anomalies' 141 91 IF(lwp)WRITE(numout,*) ' ------------- ' 142 IF(lwp)WRITE(numout,*) ' Read MDT for referencing altimeter', & 143 & ' anomalies' 144 145 ! Open the file 146 147 CALL iom_open( mdtname, nummdt ) 148 149 ! Get the MDT data 150 151 CALL iom_get( nummdt, jpdom_data, 'sossheig', z_mdt(:,:), 1 ) 152 153 ! Close the file 154 155 CALL iom_close(nummdt) 92 93 CALL iom_open( mdtname, nummdt ) ! Open the file 94 ! ! Get the MDT data 95 CALL iom_get ( nummdt, jpdom_data, 'sossheig', z_mdt(:,:), 1 ) 96 CALL iom_close(nummdt) ! Close the file 156 97 157 98 ! Read in the fill value … … 163 104 i_stat = nf90_close( nummdt ) 164 105 165 ! setup mask based on tmask and MDT mask 166 ! set mask to 0 where the MDT is set to fillvalue 167 168 WHERE(z_mdt(:,:) /= zfill) 169 mdtmask(:,:)=tmask(:,:,1) 170 ELSEWHERE 171 mdtmask(:,:)=0 106 ! setup mask based on tmask and MDT mask 107 ! set mask to 0 where the MDT is set to fillvalue 108 WHERE(z_mdt(:,:) /= zfill) ; mdtmask(:,:) = tmask(:,:,1) 109 ELSE WHERE ; mdtmask(:,:) = 0 172 110 END WHERE 173 111 174 112 ! Remove the offset between the MDT used with the sla and the model MDT 175 176 IF ( nmsshc == 1 .OR. nmsshc == 2 ) CALL obs_offset_mdt( z_mdt, zfill ) 113 IF( nmsshc == 1 .OR. nmsshc == 2 ) CALL obs_offset_mdt( z_mdt, zfill ) 177 114 178 115 ! Intepolate the MDT already on the model grid at the observation point 179 116 180 117 DO jslano = 1, kslano 181 182 118 ALLOCATE( & 183 119 & igrdi(2,2,sladata(jslano)%nsurf), & … … 202 138 END DO 203 139 204 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, & 205 & igrdi, igrdj, glamt, zglam ) 206 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, & 207 & igrdi, igrdj, gphit, zgphi ) 208 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, & 209 & igrdi, igrdj, mdtmask, zmask ) 210 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, & 211 & igrdi, igrdj, z_mdt, zmdtl ) 140 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, igrdi, igrdj, glamt , zglam ) 141 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, igrdi, igrdj, gphit , zgphi ) 142 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, igrdi, igrdj, mdtmask, zmask ) 143 CALL obs_int_comm_2d( 2, 2, sladata(jslano)%nsurf, igrdi, igrdj, z_mdt , zmdtl ) 212 144 213 145 DO jobs = 1, sladata(jslano)%nsurf … … 220 152 & zmask(:,:,jobs), zweig, zobsmask ) 221 153 222 CALL obs_int_h2d( 1, 1, & 223 & zweig, zmdtl(:,:,jobs), zext ) 154 CALL obs_int_h2d( 1, 1, zweig, zmdtl(:,:,jobs), zext ) 224 155 225 156 sladata(jslano)%rext(jobs,2) = zext(1) 226 157 227 158 ! mark any masked data with a QC flag 228 IF ( zobsmask(1) == 0 )sladata(jslano)%nqc(jobs) = 11159 IF( zobsmask(1) == 0 ) sladata(jslano)%nqc(jobs) = 11 229 160 230 161 END DO … … 241 172 END DO 242 173 243 IF(wrk_not_released(2, 1,2))THEN 244 CALL ctl_stop('obs_rea_mdt : failed to release workspace arrays.') 245 END IF 246 174 IF( wrk_not_released(2, 1,2) ) CALL ctl_stop('obs_rea_mdt: failed to release workspace arrays') 175 ! 247 176 END SUBROUTINE obs_rea_mdt 248 177 178 249 179 SUBROUTINE obs_offset_mdt( mdt, zfill ) 250 251 180 !!--------------------------------------------------------------------- 252 181 !! … … 260 189 !! 261 190 !! ** Action : 262 !! 263 !! References : 264 !! 265 !! History : 266 !! ! : 2007-04 (E. Remy) migration from OPAVAR 267 !!---------------------------------------------------------------------- 268 !! * Modules used 269 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 270 USE wrk_nemo, ONLY: zpromsk => wrk_2d_3 271 !! 272 !! * Arguments 273 REAL(wp), DIMENSION(jpi,jpj), INTENT(INOUT) :: & 274 & mdt ! MDT used on the model grid 275 REAL(wp), INTENT(IN) :: zfill 276 277 !! * Local declarations 278 REAL(wp) :: zdxdy 279 REAL(wp) :: zarea 280 REAL(wp) :: zeta1 281 REAL(wp) :: zeta2 282 REAL(wp) :: zcorr_mdt 283 REAL(wp) :: zcorr_bcketa 284 REAL(wp) :: zcorr 285 INTEGER :: jj 286 INTEGER :: ji 287 CHARACTER(LEN=14), PARAMETER :: & 288 & cpname = 'obs_offset_mdt' 289 !!---------------------------------------------------------------------- 290 291 IF(wrk_in_use(2, 3))THEN 292 CALL ctl_stop('obs_offset_mdt : requested workspace array unavailable.') 293 RETURN 294 END IF 191 !!---------------------------------------------------------------------- 192 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 193 USE wrk_nemo, ONLY: zpromsk => wrk_2d_3 194 ! 195 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: mdt ! MDT used on the model grid 196 REAL(wp) , INTENT(in ) :: zfill 197 ! 198 INTEGER :: ji, jj 199 REAL(wp) :: zdxdy, zarea, zeta1, zeta2, zcorr_mdt, zcorr_bcketa, zcorr ! local scalar 200 CHARACTER(LEN=14), PARAMETER :: cpname = 'obs_offset_mdt' 201 !!---------------------------------------------------------------------- 202 203 IF( wrk_in_use(2, 3) ) THEN 204 CALL ctl_stop('obs_offset_mdt: requested workspace array unavailable') ; RETURN 205 ENDIF 295 206 296 207 ! Initialize the local mask, for domain projection … … 322 233 END DO 323 234 324 IF( lk_mpp) CALL mpp_sum( zeta1 )325 IF( lk_mpp) CALL mpp_sum( zeta2 )326 IF( lk_mpp) CALL mpp_sum( zarea )235 IF( lk_mpp) CALL mpp_sum( zeta1 ) 236 IF( lk_mpp) CALL mpp_sum( zeta2 ) 237 IF( lk_mpp) CALL mpp_sum( zarea ) 327 238 328 zcorr_mdt = zeta1 / zarea329 zcorr_bcketa 239 zcorr_mdt = zeta1 / zarea 240 zcorr_bcketa = zeta2 / zarea 330 241 331 242 ! Define correction term … … 335 246 ! Correct spatial mean of the MSSH 336 247 337 IF ( nmsshc == 1 )mdt(:,:) = mdt(:,:) - zcorr248 IF( nmsshc == 1 ) mdt(:,:) = mdt(:,:) - zcorr 338 249 339 250 ! User defined value : 1.6 m for the Rio MDT compared to ORCA2 MDT 340 251 341 IF ( nmsshc == 2 )mdt(:,:) = mdt(:,:) - mdtcorr252 IF( nmsshc == 2 ) mdt(:,:) = mdt(:,:) - mdtcorr 342 253 343 254 IF(lwp) THEN … … 348 259 WRITE(numout,*) ' zcorr = ', zcorr 349 260 WRITE(numout,*) ' nmsshc = ', nmsshc 350 WRITE(numout,*)351 261 ENDIF 352 262 353 IF ( nmsshc == 0 ) WRITE(numout,*) & 354 & ' MSSH correction is not applied' 355 IF ( nmsshc == 1 ) WRITE(numout,*) & 356 & ' MSSH correction is applied' 357 IF ( nmsshc == 2 ) WRITE(numout,*) & 358 & ' User defined MSSH correction' 359 360 361 IF(wrk_not_released(2, 3))THEN 362 CALL ctl_stop('obs_offset_mdt : failed to release workspace array.') 363 END IF 364 263 IF ( nmsshc == 0 ) WRITE(numout,*) ' MSSH correction is not applied' 264 IF ( nmsshc == 1 ) WRITE(numout,*) ' MSSH correction is applied' 265 IF ( nmsshc == 2 ) WRITE(numout,*) ' User defined MSSH correction' 266 267 IF( wrk_not_released(2, 3) ) CALL ctl_stop('obs_offset_mdt: failed to release workspace array') 268 ! 365 269 END SUBROUTINE obs_offset_mdt 366 270 271 !!====================================================================== 367 272 END MODULE obs_readmdt -
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90
r2647 r2651 502 502 !! ** Method : 503 503 !!---------------------------------------------------------------------- 504 USE par_oce505 504 INTEGER, INTENT(in) :: num_pes ! The number of MPI processes we have 506 ! Local variables505 ! 507 506 INTEGER, PARAMETER :: nfactmax = 20 508 507 INTEGER :: nfact ! The no. of factors returned 509 508 INTEGER :: ierr ! Error flag 510 INTEGER :: i509 INTEGER :: ji 511 510 INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are closest in value 512 511 INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors … … 526 525 mindiff = 1000000 527 526 imin = 1 528 DO i=1,nfact-1,2529 idiff = ABS( ifact(i) - ifact(i+1))530 IF( idiff < mindiff)THEN527 DO ji = 1, nfact-1, 2 528 idiff = ABS( ifact(ji) - ifact(ji+1) ) 529 IF( idiff < mindiff ) THEN 531 530 mindiff = idiff 532 imin = i533 END 531 imin = ji 532 ENDIF 534 533 END DO 535 534 jpnj = ifact(imin) … … 543 542 544 543 545 SUBROUTINE factorise( ifax, maxfax, nfax, n, ierr )544 SUBROUTINE factorise( kfax, kmaxfax, knfax, kn, kerr ) 546 545 !!---------------------------------------------------------------------- 547 546 !! *** ROUTINE factorise *** 548 547 !! 549 548 !! ** Purpose : return the prime factors of n. 550 !! nfax factors are returned in array ifax which is of551 !! maximum dimension maxfax.549 !! knfax factors are returned in array kfax which is of 550 !! maximum dimension kmaxfax. 552 551 !! ** Method : 553 552 !!---------------------------------------------------------------------- 554 INTEGER , INTENT(in) :: n,maxfax555 INTEGER , INTENT(Out) :: ierr,nfax556 INTEGER, INTENT(out) :: ifax(maxfax)557 ! Local variables.558 INTEGER :: i , ifac, l,nu553 INTEGER , INTENT(in ) :: kn, kmaxfax 554 INTEGER , INTENT( out) :: kerr, knfax 555 INTEGER, DIMENSION(kmaxfax), INTENT( out) :: kfax 556 ! 557 INTEGER :: ifac, jl, inu 559 558 INTEGER, PARAMETER :: ntest = 14 560 INTEGER :: lfax(ntest)559 INTEGER :: ilfax(ntest) 561 560 562 561 ! lfax contains the set of allowed factors. 563 data ( lfax(i),i=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256, &564 & 128, 64, 32, 16, 8, 4, 2 /562 data (ilfax(jl),jl=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256, & 563 & 128, 64, 32, 16, 8, 4, 2 / 565 564 !!---------------------------------------------------------------------- 566 565 567 566 ! Clear the error flag and initialise output vars 568 ierr = 0569 ifax = 1570 nfax = 0567 kerr = 0 568 kfax = 1 569 knfax = 0 571 570 572 571 ! Find the factors of n. 573 IF( n == 1 )GOTO 20572 IF( kn == 1 ) GOTO 20 574 573 575 574 ! nu holds the unfactorised part of the number. 576 ! nfax holds the number of factors found.575 ! knfax holds the number of factors found. 577 576 ! l points to the allowed factor list. 578 577 ! ifac holds the current factor. 579 578 580 nu =n581 nfax = 0582 583 DO l = ntest, 1, -1579 inu = kn 580 knfax = 0 581 582 DO jl = ntest, 1, -1 584 583 ! 585 ifac = lfax(l)586 IF( ifac > nu)CYCLE584 ifac = ilfax(jl) 585 IF( ifac > inu ) CYCLE 587 586 588 587 ! Test whether the factor will divide. 589 588 590 IF( MOD( nu,ifac) == 0 ) THEN589 IF( MOD(inu,ifac) == 0 ) THEN 591 590 ! 592 nfax = nfax+1 ! Add the factor to the list593 IF( nfax >maxfax ) THEN594 ierr = 6595 write (*,*) 'FACTOR: insufficient space in factor array ', nfax591 knfax = knfax + 1 ! Add the factor to the list 592 IF( knfax > kmaxfax ) THEN 593 kerr = 6 594 write (*,*) 'FACTOR: insufficient space in factor array ', knfax 596 595 return 597 596 ENDIF 598 ifax(nfax) = ifac597 kfax(knfax) = ifac 599 598 ! Store the other factor that goes with this one 600 nfax = nfax + 1 601 ifax(nfax) = nu / ifac 602 !WRITE (*,*) 'ARPDBG, factors ',nfax-1,' & ',nfax,' are ', & 603 ! ifax(nfax-1),' and ',ifax(nfax) 599 knfax = knfax + 1 600 kfax(knfax) = inu / ifac 601 !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 604 602 ENDIF 605 603 ! … … 608 606 20 CONTINUE ! Label 20 is the exit point from the factor search loop. 609 607 ! 610 RETURN611 !612 608 END SUBROUTINE factorise 613 609 -
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbc.F90
r2643 r2651 57 57 !! 58 58 !!---------------------------------------------------------------------- 59 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released60 USE wrk_nemo, zemps => wrk_2d_161 USE wrk_nemo, ztrtrd => wrk_3d_159 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 60 USE wrk_nemo, ONLY: zemps => wrk_2d_1 61 USE wrk_nemo, ONLY: ztrtrd => wrk_3d_1 62 62 ! 63 63 INTEGER, INTENT( in ) :: kt ! ocean time-step index
Note: See TracChangeset
for help on using the changeset viewer.
