Changeset 3294 for trunk/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90
- Timestamp:
- 2012-01-28T17:44:18+01:00 (12 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90
r2715 r3294 10 10 !! 3.3 ! 2010-09 (E.O'Dea) modifications for Shelf configurations 11 11 !! 3.3 ! 2010-09 (D.Storkey) add ice boundary conditions 12 !! 3.4 ! 2011 (D. Storkey) rewrite in preparation for OBC-BDY merge 12 13 !!---------------------------------------------------------------------- 13 14 #if defined key_bdy 14 15 !!---------------------------------------------------------------------- 15 !! 'key_bdy' Unstructured Open Boundary Conditions 16 !!---------------------------------------------------------------------- 17 !! bdy_dta_frs : read u, v, t, s data along open boundaries 18 !! bdy_dta_fla : read depth-mean velocities and elevation along open boundaries 19 !!---------------------------------------------------------------------- 16 !! 'key_bdy' Open Boundary Conditions 17 !!---------------------------------------------------------------------- 18 !! bdy_dta : read external data along open boundaries from file 19 !! bdy_dta_init : initialise arrays etc for reading of external data 20 !!---------------------------------------------------------------------- 21 USE wrk_nemo ! Memory Allocation 22 USE timing ! Timing 20 23 USE oce ! ocean dynamics and tracers 21 24 USE dom_oce ! ocean space and time domain 22 25 USE phycst ! physical constants 23 USE bdy_oce ! ocean open boundary conditions 26 USE bdy_oce ! ocean open boundary conditions 24 27 USE bdytides ! tidal forcing at boundaries 25 USE iom26 USE io ipsl28 USE fldread ! read input fields 29 USE iom ! IOM library 27 30 USE in_out_manager ! I/O logical units 28 31 #if defined key_lim2 … … 33 36 PRIVATE 34 37 35 PUBLIC bdy_dta_frs ! routines called by step.F90 36 PUBLIC bdy_dta_fla 37 PUBLIC bdy_dta_alloc ! routine called by bdy_init.F90 38 39 INTEGER :: numbdyt, numbdyu, numbdyv ! logical units for T-, U-, & V-points data file, resp. 40 INTEGER :: ntimes_bdy ! exact number of time dumps in data files 41 INTEGER :: nbdy_b, nbdy_a ! record of bdy data file for before and after time step 42 INTEGER :: numbdyt_bt, numbdyu_bt, numbdyv_bt ! logical unit for T-, U- & V-points data file, resp. 43 INTEGER :: ntimes_bdy_bt ! exact number of time dumps in data files 44 INTEGER :: nbdy_b_bt, nbdy_a_bt ! record of bdy data file for before and after time step 45 46 INTEGER, DIMENSION (jpbtime) :: istep, istep_bt ! time array in seconds in each data file 47 48 REAL(wp) :: zoffset ! time offset between time origin in file & start time of model run 49 50 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tbdydta, sbdydta ! time interpolated values of T and S bdy data 51 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ubdydta, vbdydta ! time interpolated values of U and V bdy data 52 REAL(wp), DIMENSION(jpbdim,2) :: ubtbdydta, vbtbdydta ! Arrays used for time interpolation of bdy data 53 REAL(wp), DIMENSION(jpbdim,2) :: sshbdydta ! bdy data of ssh 54 55 #if defined key_lim2 56 REAL(wp), DIMENSION(jpbdim,2) :: frld_bdydta ! } 57 REAL(wp), DIMENSION(jpbdim,2) :: hicif_bdydta ! } Arrays used for time interp. of ice bdy data 58 REAL(wp), DIMENSION(jpbdim,2) :: hsnif_bdydta ! } 59 #endif 60 38 PUBLIC bdy_dta ! routine called by step.F90 and dynspg_ts.F90 39 PUBLIC bdy_dta_init ! routine called by nemogcm.F90 40 41 INTEGER, ALLOCATABLE, DIMENSION(:) :: nb_bdy_fld ! Number of fields to update for each boundary set. 42 INTEGER :: nb_bdy_fld_sum ! Total number of fields to update for all boundary sets. 43 44 LOGICAL, DIMENSION(jp_bdy) :: ln_full_vel_array ! =T => full velocities in 3D boundary conditions 45 ! =F => baroclinic velocities in 3D boundary conditions 46 47 TYPE(FLD), PUBLIC, ALLOCATABLE, DIMENSION(:), TARGET :: bf ! structure of input fields (file informations, fields read) 48 49 TYPE(MAP_POINTER), ALLOCATABLE, DIMENSION(:) :: nbmap_ptr ! array of pointers to nbmap 50 51 # include "domzgr_substitute.h90" 61 52 !!---------------------------------------------------------------------- 62 53 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 66 57 CONTAINS 67 58 68 FUNCTION bdy_dta_alloc() 69 !!---------------------------------------------------------------------- 70 USE lib_mpp, ONLY: ctl_warn, mpp_sum 71 ! 72 INTEGER :: bdy_dta_alloc 73 !!---------------------------------------------------------------------- 74 ! 75 ALLOCATE(tbdydta(jpbdim,jpk,2), sbdydta(jpbdim,jpk,2), & 76 ubdydta(jpbdim,jpk,2), vbdydta(jpbdim,jpk,2), Stat=bdy_dta_alloc) 77 78 IF( lk_mpp ) CALL mpp_sum ( bdy_dta_alloc ) 79 IF(bdy_dta_alloc /= 0) CALL ctl_warn('bdy_dta_alloc: failed to allocate arrays') 80 81 END FUNCTION bdy_dta_alloc 82 83 84 SUBROUTINE bdy_dta_frs( kt ) 59 SUBROUTINE bdy_dta( kt, jit, time_offset ) 85 60 !!---------------------------------------------------------------------- 86 !! *** SUBROUTINE bdy_dta _frs***61 !! *** SUBROUTINE bdy_dta *** 87 62 !! 88 !! ** Purpose : Read unstructured boundary data for FRS condition. 89 !! 90 !! ** Method : At the first timestep, read in boundary data for two 91 !! times from the file and time-interpolate. At other 92 !! timesteps, check to see if we need another time from 93 !! the file. If so read it in. Time interpolate. 63 !! ** Purpose : Update external data for open boundary conditions 64 !! 65 !! ** Method : Use fldread.F90 66 !! 94 67 !!---------------------------------------------------------------------- 95 INTEGER, INTENT( in ) :: kt ! ocean time-step index (for timesplitting option, otherwise zero) 96 !! 97 CHARACTER(LEN=80), DIMENSION(3) :: clfile ! names of input files 98 CHARACTER(LEN=70 ) :: clunits ! units attribute of time coordinate 99 LOGICAL :: lect ! flag for reading 100 INTEGER :: it, ib, ik, igrd ! dummy loop indices 101 INTEGER :: igrd_start, igrd_end ! start and end of loops on igrd 102 INTEGER :: idvar ! netcdf var ID 103 INTEGER :: iman, i15, imois ! Time variables for monthly clim forcing 104 INTEGER :: ntimes_bdyt, ntimes_bdyu, ntimes_bdyv 105 INTEGER :: itimer, totime 106 INTEGER :: ii, ij ! array addresses 107 INTEGER :: ipi, ipj, ipk, inum ! local integers (NetCDF read) 108 INTEGER :: iyear0, imonth0, iday0 109 INTEGER :: ihours0, iminutes0, isec0 110 INTEGER :: iyear, imonth, iday, isecs 111 INTEGER, DIMENSION(jpbtime) :: istept, istepu, istepv ! time arrays from data files 112 REAL(wp) :: dayfrac, zxy, zoffsett 113 REAL(wp) :: zoffsetu, zoffsetv 114 REAL(wp) :: dayjul0, zdayjulini 115 REAL(wp), DIMENSION(jpbtime) :: zstepr ! REAL time array from data files 116 REAL(wp), DIMENSION(jpbdta,1,jpk) :: zdta ! temporary array for data fields 68 !! 69 INTEGER, INTENT( in ) :: kt ! ocean time-step index 70 INTEGER, INTENT( in ), OPTIONAL :: jit ! subcycle time-step index (for timesplitting option) 71 INTEGER, INTENT( in ), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit 72 ! is present then units = subcycle timesteps. 73 ! time_offset = 0 => get data at "now" time level 74 ! time_offset = -1 => get data at "before" time level 75 ! time_offset = +1 => get data at "after" time level 76 ! etc. 77 !! 78 INTEGER :: ib_bdy, jfld, jstart, jend, ib, ii, ij, ik, igrd ! local indices 79 INTEGER, DIMENSION(jpbgrd) :: ilen1 80 INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts 81 !! 117 82 !!--------------------------------------------------------------------------- 118 119 120 IF( ln_dyn_frs .OR. ln_tra_frs & 121 & .OR. ln_ice_frs ) THEN ! If these are both false then this routine does nothing 122 123 ! -------------------- ! 124 ! Initialization ! 125 ! -------------------- ! 126 127 lect = .false. ! If true, read a time record 128 129 ! Some time variables for monthly climatological forcing: 130 ! ******************************************************* 131 132 !!gm here use directely daymod calendar variables 133 134 iman = INT( raamo ) ! Number of months in a year 135 136 i15 = INT( 2*REAL( nday, wp ) / ( REAL( nmonth_len(nmonth), wp ) + 0.5 ) ) 137 ! i15=0 if the current day is in the first half of the month, else i15=1 138 139 imois = nmonth + i15 - 1 ! imois is the first month record 140 IF( imois == 0 ) imois = iman 141 142 ! Time variable for non-climatological forcing: 143 ! ********************************************* 144 itimer = (kt-nit000+1)*rdt ! current time in seconds for interpolation 145 146 147 ! !-------------------! 148 IF( kt == nit000 ) THEN ! First call only ! 149 ! !-------------------! 150 istep(:) = 0 151 nbdy_b = 0 152 nbdy_a = 0 153 154 ! Get time information from bdy data file 155 ! *************************************** 156 157 IF(lwp) WRITE(numout,*) 158 IF(lwp) WRITE(numout,*) 'bdy_dta_frs : Initialize unstructured boundary data' 159 IF(lwp) WRITE(numout,*) '~~~~~~~' 160 161 IF ( nn_dtactl == 0 ) THEN 162 ! 163 IF(lwp) WRITE(numout,*) ' Bdy data are taken from initial conditions' 164 ! 165 ELSEIF (nn_dtactl == 1) THEN 166 ! 167 IF(lwp) WRITE(numout,*) ' Bdy data are read in netcdf files' 168 ! 169 dayfrac = adatrj - REAL( itimer, wp ) / 86400. ! day fraction at time step kt-1 170 dayfrac = dayfrac - INT ( dayfrac ) ! 171 totime = ( nitend - nit000 + 1 ) * rdt ! Total time of the run to verify that all the 172 ! ! necessary time dumps in file are included 173 ! 174 clfile(1) = cn_dta_frs_T 175 clfile(2) = cn_dta_frs_U 176 clfile(3) = cn_dta_frs_V 177 ! 178 ! how many files are we to read in? 179 igrd_start = 1 180 igrd_end = 3 181 IF(.NOT. ln_tra_frs .AND. .NOT. ln_ice_frs) THEN ! No T-grid file. 182 igrd_start = 2 183 ELSEIF ( .NOT. ln_dyn_frs ) THEN ! No U-grid or V-grid file. 184 igrd_end = 1 185 ENDIF 186 187 DO igrd = igrd_start, igrd_end ! loop over T, U & V grid ! 188 ! !---------------------------! 189 CALL iom_open( clfile(igrd), inum ) 190 CALL iom_gettime( inum, zstepr, kntime=ntimes_bdy, cdunits=clunits ) 191 192 SELECT CASE( igrd ) 193 CASE (1) ; numbdyt = inum 194 CASE (2) ; numbdyu = inum 195 CASE (3) ; numbdyv = inum 196 END SELECT 197 198 ! Calculate time offset 199 READ(clunits,7000) iyear0, imonth0, iday0, ihours0, iminutes0, isec0 200 ! Convert time origin in file to julian days 201 isec0 = isec0 + ihours0*60.*60. + iminutes0*60. 202 CALL ymds2ju(iyear0, imonth0, iday0, REAL(isec0, wp), dayjul0) 203 ! Compute model initialization time 204 iyear = ndastp / 10000 205 imonth = ( ndastp - iyear * 10000 ) / 100 206 iday = ndastp - iyear * 10000 - imonth * 100 207 isecs = dayfrac * 86400 208 CALL ymds2ju(iyear, imonth, iday, REAL(isecs, wp) , zdayjulini) 209 ! offset from initialization date: 210 zoffset = (dayjul0-zdayjulini)*86400 211 ! 212 7000 FORMAT('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2) 213 214 !! TO BE DONE... Check consistency between calendar from file 215 !! (available optionally from iom_gettime) and calendar in model 216 !! when calendar in model available outside of IOIPSL. 217 218 IF(lwp) WRITE(numout,*) 'number of times: ',ntimes_bdy 219 IF(lwp) WRITE(numout,*) 'offset: ',zoffset 220 IF(lwp) WRITE(numout,*) 'totime: ',totime 221 IF(lwp) WRITE(numout,*) 'zstepr: ',zstepr(1:ntimes_bdy) 222 223 ! Check that there are not too many times in the file. 224 IF( ntimes_bdy > jpbtime ) THEN 225 WRITE(ctmp1,*) 'Check file: ', clfile(igrd), 'jpbtime= ', jpbtime, ' ntimes_bdy= ', ntimes_bdy 226 CALL ctl_stop( 'Number of time dumps in files exceed jpbtime parameter', ctmp1 ) 227 ENDIF 228 229 ! Check that time array increases: 230 it = 1 231 DO WHILE( zstepr(it+1) > zstepr(it) .AND. it /= ntimes_bdy - 1 ) 232 it = it + 1 233 END DO 234 ! 235 IF( it /= ntimes_bdy-1 .AND. ntimes_bdy > 1 ) THEN 236 WRITE(ctmp1,*) 'Check file: ', clfile(igrd) 237 CALL ctl_stop( 'Time array in unstructured boundary data files', & 238 & 'does not continuously increase.' , ctmp1 ) 239 ENDIF 240 ! 241 ! Check that times in file span model run time: 242 IF( zstepr(1) + zoffset > 0 ) THEN 243 WRITE(ctmp1,*) 'Check file: ', clfile(igrd) 244 CALL ctl_stop( 'First time dump in bdy file is after model initial time', ctmp1 ) 245 END IF 246 IF( zstepr(ntimes_bdy) + zoffset < totime ) THEN 247 WRITE(ctmp1,*) 'Check file: ', clfile(igrd) 248 CALL ctl_stop( 'Last time dump in bdy file is before model final time', ctmp1 ) 249 END IF 250 ! 251 SELECT CASE( igrd ) 252 CASE (1) 253 ntimes_bdyt = ntimes_bdy 254 zoffsett = zoffset 255 istept(:) = INT( zstepr(:) + zoffset ) 256 numbdyt = inum 257 CASE (2) 258 ntimes_bdyu = ntimes_bdy 259 zoffsetu = zoffset 260 istepu(:) = INT( zstepr(:) + zoffset ) 261 numbdyu = inum 262 CASE (3) 263 ntimes_bdyv = ntimes_bdy 264 zoffsetv = zoffset 265 istepv(:) = INT( zstepr(:) + zoffset ) 266 numbdyv = inum 267 END SELECT 268 ! 269 END DO ! end loop over T, U & V grid 270 271 IF (igrd_start == 1 .and. igrd_end == 3) THEN 272 ! Only test differences if we are reading in 3 files 273 ! Verify time consistency between files 274 IF( ntimes_bdyu /= ntimes_bdyt .OR. ntimes_bdyv /= ntimes_bdyt ) THEN 275 CALL ctl_stop( 'Bdy data files must have the same number of time dumps', & 276 & 'Multiple time frequencies not implemented yet' ) 277 ENDIF 278 ntimes_bdy = ntimes_bdyt 279 ! 280 IF( zoffsetu /= zoffsett .OR. zoffsetv /= zoffsett ) THEN 281 CALL ctl_stop( 'Bdy data files must have the same time origin', & 282 & 'Multiple time frequencies not implemented yet' ) 283 ENDIF 284 zoffset = zoffsett 285 ENDIF 286 287 IF( igrd_start == 1 ) THEN ; istep(:) = istept(:) 288 ELSE ; istep(:) = istepu(:) 289 ENDIF 290 291 ! Check number of time dumps: 292 IF( ntimes_bdy == 1 .AND. .NOT. ln_clim ) THEN 293 CALL ctl_stop( 'There is only one time dump in data files', & 294 & 'Choose ln_clim=.true. in namelist for constant bdy forcing.' ) 295 ENDIF 296 297 IF( ln_clim ) THEN 298 IF( ntimes_bdy /= 1 .AND. ntimes_bdy /= 12 ) THEN 299 CALL ctl_stop( 'For climatological boundary forcing (ln_clim=.true.),', & 300 & 'bdy data files must contain 1 or 12 time dumps.' ) 301 ELSEIF( ntimes_bdy == 1 ) THEN 302 IF(lwp) WRITE(numout,*) 303 IF(lwp) WRITE(numout,*) 'We assume constant boundary forcing from bdy data files' 304 ELSEIF( ntimes_bdy == 12 ) THEN 305 IF(lwp) WRITE(numout,*) 306 IF(lwp) WRITE(numout,*) 'We assume monthly (and cyclic) boundary forcing from bdy data files' 307 ENDIF 308 ENDIF 309 310 ! Find index of first record to read (before first model time). 311 it = 1 312 DO WHILE( istep(it+1) <= 0 .AND. it <= ntimes_bdy - 1 ) 313 it = it + 1 314 END DO 315 nbdy_b = it 316 ! 317 IF(lwp) WRITE(numout,*) 'Time offset is ',zoffset 318 IF(lwp) WRITE(numout,*) 'First record to read is ',nbdy_b 319 320 ENDIF ! endif (nn_dtactl == 1) 321 322 323 ! 1.2 Read first record in file if necessary (ie if nn_dtactl == 1) 324 ! ***************************************************************** 325 326 IF( nn_dtactl == 0 ) THEN ! boundary data arrays are filled with initial conditions 327 ! 328 IF (ln_tra_frs) THEN 329 igrd = 1 ! T-points data 330 DO ib = 1, nblen(igrd) 331 ii = nbi(ib,igrd) 332 ij = nbj(ib,igrd) 83 !! 84 IF( nn_timing == 1 ) CALL timing_start('bdy_dta') 85 86 ! Initialise data arrays once for all from initial conditions where required 87 !--------------------------------------------------------------------------- 88 IF( kt .eq. nit000 .and. .not. PRESENT(jit) ) THEN 89 90 ! Calculate depth-mean currents 91 !----------------------------- 92 CALL wrk_alloc(jpi,jpj,pu2d,pv2d) 93 94 pu2d(:,:) = 0.e0 95 pv2d(:,:) = 0.e0 96 97 DO ik = 1, jpkm1 !! Vertically integrated momentum trends 98 pu2d(:,:) = pu2d(:,:) + fse3u(:,:,ik) * umask(:,:,ik) * un(:,:,ik) 99 pv2d(:,:) = pv2d(:,:) + fse3v(:,:,ik) * vmask(:,:,ik) * vn(:,:,ik) 100 END DO 101 pu2d(:,:) = pu2d(:,:) * hur(:,:) 102 pv2d(:,:) = pv2d(:,:) * hvr(:,:) 103 104 DO ib_bdy = 1, nb_bdy 105 106 nblen => idx_bdy(ib_bdy)%nblen 107 nblenrim => idx_bdy(ib_bdy)%nblenrim 108 109 IF( nn_dyn2d(ib_bdy) .gt. 0 .and. nn_dyn2d_dta(ib_bdy) .eq. 0 ) THEN 110 IF( nn_dyn2d(ib_bdy) .eq. jp_frs ) THEN 111 ilen1(:) = nblen(:) 112 ELSE 113 ilen1(:) = nblenrim(:) 114 ENDIF 115 igrd = 1 116 DO ib = 1, ilen1(igrd) 117 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 118 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 119 dta_bdy(ib_bdy)%ssh(ib) = sshn(ii,ij) * tmask(ii,ij,1) 120 END DO 121 igrd = 2 122 DO ib = 1, ilen1(igrd) 123 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 124 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 125 dta_bdy(ib_bdy)%u2d(ib) = pu2d(ii,ij) * umask(ii,ij,1) 126 END DO 127 igrd = 3 128 DO ib = 1, ilen1(igrd) 129 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 130 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 131 dta_bdy(ib_bdy)%v2d(ib) = pv2d(ii,ij) * vmask(ii,ij,1) 132 END DO 133 ENDIF 134 135 IF( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 0 ) THEN 136 IF( nn_dyn3d(ib_bdy) .eq. jp_frs ) THEN 137 ilen1(:) = nblen(:) 138 ELSE 139 ilen1(:) = nblenrim(:) 140 ENDIF 141 igrd = 2 142 DO ib = 1, ilen1(igrd) 333 143 DO ik = 1, jpkm1 334 tbdy(ib,ik) = tn(ii,ij,ik) 335 sbdy(ib,ik) = sn(ii,ij,ik) 144 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 145 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 146 dta_bdy(ib_bdy)%u3d(ib,ik) = ( un(ii,ij,ik) - pu2d(ii,ij) ) * umask(ii,ij,ik) 147 END DO 148 END DO 149 igrd = 3 150 DO ib = 1, ilen1(igrd) 151 DO ik = 1, jpkm1 152 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 153 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 154 dta_bdy(ib_bdy)%v3d(ib,ik) = ( vn(ii,ij,ik) - pv2d(ii,ij) ) * vmask(ii,ij,ik) 155 END DO 156 END DO 157 ENDIF 158 159 IF( nn_tra(ib_bdy) .gt. 0 .and. nn_tra_dta(ib_bdy) .eq. 0 ) THEN 160 IF( nn_tra(ib_bdy) .eq. jp_frs ) THEN 161 ilen1(:) = nblen(:) 162 ELSE 163 ilen1(:) = nblenrim(:) 164 ENDIF 165 igrd = 1 ! Everything is at T-points here 166 DO ib = 1, ilen1(igrd) 167 DO ik = 1, jpkm1 168 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 169 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 170 dta_bdy(ib_bdy)%tem(ib,ik) = tsn(ii,ij,ik,jp_tem) * tmask(ii,ij,ik) 171 dta_bdy(ib_bdy)%sal(ib,ik) = tsn(ii,ij,ik,jp_sal) * tmask(ii,ij,ik) 172 END DO 173 END DO 174 ENDIF 175 176 #if defined key_lim2 177 IF( nn_ice_lim2(ib_bdy) .gt. 0 .and. nn_ice_lim2_dta(ib_bdy) .eq. 0 ) THEN 178 IF( nn_ice_lim2(ib_bdy) .eq. jp_frs ) THEN 179 ilen1(:) = nblen(:) 180 ELSE 181 ilen1(:) = nblenrim(:) 182 ENDIF 183 igrd = 1 ! Everything is at T-points here 184 DO ib = 1, ilen1(igrd) 185 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 186 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 187 dta_bdy(ib_bdy)%frld(ib) = frld(ii,ij) * tmask(ii,ij,1) 188 dta_bdy(ib_bdy)%hicif(ib) = hicif(ii,ij) * tmask(ii,ij,1) 189 dta_bdy(ib_bdy)%hsnif(ib) = hsnif(ii,ij) * tmask(ii,ij,1) 190 END DO 191 ENDIF 192 #endif 193 194 ENDDO ! ib_bdy 195 196 CALL wrk_dealloc(jpi,jpj,pu2d,pv2d) 197 198 ENDIF ! kt .eq. nit000 199 200 ! update external data from files 201 !-------------------------------- 202 203 jstart = 1 204 DO ib_bdy = 1, nb_bdy 205 IF( nn_dta(ib_bdy) .eq. 1 ) THEN ! skip this bit if no external data required 206 207 IF( PRESENT(jit) ) THEN 208 ! Update barotropic boundary conditions only 209 ! jit is optional argument for fld_read and tide_update 210 IF( nn_dyn2d(ib_bdy) .gt. 0 ) THEN 211 IF( nn_dyn2d_dta(ib_bdy) .eq. 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 212 dta_bdy(ib_bdy)%ssh(:) = 0.0 213 dta_bdy(ib_bdy)%u2d(:) = 0.0 214 dta_bdy(ib_bdy)%v2d(:) = 0.0 215 ENDIF 216 IF( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) THEN ! update external data 217 jend = jstart + 2 218 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), & 219 & jit=jit, time_offset=time_offset ) 220 ENDIF 221 IF( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN ! update tidal harmonic forcing 222 CALL tide_update( kt=kt, idx=idx_bdy(ib_bdy), dta=dta_bdy(ib_bdy), td=tides(ib_bdy), & 223 & jit=jit, time_offset=time_offset ) 224 ENDIF 225 ENDIF 226 ELSE 227 IF( nn_dyn2d(ib_bdy) .gt. 0 .and. nn_dyn2d_dta(ib_bdy) .eq. 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 228 dta_bdy(ib_bdy)%ssh(:) = 0.0 229 dta_bdy(ib_bdy)%u2d(:) = 0.0 230 dta_bdy(ib_bdy)%v2d(:) = 0.0 231 ENDIF 232 IF( nb_bdy_fld(ib_bdy) .gt. 0 ) THEN ! update external data 233 jend = jstart + nb_bdy_fld(ib_bdy) - 1 234 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), time_offset=time_offset ) 235 ENDIF 236 IF( nn_dyn2d(ib_bdy) .gt. 0 .and. nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN ! update tidal harmonic forcing 237 CALL tide_update( kt=kt, idx=idx_bdy(ib_bdy), dta=dta_bdy(ib_bdy), td=tides(ib_bdy), time_offset=time_offset ) 238 ENDIF 239 ENDIF 240 jstart = jend+1 241 242 ! If full velocities in boundary data then split into barotropic and baroclinic data 243 ! (Note that we have already made sure that you can't use ln_full_vel = .true. at the same 244 ! time as the dynspg_ts option). 245 246 IF( ln_full_vel_array(ib_bdy) .and. & 247 & ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 .or. nn_dyn3d_dta(ib_bdy) .eq. 1 ) ) THEN 248 249 igrd = 2 ! zonal velocity 250 dta_bdy(ib_bdy)%u2d(:) = 0.0 251 DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 252 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 253 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 254 DO ik = 1, jpkm1 255 dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) & 256 & + fse3u(ii,ij,ik) * umask(ii,ij,ik) * dta_bdy(ib_bdy)%u3d(ib,ik) 257 END DO 258 dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) * hur(ii,ij) 259 DO ik = 1, jpkm1 260 dta_bdy(ib_bdy)%u3d(ib,ik) = dta_bdy(ib_bdy)%u3d(ib,ik) - dta_bdy(ib_bdy)%u2d(ib) 336 261 END DO 337 262 END DO 338 ENDIF 339 340 IF(ln_dyn_frs) THEN 341 igrd = 2 ! U-points data 342 DO ib = 1, nblen(igrd) 343 ii = nbi(ib,igrd) 344 ij = nbj(ib,igrd) 263 264 igrd = 3 ! meridional velocity 265 dta_bdy(ib_bdy)%v2d(:) = 0.0 266 DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 267 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 268 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 345 269 DO ik = 1, jpkm1 346 ubdy(ib,ik) = un(ii, ij, ik) 270 dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) & 271 & + fse3v(ii,ij,ik) * vmask(ii,ij,ik) * dta_bdy(ib_bdy)%v3d(ib,ik) 272 END DO 273 dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) * hvr(ii,ij) 274 DO ik = 1, jpkm1 275 dta_bdy(ib_bdy)%v3d(ib,ik) = dta_bdy(ib_bdy)%v3d(ib,ik) - dta_bdy(ib_bdy)%v2d(ib) 347 276 END DO 348 277 END DO 349 ! 350 igrd = 3 ! V-points data 351 DO ib = 1, nblen(igrd) 352 ii = nbi(ib,igrd) 353 ij = nbj(ib,igrd) 354 DO ik = 1, jpkm1 355 vbdy(ib,ik) = vn(ii, ij, ik) 356 END DO 357 END DO 358 ENDIF 359 ! 360 #if defined key_lim2 361 IF( ln_ice_frs ) THEN 362 igrd = 1 ! T-points data 363 DO ib = 1, nblen(igrd) 364 frld_bdy (ib) = frld(nbi(ib,igrd), nbj(ib,igrd)) 365 hicif_bdy(ib) = hicif(nbi(ib,igrd), nbj(ib,igrd)) 366 hsnif_bdy(ib) = hsnif(nbi(ib,igrd), nbj(ib,igrd)) 367 END DO 368 ENDIF 369 #endif 370 ELSEIF( nn_dtactl == 1 ) THEN ! Set first record in the climatological case: 371 ! 372 IF( ln_clim .AND. ntimes_bdy == 1 ) THEN 373 nbdy_a = 1 374 ELSEIF( ln_clim .AND. ntimes_bdy == iman ) THEN 375 nbdy_b = 0 376 nbdy_a = imois 278 279 ENDIF 280 281 END IF ! nn_dta(ib_bdy) = 1 282 END DO ! ib_bdy 283 284 IF( nn_timing == 1 ) CALL timing_stop('bdy_dta') 285 286 END SUBROUTINE bdy_dta 287 288 289 SUBROUTINE bdy_dta_init 290 !!---------------------------------------------------------------------- 291 !! *** SUBROUTINE bdy_dta_init *** 292 !! 293 !! ** Purpose : Initialise arrays for reading of external data 294 !! for open boundary conditions 295 !! 296 !! ** Method : Use fldread.F90 297 !! 298 !!---------------------------------------------------------------------- 299 USE dynspg_oce, ONLY: lk_dynspg_ts 300 !! 301 INTEGER :: ib_bdy, jfld, jstart, jend, ierror ! local indices 302 !! 303 CHARACTER(len=100) :: cn_dir ! Root directory for location of data files 304 CHARACTER(len=100), DIMENSION(nb_bdy) :: cn_dir_array ! Root directory for location of data files 305 LOGICAL :: ln_full_vel ! =T => full velocities in 3D boundary data 306 ! =F => baroclinic velocities in 3D boundary data 307 INTEGER :: ilen_global ! Max length required for global bdy dta arrays 308 INTEGER, DIMENSION(jpbgrd) :: ilen0 ! size of local arrays 309 INTEGER, ALLOCATABLE, DIMENSION(:) :: ilen1, ilen3 ! size of 1st and 3rd dimensions of local arrays 310 INTEGER, ALLOCATABLE, DIMENSION(:) :: ibdy ! bdy set for a particular jfld 311 INTEGER, ALLOCATABLE, DIMENSION(:) :: igrid ! index for grid type (1,2,3 = T,U,V) 312 INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts 313 TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: blf_i ! array of namelist information structures 314 TYPE(FLD_N) :: bn_tem, bn_sal, bn_u3d, bn_v3d ! 315 TYPE(FLD_N) :: bn_ssh, bn_u2d, bn_v2d ! informations about the fields to be read 316 #if defined key_lim2 317 TYPE(FLD_N) :: bn_frld, bn_hicif, bn_hsnif ! 318 #endif 319 NAMELIST/nambdy_dta/ cn_dir, bn_tem, bn_sal, bn_u3d, bn_v3d, bn_ssh, bn_u2d, bn_v2d 320 #if defined key_lim2 321 NAMELIST/nambdy_dta/ bn_frld, bn_hicif, bn_hsnif 322 #endif 323 NAMELIST/nambdy_dta/ ln_full_vel 324 !!--------------------------------------------------------------------------- 325 326 IF( nn_timing == 1 ) CALL timing_start('bdy_dta_init') 327 328 ! Set nn_dta 329 DO ib_bdy = 1, nb_bdy 330 nn_dta(ib_bdy) = MAX( nn_dyn2d_dta(ib_bdy) & 331 ,nn_dyn3d_dta(ib_bdy) & 332 ,nn_tra_dta(ib_bdy) & 333 #if defined key_lim2 334 ,nn_ice_lim2_dta(ib_bdy) & 335 #endif 336 ) 337 IF(nn_dta(ib_bdy) .gt. 1) nn_dta(ib_bdy) = 1 338 END DO 339 340 ! Work out upper bound of how many fields there are to read in and allocate arrays 341 ! --------------------------------------------------------------------------- 342 ALLOCATE( nb_bdy_fld(nb_bdy) ) 343 nb_bdy_fld(:) = 0 344 DO ib_bdy = 1, nb_bdy 345 IF( nn_dyn2d(ib_bdy) .gt. 0 .and. ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) THEN 346 nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 3 347 ENDIF 348 IF( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 1 ) THEN 349 nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 2 350 ENDIF 351 IF( nn_tra(ib_bdy) .gt. 0 .and. nn_tra_dta(ib_bdy) .eq. 1 ) THEN 352 nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 2 353 ENDIF 354 #if defined key_lim2 355 IF( nn_ice_lim2(ib_bdy) .gt. 0 .and. nn_ice_lim2_dta(ib_bdy) .eq. 1 ) THEN 356 nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 3 357 ENDIF 358 #endif 359 ENDDO 360 361 nb_bdy_fld_sum = SUM( nb_bdy_fld ) 362 363 ALLOCATE( bf(nb_bdy_fld_sum), STAT=ierror ) 364 IF( ierror > 0 ) THEN 365 CALL ctl_stop( 'bdy_dta: unable to allocate bf structure' ) ; RETURN 366 ENDIF 367 ALLOCATE( blf_i(nb_bdy_fld_sum), STAT=ierror ) 368 IF( ierror > 0 ) THEN 369 CALL ctl_stop( 'bdy_dta: unable to allocate blf_i structure' ) ; RETURN 370 ENDIF 371 ALLOCATE( nbmap_ptr(nb_bdy_fld_sum), STAT=ierror ) 372 IF( ierror > 0 ) THEN 373 CALL ctl_stop( 'bdy_dta: unable to allocate nbmap_ptr structure' ) ; RETURN 374 ENDIF 375 ALLOCATE( ilen1(nb_bdy_fld_sum), ilen3(nb_bdy_fld_sum) ) 376 ALLOCATE( ibdy(nb_bdy_fld_sum) ) 377 ALLOCATE( igrid(nb_bdy_fld_sum) ) 378 379 ! Read namelists 380 ! -------------- 381 REWIND(numnam) 382 jfld = 0 383 DO ib_bdy = 1, nb_bdy 384 IF( nn_dta(ib_bdy) .eq. 1 ) THEN 385 ! set file information 386 cn_dir = './' ! directory in which the model is executed 387 ln_full_vel = .false. 388 ! ... default values (NB: frequency positive => hours, negative => months) 389 ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation ! 390 ! ! name ! hours ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs ! 391 bn_ssh = FLD_N( 'bdy_ssh' , 24 , 'sossheig' , .false. , .false. , 'yearly' , '' , '' ) 392 bn_u2d = FLD_N( 'bdy_vel2d_u' , 24 , 'vobtcrtx' , .false. , .false. , 'yearly' , '' , '' ) 393 bn_v2d = FLD_N( 'bdy_vel2d_v' , 24 , 'vobtcrty' , .false. , .false. , 'yearly' , '' , '' ) 394 bn_u3d = FLD_N( 'bdy_vel3d_u' , 24 , 'vozocrtx' , .false. , .false. , 'yearly' , '' , '' ) 395 bn_v3d = FLD_N( 'bdy_vel3d_v' , 24 , 'vomecrty' , .false. , .false. , 'yearly' , '' , '' ) 396 bn_tem = FLD_N( 'bdy_tem' , 24 , 'votemper' , .false. , .false. , 'yearly' , '' , '' ) 397 bn_sal = FLD_N( 'bdy_sal' , 24 , 'vosaline' , .false. , .false. , 'yearly' , '' , '' ) 398 #if defined key_lim2 399 bn_frld = FLD_N( 'bdy_frld' , 24 , 'ildsconc' , .false. , .false. , 'yearly' , '' , '' ) 400 bn_hicif = FLD_N( 'bdy_hicif' , 24 , 'iicethic' , .false. , .false. , 'yearly' , '' , '' ) 401 bn_hsnif = FLD_N( 'bdy_hsnif' , 24 , 'isnothic' , .false. , .false. , 'yearly' , '' , '' ) 402 #endif 403 404 ! Important NOT to rewind here. 405 READ( numnam, nambdy_dta ) 406 407 cn_dir_array(ib_bdy) = cn_dir 408 ln_full_vel_array(ib_bdy) = ln_full_vel 409 410 IF( ln_full_vel_array(ib_bdy) .and. lk_dynspg_ts ) THEN 411 CALL ctl_stop( 'bdy_dta_init: ERROR, cannot specify full velocities in boundary data',& 412 & 'with dynspg_ts option' ) ; RETURN 413 ENDIF 414 415 nblen => idx_bdy(ib_bdy)%nblen 416 nblenrim => idx_bdy(ib_bdy)%nblenrim 417 418 ! Only read in necessary fields for this set. 419 ! Important that barotropic variables come first. 420 IF( nn_dyn2d(ib_bdy) .gt. 0 .and. ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) THEN 421 422 IF( nn_dyn2d(ib_bdy) .ne. jp_frs ) THEN 423 jfld = jfld + 1 424 blf_i(jfld) = bn_ssh 425 ibdy(jfld) = ib_bdy 426 igrid(jfld) = 1 427 ilen1(jfld) = nblenrim(igrid(jfld)) 428 ilen3(jfld) = 1 429 ENDIF 430 431 IF( .not. ln_full_vel_array(ib_bdy) ) THEN 432 433 jfld = jfld + 1 434 blf_i(jfld) = bn_u2d 435 ibdy(jfld) = ib_bdy 436 igrid(jfld) = 2 437 IF( nn_dyn2d(ib_bdy) .eq. jp_frs ) THEN 438 ilen1(jfld) = nblen(igrid(jfld)) 439 ELSE 440 ilen1(jfld) = nblenrim(igrid(jfld)) 441 ENDIF 442 ilen3(jfld) = 1 443 444 jfld = jfld + 1 445 blf_i(jfld) = bn_v2d 446 ibdy(jfld) = ib_bdy 447 igrid(jfld) = 3 448 IF( nn_dyn2d(ib_bdy) .eq. jp_frs ) THEN 449 ilen1(jfld) = nblen(igrid(jfld)) 450 ELSE 451 ilen1(jfld) = nblenrim(igrid(jfld)) 452 ENDIF 453 ilen3(jfld) = 1 454 455 ENDIF 456 457 ENDIF 458 459 ! baroclinic velocities 460 IF( ( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 1 ) .or. & 461 & ( ln_full_vel_array(ib_bdy) .and. nn_dyn2d(ib_bdy) .gt. 0 .and. & 462 & ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) ) THEN 463 464 jfld = jfld + 1 465 blf_i(jfld) = bn_u3d 466 ibdy(jfld) = ib_bdy 467 igrid(jfld) = 2 468 IF( nn_dyn3d(ib_bdy) .eq. jp_frs ) THEN 469 ilen1(jfld) = nblen(igrid(jfld)) 470 ELSE 471 ilen1(jfld) = nblenrim(igrid(jfld)) 472 ENDIF 473 ilen3(jfld) = jpk 474 475 jfld = jfld + 1 476 blf_i(jfld) = bn_v3d 477 ibdy(jfld) = ib_bdy 478 igrid(jfld) = 3 479 IF( nn_dyn3d(ib_bdy) .eq. jp_frs ) THEN 480 ilen1(jfld) = nblen(igrid(jfld)) 481 ELSE 482 ilen1(jfld) = nblenrim(igrid(jfld)) 483 ENDIF 484 ilen3(jfld) = jpk 485 486 ENDIF 487 488 ! temperature and salinity 489 IF( nn_tra(ib_bdy) .gt. 0 .and. nn_tra_dta(ib_bdy) .eq. 1 ) THEN 490 491 jfld = jfld + 1 492 blf_i(jfld) = bn_tem 493 ibdy(jfld) = ib_bdy 494 igrid(jfld) = 1 495 IF( nn_tra(ib_bdy) .eq. jp_frs ) THEN 496 ilen1(jfld) = nblen(igrid(jfld)) 497 ELSE 498 ilen1(jfld) = nblenrim(igrid(jfld)) 499 ENDIF 500 ilen3(jfld) = jpk 501 502 jfld = jfld + 1 503 blf_i(jfld) = bn_sal 504 ibdy(jfld) = ib_bdy 505 igrid(jfld) = 1 506 IF( nn_tra(ib_bdy) .eq. jp_frs ) THEN 507 ilen1(jfld) = nblen(igrid(jfld)) 508 ELSE 509 ilen1(jfld) = nblenrim(igrid(jfld)) 510 ENDIF 511 ilen3(jfld) = jpk 512 513 ENDIF 514 515 #if defined key_lim2 516 ! sea ice 517 IF( nn_ice_lim2(ib_bdy) .gt. 0 .and. nn_ice_lim2_dta(ib_bdy) .eq. 1 ) THEN 518 519 jfld = jfld + 1 520 blf_i(jfld) = bn_frld 521 ibdy(jfld) = ib_bdy 522 igrid(jfld) = 1 523 IF( nn_ice_lim2(ib_bdy) .eq. jp_frs ) THEN 524 ilen1(jfld) = nblen(igrid(jfld)) 525 ELSE 526 ilen1(jfld) = nblenrim(igrid(jfld)) 527 ENDIF 528 ilen3(jfld) = 1 529 530 jfld = jfld + 1 531 blf_i(jfld) = bn_hicif 532 ibdy(jfld) = ib_bdy 533 igrid(jfld) = 1 534 IF( nn_ice_lim2(ib_bdy) .eq. jp_frs ) THEN 535 ilen1(jfld) = nblen(igrid(jfld)) 536 ELSE 537 ilen1(jfld) = nblenrim(igrid(jfld)) 538 ENDIF 539 ilen3(jfld) = 1 540 541 jfld = jfld + 1 542 blf_i(jfld) = bn_hsnif 543 ibdy(jfld) = ib_bdy 544 igrid(jfld) = 1 545 IF( nn_ice_lim2(ib_bdy) .eq. jp_frs ) THEN 546 ilen1(jfld) = nblen(igrid(jfld)) 547 ELSE 548 ilen1(jfld) = nblenrim(igrid(jfld)) 549 ENDIF 550 ilen3(jfld) = 1 551 552 ENDIF 553 #endif 554 ! Recalculate field counts 555 !------------------------- 556 nb_bdy_fld_sum = 0 557 IF( ib_bdy .eq. 1 ) THEN 558 nb_bdy_fld(ib_bdy) = jfld 559 nb_bdy_fld_sum = jfld 377 560 ELSE 378 nbdy_a = nbdy_b 379 ENDIF 380 381 ! Read first record: 382 ipj = 1 383 ipk = jpk 384 igrd = 1 385 ipi = nblendta(igrd) 386 387 IF(ln_tra_frs) THEN 388 ! 389 igrd = 1 ! Temperature 390 IF( nblendta(igrd) <= 0 ) THEN 391 idvar = iom_varid( numbdyt, 'votemper' ) 392 nblendta(igrd) = iom_file(numbdyt)%dimsz(1,idvar) 393 ENDIF 394 IF(lwp) WRITE(numout,*) 'Dim size for votemper is ', nblendta(igrd) 395 ipi = nblendta(igrd) 396 CALL iom_get ( numbdyt, jpdom_unknown, 'votemper', zdta(1:ipi,1:ipj,1:ipk), nbdy_a ) 397 ! 398 DO ib = 1, nblen(igrd) 399 DO ik = 1, jpkm1 400 tbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 401 END DO 402 END DO 403 ! 404 igrd = 1 ! salinity 405 IF( nblendta(igrd) .le. 0 ) THEN 406 idvar = iom_varid( numbdyt, 'vosaline' ) 407 nblendta(igrd) = iom_file(numbdyt)%dimsz(1,idvar) 408 ENDIF 409 IF(lwp) WRITE(numout,*) 'Dim size for vosaline is ', nblendta(igrd) 410 ipi = nblendta(igrd) 411 CALL iom_get ( numbdyt, jpdom_unknown, 'vosaline', zdta(1:ipi,1:ipj,1:ipk), nbdy_a ) 412 ! 413 DO ib = 1, nblen(igrd) 414 DO ik = 1, jpkm1 415 sbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 416 END DO 417 END DO 418 ENDIF ! ln_tra_frs 419 420 IF( ln_dyn_frs ) THEN 421 ! 422 igrd = 2 ! u-velocity 423 IF ( nblendta(igrd) .le. 0 ) THEN 424 idvar = iom_varid( numbdyu,'vozocrtx' ) 425 nblendta(igrd) = iom_file(numbdyu)%dimsz(1,idvar) 426 ENDIF 427 IF(lwp) WRITE(numout,*) 'Dim size for vozocrtx is ', nblendta(igrd) 428 ipi = nblendta(igrd) 429 CALL iom_get ( numbdyu, jpdom_unknown,'vozocrtx',zdta(1:ipi,1:ipj,1:ipk),nbdy_a ) 430 DO ib = 1, nblen(igrd) 431 DO ik = 1, jpkm1 432 ubdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 433 END DO 434 END DO 435 ! 436 igrd = 3 ! v-velocity 437 IF ( nblendta(igrd) .le. 0 ) THEN 438 idvar = iom_varid( numbdyv,'vomecrty' ) 439 nblendta(igrd) = iom_file(numbdyv)%dimsz(1,idvar) 440 ENDIF 441 IF(lwp) WRITE(numout,*) 'Dim size for vomecrty is ', nblendta(igrd) 442 ipi = nblendta(igrd) 443 CALL iom_get ( numbdyv, jpdom_unknown,'vomecrty',zdta(1:ipi,1:ipj,1:ipk),nbdy_a ) 444 DO ib = 1, nblen(igrd) 445 DO ik = 1, jpkm1 446 vbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 447 END DO 448 END DO 449 ENDIF ! ln_dyn_frs 450 451 #if defined key_lim2 452 IF( ln_ice_frs ) THEN 453 ! 454 igrd=1 ! leads fraction 455 IF(lwp) WRITE(numout,*) 'Dim size for ildsconc is ',nblendta(igrd) 456 ipi=nblendta(igrd) 457 CALL iom_get ( numbdyt, jpdom_unknown,'ildsconc',zdta(1:ipi,:,1),nbdy_a ) 458 DO ib=1, nblen(igrd) 459 frld_bdydta(ib,2) = zdta(nbmap(ib,igrd),1,1) 460 END DO 461 ! 462 igrd=1 ! ice thickness 463 IF(lwp) WRITE(numout,*) 'Dim size for iicethic is ',nblendta(igrd) 464 ipi=nblendta(igrd) 465 CALL iom_get ( numbdyt, jpdom_unknown,'iicethic',zdta(1:ipi,:,1),nbdy_a ) 466 DO ib=1, nblen(igrd) 467 hicif_bdydta(ib,2) = zdta(nbmap(ib,igrd),1,1) 468 END DO 469 ! 470 igrd=1 ! snow thickness 471 IF(lwp) WRITE(numout,*) 'Dim size for isnowthi is ',nblendta(igrd) 472 ipi=nblendta(igrd) 473 CALL iom_get ( numbdyt, jpdom_unknown,'isnowthi',zdta(1:ipi,:,1),nbdy_a ) 474 DO ib=1, nblen(igrd) 475 hsnif_bdydta(ib,2) = zdta(nbmap(ib,igrd),1,1) 476 END DO 477 ENDIF ! just if ln_ice_frs is set 478 #endif 479 480 IF( .NOT.ln_clim .AND. istep(1) > 0 ) THEN ! First data time is after start of run 481 nbdy_b = nbdy_a ! Put first value in both time levels 482 IF( ln_tra_frs ) THEN 483 tbdydta(:,:,1) = tbdydta(:,:,2) 484 sbdydta(:,:,1) = sbdydta(:,:,2) 485 ENDIF 486 IF( ln_dyn_frs ) THEN 487 ubdydta(:,:,1) = ubdydta(:,:,2) 488 vbdydta(:,:,1) = vbdydta(:,:,2) 489 ENDIF 490 #if defined key_lim2 491 IF( ln_ice_frs ) THEN 492 frld_bdydta (:,1) = frld_bdydta(:,2) 493 hicif_bdydta(:,1) = hicif_bdydta(:,2) 494 hsnif_bdydta(:,1) = hsnif_bdydta(:,2) 495 ENDIF 496 #endif 497 END IF 498 ! 499 END IF ! nn_dtactl == 0/1 500 501 ! In the case of constant boundary forcing fill bdy arrays once for all 502 IF( ln_clim .AND. ntimes_bdy == 1 ) THEN 503 IF( ln_tra_frs ) THEN 504 tbdy (:,:) = tbdydta (:,:,2) 505 sbdy (:,:) = sbdydta (:,:,2) 506 ENDIF 507 IF( ln_dyn_frs) THEN 508 ubdy (:,:) = ubdydta (:,:,2) 509 vbdy (:,:) = vbdydta (:,:,2) 510 ENDIF 511 #if defined key_lim2 512 IF( ln_ice_frs ) THEN 513 frld_bdy (:) = frld_bdydta (:,2) 514 hicif_bdy(:) = hicif_bdydta(:,2) 515 hsnif_bdy(:) = hsnif_bdydta(:,2) 516 ENDIF 517 #endif 518 519 IF( ln_tra_frs .OR. ln_ice_frs) CALL iom_close( numbdyt ) 520 IF( ln_dyn_frs ) CALL iom_close( numbdyu ) 521 IF( ln_dyn_frs ) CALL iom_close( numbdyv ) 522 END IF 523 ! 524 ENDIF ! End if nit000 525 526 527 ! !---------------------! 528 IF( nn_dtactl == 1 .AND. ntimes_bdy > 1 ) THEN ! at each time step ! 529 ! !---------------------! 530 ! Read one more record if necessary 531 !********************************** 532 533 IF( ln_clim .AND. imois /= nbdy_b ) THEN ! remember that nbdy_b=0 for kt=nit000 534 nbdy_b = imois 535 nbdy_a = imois + 1 536 nbdy_b = MOD( nbdy_b, iman ) ; IF( nbdy_b == 0 ) nbdy_b = iman 537 nbdy_a = MOD( nbdy_a, iman ) ; IF( nbdy_a == 0 ) nbdy_a = iman 538 lect=.true. 539 ELSEIF( .NOT.ln_clim .AND. itimer >= istep(nbdy_a) ) THEN 540 541 IF( nbdy_a < ntimes_bdy ) THEN 542 nbdy_b = nbdy_a 543 nbdy_a = nbdy_a + 1 544 lect =.true. 561 nb_bdy_fld(ib_bdy) = jfld - nb_bdy_fld_sum 562 nb_bdy_fld_sum = nb_bdy_fld_sum + nb_bdy_fld(ib_bdy) 563 ENDIF 564 565 ENDIF ! nn_dta .eq. 1 566 ENDDO ! ib_bdy 567 568 569 DO jfld = 1, nb_bdy_fld_sum 570 ALLOCATE( bf(jfld)%fnow(ilen1(jfld),1,ilen3(jfld)) ) 571 IF( blf_i(jfld)%ln_tint ) ALLOCATE( bf(jfld)%fdta(ilen1(jfld),1,ilen3(jfld),2) ) 572 nbmap_ptr(jfld)%ptr => idx_bdy(ibdy(jfld))%nbmap(:,igrid(jfld)) 573 ENDDO 574 575 ! fill bf with blf_i and control print 576 !------------------------------------- 577 jstart = 1 578 DO ib_bdy = 1, nb_bdy 579 jend = jstart + nb_bdy_fld(ib_bdy) - 1 580 CALL fld_fill( bf(jstart:jend), blf_i(jstart:jend), cn_dir_array(ib_bdy), 'bdy_dta', & 581 & 'open boundary conditions', 'nambdy_dta' ) 582 jstart = jend + 1 583 ENDDO 584 585 ! Initialise local boundary data arrays 586 ! nn_xxx_dta=0 : allocate space - will be filled from initial conditions later 587 ! nn_xxx_dta=1 : point to "fnow" arrays 588 !------------------------------------- 589 590 jfld = 0 591 DO ib_bdy=1, nb_bdy 592 593 nblen => idx_bdy(ib_bdy)%nblen 594 nblenrim => idx_bdy(ib_bdy)%nblenrim 595 596 IF (nn_dyn2d(ib_bdy) .gt. 0) THEN 597 IF( nn_dyn2d_dta(ib_bdy) .eq. 0 .or. nn_dyn2d_dta(ib_bdy) .eq. 2 .or. ln_full_vel_array(ib_bdy) ) THEN 598 IF( nn_dyn2d(ib_bdy) .eq. jp_frs ) THEN 599 ilen0(1:3) = nblen(1:3) 600 ELSE 601 ilen0(1:3) = nblenrim(1:3) 602 ENDIF 603 ALLOCATE( dta_bdy(ib_bdy)%ssh(ilen0(1)) ) 604 ALLOCATE( dta_bdy(ib_bdy)%u2d(ilen0(2)) ) 605 ALLOCATE( dta_bdy(ib_bdy)%v2d(ilen0(3)) ) 545 606 ELSE 546 ! We have reached the end of the file 547 ! put the last data time into both time levels 548 nbdy_b = nbdy_a 549 IF(ln_tra_frs) THEN 550 tbdydta(:,:,1) = tbdydta(:,:,2) 551 sbdydta(:,:,1) = sbdydta(:,:,2) 552 ENDIF 553 IF(ln_dyn_frs) THEN 554 ubdydta(:,:,1) = ubdydta(:,:,2) 555 vbdydta(:,:,1) = vbdydta(:,:,2) 556 ENDIF 557 #if defined key_lim2 558 IF(ln_ice_frs) THEN 559 frld_bdydta (:,1) = frld_bdydta (:,2) 560 hicif_bdydta(:,1) = hicif_bdydta(:,2) 561 hsnif_bdydta(:,1) = hsnif_bdydta(:,2) 562 ENDIF 563 #endif 564 END IF ! nbdy_a < ntimes_bdy 565 ! 566 END IF 567 568 IF( lect ) THEN ! Swap arrays 569 IF( ln_tra_frs ) THEN 570 tbdydta(:,:,1) = tbdydta(:,:,2) 571 sbdydta(:,:,1) = sbdydta(:,:,2) 572 ENDIF 573 IF( ln_dyn_frs ) THEN 574 ubdydta(:,:,1) = ubdydta(:,:,2) 575 vbdydta(:,:,1) = vbdydta(:,:,2) 576 ENDIF 577 #if defined key_lim2 578 IF( ln_ice_frs ) THEN 579 frld_bdydta (:,1) = frld_bdydta (:,2) 580 hicif_bdydta(:,1) = hicif_bdydta(:,2) 581 hsnif_bdydta(:,1) = hsnif_bdydta(:,2) 582 ENDIF 583 #endif 584 ! read another set 585 ipj = 1 586 ipk = jpk 587 588 IF( ln_tra_frs ) THEN 589 ! 590 igrd = 1 ! temperature 591 ipi = nblendta(igrd) 592 CALL iom_get ( numbdyt, jpdom_unknown, 'votemper', zdta(1:ipi,1:ipj,1:ipk), nbdy_a ) 593 DO ib = 1, nblen(igrd) 594 DO ik = 1, jpkm1 595 tbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 596 END DO 597 END DO 598 ! 599 igrd = 1 ! salinity 600 ipi = nblendta(igrd) 601 CALL iom_get ( numbdyt, jpdom_unknown, 'vosaline', zdta(1:ipi,1:ipj,1:ipk), nbdy_a ) 602 DO ib = 1, nblen(igrd) 603 DO ik = 1, jpkm1 604 sbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 605 END DO 606 END DO 607 ENDIF ! ln_tra_frs 608 609 IF(ln_dyn_frs) THEN 610 ! 611 igrd = 2 ! u-velocity 612 ipi = nblendta(igrd) 613 CALL iom_get ( numbdyu, jpdom_unknown,'vozocrtx',zdta(1:ipi,1:ipj,1:ipk),nbdy_a ) 614 DO ib = 1, nblen(igrd) 615 DO ik = 1, jpkm1 616 ubdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 617 END DO 618 END DO 619 ! 620 igrd = 3 ! v-velocity 621 ipi = nblendta(igrd) 622 CALL iom_get ( numbdyv, jpdom_unknown,'vomecrty',zdta(1:ipi,1:ipj,1:ipk),nbdy_a ) 623 DO ib = 1, nblen(igrd) 624 DO ik = 1, jpkm1 625 vbdydta(ib,ik,2) = zdta(nbmap(ib,igrd),1,ik) 626 END DO 627 END DO 628 ENDIF ! ln_dyn_frs 629 ! 630 #if defined key_lim2 631 IF(ln_ice_frs) THEN 632 ! 633 igrd = 1 ! ice concentration 634 ipi=nblendta(igrd) 635 CALL iom_get ( numbdyt, jpdom_unknown,'ildsconc',zdta(1:ipi,:,1),nbdy_a ) 636 DO ib=1, nblen(igrd) 637 frld_bdydta(ib,2) = zdta( nbmap(ib,igrd), 1, 1 ) 638 END DO 639 ! 640 igrd=1 ! ice thickness 641 ipi=nblendta(igrd) 642 CALL iom_get ( numbdyt, jpdom_unknown,'iicethic',zdta(1:ipi,:,1),nbdy_a ) 643 DO ib=1, nblen(igrd) 644 hicif_bdydta(ib,2) = zdta( nbmap(ib,igrd), 1, 1 ) 645 END DO 646 ! 647 igrd=1 ! snow thickness 648 ipi=nblendta(igrd) 649 CALL iom_get ( numbdyt, jpdom_unknown,'isnowthi',zdta(1:ipi,:,1),nbdy_a ) 650 DO ib=1, nblen(igrd) 651 hsnif_bdydta(ib,2) = zdta( nbmap(ib,igrd), 1, 1 ) 652 END DO 653 ENDIF ! ln_ice_frs 654 #endif 655 ! 656 IF(lwp) WRITE(numout,*) 'bdy_dta_frs : first record file used nbdy_b ',nbdy_b 657 IF(lwp) WRITE(numout,*) '~~~~~~~~ last record file used nbdy_a ',nbdy_a 658 IF (.NOT.ln_clim) THEN 659 IF(lwp) WRITE(numout,*) 'first record time (s): ', istep(nbdy_b) 660 IF(lwp) WRITE(numout,*) 'model time (s) : ', itimer 661 IF(lwp) WRITE(numout,*) 'second record time (s): ', istep(nbdy_a) 662 ENDIF 663 ! 664 ENDIF ! end lect=.true. 665 666 667 ! Interpolate linearly 668 ! ******************** 669 ! 670 IF( ln_clim ) THEN ; zxy = REAL( nday ) / REAL( nmonth_len(nbdy_b) ) + 0.5 - i15 671 ELSEIF( istep(nbdy_b) == istep(nbdy_a) ) THEN 672 zxy = 0.0_wp 673 ELSE ; zxy = REAL( istep(nbdy_b) - itimer ) / REAL( istep(nbdy_b) - istep(nbdy_a) ) 674 END IF 675 676 IF(ln_tra_frs) THEN 677 igrd = 1 ! temperature & salinity 678 DO ib = 1, nblen(igrd) 679 DO ik = 1, jpkm1 680 tbdy(ib,ik) = zxy * tbdydta(ib,ik,2) + (1.-zxy) * tbdydta(ib,ik,1) 681 sbdy(ib,ik) = zxy * sbdydta(ib,ik,2) + (1.-zxy) * sbdydta(ib,ik,1) 682 END DO 683 END DO 684 ENDIF 685 686 IF(ln_dyn_frs) THEN 687 igrd = 2 ! u-velocity 688 DO ib = 1, nblen(igrd) 689 DO ik = 1, jpkm1 690 ubdy(ib,ik) = zxy * ubdydta(ib,ik,2) + (1.-zxy) * ubdydta(ib,ik,1) 691 END DO 692 END DO 693 ! 694 igrd = 3 ! v-velocity 695 DO ib = 1, nblen(igrd) 696 DO ik = 1, jpkm1 697 vbdy(ib,ik) = zxy * vbdydta(ib,ik,2) + (1.-zxy) * vbdydta(ib,ik,1) 698 END DO 699 END DO 700 ENDIF 701 702 #if defined key_lim2 703 IF(ln_ice_frs) THEN 704 igrd=1 705 DO ib=1, nblen(igrd) 706 frld_bdy(ib) = zxy * frld_bdydta(ib,2) + (1.-zxy) * frld_bdydta(ib,1) 707 hicif_bdy(ib) = zxy * hicif_bdydta(ib,2) + (1.-zxy) * hicif_bdydta(ib,1) 708 hsnif_bdy(ib) = zxy * hsnif_bdydta(ib,2) + (1.-zxy) * hsnif_bdydta(ib,1) 709 END DO 710 ENDIF ! just if ln_ice_frs is true 711 #endif 712 713 END IF !end if ((nn_dtactl==1).AND.(ntimes_bdy>1)) 714 715 716 ! !---------------------! 717 ! ! last call ! 718 ! !---------------------! 719 IF( kt == nitend ) THEN 720 IF(ln_tra_frs .or. ln_ice_frs) CALL iom_close( numbdyt ) ! Closing of the 3 files 721 IF(ln_dyn_frs) CALL iom_close( numbdyu ) 722 IF(ln_dyn_frs) CALL iom_close( numbdyv ) 723 ENDIF 724 ! 725 ENDIF ! ln_dyn_frs .OR. ln_tra_frs 726 ! 727 END SUBROUTINE bdy_dta_frs 728 729 730 SUBROUTINE bdy_dta_fla( kt, jit, icycl ) 731 !!--------------------------------------------------------------------------- 732 !! *** SUBROUTINE bdy_dta_fla *** 733 !! 734 !! ** Purpose : Read unstructured boundary data for Flather condition 735 !! 736 !! ** Method : At the first timestep, read in boundary data for two 737 !! times from the file and time-interpolate. At other 738 !! timesteps, check to see if we need another time from 739 !! the file. If so read it in. Time interpolate. 740 !!--------------------------------------------------------------------------- 741 !!gm DOCTOR names : argument integer : start with "k" 742 INTEGER, INTENT( in ) :: kt ! ocean time-step index 743 INTEGER, INTENT( in ) :: jit ! barotropic time step index 744 INTEGER, INTENT( in ) :: icycl ! number of cycles need for final file close 745 ! ! (for timesplitting option, otherwise zero) 746 !! 747 LOGICAL :: lect ! flag for reading 748 INTEGER :: it, ib, igrd ! dummy loop indices 749 INTEGER :: idvar ! netcdf var ID 750 INTEGER :: iman, i15, imois ! Time variables for monthly clim forcing 751 INTEGER :: ntimes_bdyt, ntimes_bdyu, ntimes_bdyv 752 INTEGER :: itimer, totime 753 INTEGER :: ipi, ipj, ipk, inum ! temporary integers (NetCDF read) 754 INTEGER :: iyear0, imonth0, iday0 755 INTEGER :: ihours0, iminutes0, isec0 756 INTEGER :: iyear, imonth, iday, isecs 757 INTEGER, DIMENSION(jpbtime) :: istept, istepu, istepv ! time arrays from data files 758 REAL(wp) :: dayfrac, zxy, zoffsett 759 REAL(wp) :: zoffsetu, zoffsetv 760 REAL(wp) :: dayjul0, zdayjulini 761 REAL(wp) :: zinterval_s, zinterval_e ! First and last interval in time axis 762 REAL(wp), DIMENSION(jpbtime) :: zstepr ! REAL time array from data files 763 REAL(wp), DIMENSION(jpbdta,1) :: zdta ! temporary array for data fields 764 CHARACTER(LEN=80), DIMENSION(6) :: clfile 765 CHARACTER(LEN=70 ) :: clunits ! units attribute of time coordinate 766 !!--------------------------------------------------------------------------- 767 768 !!gm add here the same style as in bdy_dta_frs 769 !!gm clearly bdy_dta_fla and bdy_dta_frs can be combined... 770 !!gm too many things duplicated in the read of data... simplification can be done 771 772 ! -------------------- ! 773 ! Initialization ! 774 ! -------------------- ! 775 776 lect = .false. ! If true, read a time record 777 778 ! Some time variables for monthly climatological forcing: 779 ! ******************************************************* 780 !!gm here use directely daymod variables 781 782 iman = INT( raamo ) ! Number of months in a year 783 784 i15 = INT( 2*REAL( nday, wp ) / ( REAL( nmonth_len(nmonth), wp ) + 0.5 ) ) 785 ! i15=0 if the current day is in the first half of the month, else i15=1 786 787 imois = nmonth + i15 - 1 ! imois is the first month record 788 IF( imois == 0 ) imois = iman 789 790 ! Time variable for non-climatological forcing: 791 ! ********************************************* 792 793 itimer = ((kt-1)-nit000+1)*rdt ! current time in seconds for interpolation 794 itimer = itimer + jit*rdt/REAL(nn_baro,wp) ! in non-climatological case 795 796 IF ( ln_tides ) THEN 797 798 ! -------------------------------------! 799 ! Update BDY fields with tidal forcing ! 800 ! -------------------------------------! 801 802 CALL tide_update( kt, jit ) 803 804 ENDIF 805 806 IF ( ln_dyn_fla ) THEN 807 808 ! -------------------------------------! 809 ! Update BDY fields with model data ! 810 ! -------------------------------------! 811 812 ! !-------------------! 813 IF( kt == nit000 .and. jit ==2 ) THEN ! First call only ! 814 ! !-------------------! 815 istep_bt(:) = 0 816 nbdy_b_bt = 0 817 nbdy_a_bt = 0 818 819 ! Get time information from bdy data file 820 ! *************************************** 821 822 IF(lwp) WRITE(numout,*) 823 IF(lwp) WRITE(numout,*) 'bdy_dta_fla :Initialize unstructured boundary data for barotropic variables.' 824 IF(lwp) WRITE(numout,*) '~~~~~~~' 825 826 IF( nn_dtactl == 0 ) THEN 827 IF(lwp) WRITE(numout,*) 'Bdy data are taken from initial conditions' 828 829 ELSEIF (nn_dtactl == 1) THEN 830 IF(lwp) WRITE(numout,*) 'Bdy data are read in netcdf files' 831 832 dayfrac = adatrj - REAL(itimer,wp)/86400. ! day fraction at time step kt-1 833 dayfrac = dayfrac - INT (dayfrac) ! 834 totime = (nitend-nit000+1)*rdt ! Total time of the run to verify that all the 835 ! necessary time dumps in file are included 836 837 clfile(4) = cn_dta_fla_T 838 clfile(5) = cn_dta_fla_U 839 clfile(6) = cn_dta_fla_V 840 841 DO igrd = 4,6 842 843 CALL iom_open( clfile(igrd), inum ) 844 CALL iom_gettime( inum, zstepr, kntime=ntimes_bdy_bt, cdunits=clunits ) 845 846 SELECT CASE( igrd ) 847 CASE (4) 848 numbdyt_bt = inum 849 CASE (5) 850 numbdyu_bt = inum 851 CASE (6) 852 numbdyv_bt = inum 853 END SELECT 854 855 ! Calculate time offset 856 READ(clunits,7000) iyear0, imonth0, iday0, ihours0, iminutes0, isec0 857 ! Convert time origin in file to julian days 858 isec0 = isec0 + ihours0*60.*60. + iminutes0*60. 859 CALL ymds2ju(iyear0, imonth0, iday0, REAL(isec0, wp), dayjul0) 860 ! Compute model initialization time 861 iyear = ndastp / 10000 862 imonth = ( ndastp - iyear * 10000 ) / 100 863 iday = ndastp - iyear * 10000 - imonth * 100 864 isecs = dayfrac * 86400 865 CALL ymds2ju(iyear, imonth, iday, REAL(isecs, wp) , zdayjulini) 866 ! zoffset from initialization date: 867 zoffset = (dayjul0-zdayjulini)*86400 868 ! 869 870 7000 FORMAT('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2) 871 872 !! TO BE DONE... Check consistency between calendar from file 873 !! (available optionally from iom_gettime) and calendar in model 874 !! when calendar in model available outside of IOIPSL. 875 876 ! Check that there are not too many times in the file. 877 IF (ntimes_bdy_bt > jpbtime) CALL ctl_stop( & 878 'Number of time dumps in bdy file exceed jpbtime parameter', & 879 'Check file:' // TRIM(clfile(igrd)) ) 880 881 ! Check that time array increases (or interp will fail): 882 DO it = 2, ntimes_bdy_bt 883 IF ( zstepr(it-1) >= zstepr(it) ) THEN 884 CALL ctl_stop('Time array in unstructured boundary data file', & 885 'does not continuously increase.', & 886 'Check file:' // TRIM(clfile(igrd)) ) 887 EXIT 888 END IF 889 END DO 890 891 IF ( .NOT. ln_clim ) THEN 892 ! Check that times in file span model run time: 893 894 ! Note: the fields may be time means, so we allow nit000 to be before 895 ! first time in the file, provided that it falls inside the meaning 896 ! period of the first field. Until we can get the meaning period 897 ! from the file, use the interval between fields as a proxy. 898 ! If nit000 is before the first time, use the value at first time 899 ! instead of extrapolating. This is done by putting time 1 into 900 ! both time levels. 901 ! The same applies to the last time level: see setting of lect below. 902 903 IF ( ntimes_bdy_bt == 1 ) CALL ctl_stop( & 904 'There is only one time dump in data files', & 905 'Set ln_clim=.true. in namelist for constant bdy forcing.' ) 906 907 zinterval_s = zstepr(2) - zstepr(1) 908 zinterval_e = zstepr(ntimes_bdy_bt) - zstepr(ntimes_bdy_bt-1) 909 910 IF( zstepr(1) + zoffset > 0 ) THEN 911 WRITE(ctmp1,*) 'Check file: ', clfile(igrd) 912 CALL ctl_stop( 'First time dump in bdy file is after model initial time', ctmp1 ) 913 END IF 914 IF( zstepr(ntimes_bdy_bt) + zoffset < totime ) THEN 915 WRITE(ctmp1,*) 'Check file: ', clfile(igrd) 916 CALL ctl_stop( 'Last time dump in bdy file is before model final time', ctmp1 ) 917 END IF 918 END IF ! .NOT. ln_clim 919 920 IF ( igrd .EQ. 4) THEN 921 ntimes_bdyt = ntimes_bdy_bt 922 zoffsett = zoffset 923 istept(:) = INT( zstepr(:) + zoffset ) 924 ELSE IF (igrd .EQ. 5) THEN 925 ntimes_bdyu = ntimes_bdy_bt 926 zoffsetu = zoffset 927 istepu(:) = INT( zstepr(:) + zoffset ) 928 ELSE IF (igrd .EQ. 6) THEN 929 ntimes_bdyv = ntimes_bdy_bt 930 zoffsetv = zoffset 931 istepv(:) = INT( zstepr(:) + zoffset ) 932 ENDIF 933 934 ENDDO 935 936 ! Verify time consistency between files 937 938 IF ( ntimes_bdyu /= ntimes_bdyt .OR. ntimes_bdyv /= ntimes_bdyt ) THEN 939 CALL ctl_stop( & 940 'Time axis lengths differ between bdy data files', & 941 'Multiple time frequencies not implemented yet' ) 942 ELSE 943 ntimes_bdy_bt = ntimes_bdyt 944 ENDIF 945 946 IF (zoffsetu.NE.zoffsett .OR. zoffsetv.NE.zoffsett) THEN 947 CALL ctl_stop( & 948 'Bdy data files must have the same time origin', & 949 'Multiple time frequencies not implemented yet' ) 950 ENDIF 951 zoffset = zoffsett 952 953 !! Check that times are the same in the three files... HERE. 954 istep_bt(:) = istept(:) 955 956 ! Check number of time dumps: 957 IF (ln_clim) THEN 958 SELECT CASE ( ntimes_bdy_bt ) 959 CASE( 1 ) 960 IF(lwp) WRITE(numout,*) 961 IF(lwp) WRITE(numout,*) 'We assume constant boundary forcing from bdy data files' 962 IF(lwp) WRITE(numout,*) 963 CASE( 12 ) 964 IF(lwp) WRITE(numout,*) 965 IF(lwp) WRITE(numout,*) 'We assume monthly (and cyclic) boundary forcing from bdy data files' 966 IF(lwp) WRITE(numout,*) 967 CASE DEFAULT 968 CALL ctl_stop( & 969 'For climatological boundary forcing (ln_clim=.true.),',& 970 'bdy data files must contain 1 or 12 time dumps.' ) 971 END SELECT 972 ENDIF 973 974 ! Find index of first record to read (before first model time). 975 976 it=1 977 DO WHILE ( ((istep_bt(it+1)) <= 0 ).AND.(it.LE.(ntimes_bdy_bt-1))) 978 it=it+1 979 END DO 980 nbdy_b_bt = it 981 982 IF(lwp) WRITE(numout,*) 'Time offset is ',zoffset 983 IF(lwp) WRITE(numout,*) 'First record to read is ',nbdy_b_bt 984 985 ENDIF ! endif (nn_dtactl == 1) 986 987 ! 1.2 Read first record in file if necessary (ie if nn_dtactl == 1) 988 ! ***************************************************************** 989 990 IF ( nn_dtactl == 0) THEN 991 ! boundary data arrays are filled with initial conditions 992 igrd = 5 ! U-points data 993 DO ib = 1, nblen(igrd) 994 ubtbdy(ib) = un(nbi(ib,igrd), nbj(ib,igrd), 1) 995 END DO 996 997 igrd = 6 ! V-points data 998 DO ib = 1, nblen(igrd) 999 vbtbdy(ib) = vn(nbi(ib,igrd), nbj(ib,igrd), 1) 1000 END DO 1001 1002 igrd = 4 ! T-points data 1003 DO ib = 1, nblen(igrd) 1004 sshbdy(ib) = sshn(nbi(ib,igrd), nbj(ib,igrd)) 1005 END DO 1006 1007 ELSEIF (nn_dtactl == 1) THEN 1008 1009 ! Set first record in the climatological case: 1010 IF ((ln_clim).AND.(ntimes_bdy_bt==1)) THEN 1011 nbdy_a_bt = 1 1012 ELSEIF ((ln_clim).AND.(ntimes_bdy_bt==iman)) THEN 1013 nbdy_b_bt = 0 1014 nbdy_a_bt = imois 1015 ELSE 1016 nbdy_a_bt = nbdy_b_bt 1017 END IF 1018 1019 ! Open Netcdf files: 1020 1021 CALL iom_open ( cn_dta_fla_T, numbdyt_bt ) 1022 CALL iom_open ( cn_dta_fla_U, numbdyu_bt ) 1023 CALL iom_open ( cn_dta_fla_V, numbdyv_bt ) 1024 1025 ! Read first record: 1026 ipj=1 1027 igrd=4 1028 ipi=nblendta(igrd) 1029 1030 ! ssh 1031 igrd=4 1032 IF ( nblendta(igrd) .le. 0 ) THEN 1033 idvar = iom_varid( numbdyt_bt,'sossheig' ) 1034 nblendta(igrd) = iom_file(numbdyt_bt)%dimsz(1,idvar) 1035 ENDIF 1036 WRITE(numout,*) 'Dim size for sossheig is ',nblendta(igrd) 1037 ipi=nblendta(igrd) 1038 1039 CALL iom_get ( numbdyt_bt, jpdom_unknown,'sossheig',zdta(1:ipi,1:ipj),nbdy_a_bt ) 1040 1041 DO ib=1, nblen(igrd) 1042 sshbdydta(ib,2) = zdta(nbmap(ib,igrd),1) 1043 END DO 1044 1045 ! u-velocity 1046 igrd=5 1047 IF ( nblendta(igrd) .le. 0 ) THEN 1048 idvar = iom_varid( numbdyu_bt,'vobtcrtx' ) 1049 nblendta(igrd) = iom_file(numbdyu_bt)%dimsz(1,idvar) 1050 ENDIF 1051 WRITE(numout,*) 'Dim size for vobtcrtx is ',nblendta(igrd) 1052 ipi=nblendta(igrd) 1053 1054 CALL iom_get ( numbdyu_bt, jpdom_unknown,'vobtcrtx',zdta(1:ipi,1:ipj),nbdy_a_bt ) 1055 1056 DO ib=1, nblen(igrd) 1057 ubtbdydta(ib,2) = zdta(nbmap(ib,igrd),1) 1058 END DO 1059 1060 ! v-velocity 1061 igrd=6 1062 IF ( nblendta(igrd) .le. 0 ) THEN 1063 idvar = iom_varid( numbdyv_bt,'vobtcrty' ) 1064 nblendta(igrd) = iom_file(numbdyv_bt)%dimsz(1,idvar) 1065 ENDIF 1066 WRITE(numout,*) 'Dim size for vobtcrty is ',nblendta(igrd) 1067 ipi=nblendta(igrd) 1068 1069 CALL iom_get ( numbdyv_bt, jpdom_unknown,'vobtcrty',zdta(1:ipi,1:ipj),nbdy_a_bt ) 1070 1071 DO ib=1, nblen(igrd) 1072 vbtbdydta(ib,2) = zdta(nbmap(ib,igrd),1) 1073 END DO 1074 1075 END IF 1076 1077 ! In the case of constant boundary forcing fill bdy arrays once for all 1078 IF ((ln_clim).AND.(ntimes_bdy_bt==1)) THEN 1079 1080 ubtbdy (:) = ubtbdydta (:,2) 1081 vbtbdy (:) = vbtbdydta (:,2) 1082 sshbdy (:) = sshbdydta (:,2) 1083 1084 CALL iom_close( numbdyt_bt ) 1085 CALL iom_close( numbdyu_bt ) 1086 CALL iom_close( numbdyv_bt ) 1087 1088 END IF 1089 1090 ENDIF ! End if nit000 1091 1092 ! -------------------- ! 1093 ! 2. At each time step ! 1094 ! -------------------- ! 1095 1096 IF ((nn_dtactl==1).AND.(ntimes_bdy_bt>1)) THEN 1097 1098 ! 2.1 Read one more record if necessary 1099 !************************************** 1100 1101 IF ( (ln_clim).AND.(imois/=nbdy_b_bt) ) THEN ! remember that nbdy_b_bt=0 for kt=nit000 1102 nbdy_b_bt = imois 1103 nbdy_a_bt = imois+1 1104 nbdy_b_bt = MOD( nbdy_b_bt, iman ) 1105 IF( nbdy_b_bt == 0 ) nbdy_b_bt = iman 1106 nbdy_a_bt = MOD( nbdy_a_bt, iman ) 1107 IF( nbdy_a_bt == 0 ) nbdy_a_bt = iman 1108 lect=.true. 1109 1110 ELSEIF ((.NOT.ln_clim).AND.(itimer >= istep_bt(nbdy_a_bt))) THEN 1111 nbdy_b_bt=nbdy_a_bt 1112 nbdy_a_bt=nbdy_a_bt+1 1113 lect=.true. 1114 END IF 1115 1116 IF (lect) THEN 1117 1118 ! Swap arrays 1119 sshbdydta(:,1) = sshbdydta(:,2) 1120 ubtbdydta(:,1) = ubtbdydta(:,2) 1121 vbtbdydta(:,1) = vbtbdydta(:,2) 1122 1123 ! read another set 1124 1125 ipj=1 1126 ipk=jpk 1127 igrd=4 1128 ipi=nblendta(igrd) 1129 1130 1131 ! ssh 1132 igrd=4 1133 ipi=nblendta(igrd) 1134 1135 CALL iom_get ( numbdyt_bt, jpdom_unknown,'sossheig',zdta(1:ipi,1:ipj),nbdy_a_bt ) 1136 1137 DO ib=1, nblen(igrd) 1138 sshbdydta(ib,2) = zdta(nbmap(ib,igrd),1) 1139 END DO 1140 1141 ! u-velocity 1142 igrd=5 1143 ipi=nblendta(igrd) 1144 1145 CALL iom_get ( numbdyu_bt, jpdom_unknown,'vobtcrtx',zdta(1:ipi,1:ipj),nbdy_a_bt ) 1146 1147 DO ib=1, nblen(igrd) 1148 ubtbdydta(ib,2) = zdta(nbmap(ib,igrd),1) 1149 END DO 1150 1151 ! v-velocity 1152 igrd=6 1153 ipi=nblendta(igrd) 1154 1155 CALL iom_get ( numbdyv_bt, jpdom_unknown,'vobtcrty',zdta(1:ipi,1:ipj),nbdy_a_bt ) 1156 1157 DO ib=1, nblen(igrd) 1158 vbtbdydta(ib,2) = zdta(nbmap(ib,igrd),1) 1159 END DO 1160 1161 1162 IF(lwp) WRITE(numout,*) 'bdy_dta_fla : first record file used nbdy_b_bt ',nbdy_b_bt 1163 IF(lwp) WRITE(numout,*) '~~~~~~~~ last record file used nbdy_a_bt ',nbdy_a_bt 1164 IF (.NOT.ln_clim) THEN 1165 IF(lwp) WRITE(numout,*) 'first record time (s): ', istep_bt(nbdy_b_bt) 1166 IF(lwp) WRITE(numout,*) 'model time (s) : ', itimer 1167 IF(lwp) WRITE(numout,*) 'second record time (s): ', istep_bt(nbdy_a_bt) 1168 ENDIF 1169 END IF ! end lect=.true. 1170 1171 1172 ! 2.2 Interpolate linearly: 1173 ! *************************** 1174 1175 IF (ln_clim) THEN 1176 zxy = REAL( nday, wp ) / REAL( nmonth_len(nbdy_b_bt), wp ) + 0.5 - i15 1177 ELSE 1178 zxy = REAL(istep_bt(nbdy_b_bt)-itimer, wp) / REAL(istep_bt(nbdy_b_bt)-istep_bt(nbdy_a_bt), wp) 1179 END IF 1180 1181 igrd=4 1182 DO ib=1, nblen(igrd) 1183 sshbdy(ib) = zxy * sshbdydta(ib,2) + & 1184 (1.-zxy) * sshbdydta(ib,1) 1185 END DO 1186 1187 igrd=5 1188 DO ib=1, nblen(igrd) 1189 ubtbdy(ib) = zxy * ubtbdydta(ib,2) + & 1190 (1.-zxy) * ubtbdydta(ib,1) 1191 END DO 1192 1193 igrd=6 1194 DO ib=1, nblen(igrd) 1195 vbtbdy(ib) = zxy * vbtbdydta(ib,2) + & 1196 (1.-zxy) * vbtbdydta(ib,1) 1197 END DO 1198 1199 1200 END IF !end if ((nn_dtactl==1).AND.(ntimes_bdy_bt>1)) 1201 1202 ! ------------------- ! 1203 ! Last call kt=nitend ! 1204 ! ------------------- ! 1205 1206 ! Closing of the 3 files 1207 IF( kt == nitend .and. jit == icycl ) THEN 1208 CALL iom_close( numbdyt_bt ) 1209 CALL iom_close( numbdyu_bt ) 1210 CALL iom_close( numbdyv_bt ) 1211 ENDIF 1212 1213 ENDIF ! ln_dyn_frs 1214 1215 END SUBROUTINE bdy_dta_fla 1216 607 IF( nn_dyn2d(ib_bdy) .ne. jp_frs ) THEN 608 jfld = jfld + 1 609 dta_bdy(ib_bdy)%ssh => bf(jfld)%fnow(:,1,1) 610 ENDIF 611 jfld = jfld + 1 612 dta_bdy(ib_bdy)%u2d => bf(jfld)%fnow(:,1,1) 613 jfld = jfld + 1 614 dta_bdy(ib_bdy)%v2d => bf(jfld)%fnow(:,1,1) 615 ENDIF 616 ENDIF 617 618 IF ( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 0 ) THEN 619 IF( nn_dyn3d(ib_bdy) .eq. jp_frs ) THEN 620 ilen0(1:3) = nblen(1:3) 621 ELSE 622 ilen0(1:3) = nblenrim(1:3) 623 ENDIF 624 ALLOCATE( dta_bdy(ib_bdy)%u3d(ilen0(2),jpk) ) 625 ALLOCATE( dta_bdy(ib_bdy)%v3d(ilen0(3),jpk) ) 626 ENDIF 627 IF ( ( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 1 ).or. & 628 & ( ln_full_vel_array(ib_bdy) .and. nn_dyn2d(ib_bdy) .gt. 0 .and. & 629 & ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) ) THEN 630 jfld = jfld + 1 631 dta_bdy(ib_bdy)%u3d => bf(jfld)%fnow(:,1,:) 632 jfld = jfld + 1 633 dta_bdy(ib_bdy)%v3d => bf(jfld)%fnow(:,1,:) 634 ENDIF 635 636 IF (nn_tra(ib_bdy) .gt. 0) THEN 637 IF( nn_tra_dta(ib_bdy) .eq. 0 ) THEN 638 IF( nn_tra(ib_bdy) .eq. jp_frs ) THEN 639 ilen0(1:3) = nblen(1:3) 640 ELSE 641 ilen0(1:3) = nblenrim(1:3) 642 ENDIF 643 ALLOCATE( dta_bdy(ib_bdy)%tem(ilen0(1),jpk) ) 644 ALLOCATE( dta_bdy(ib_bdy)%sal(ilen0(1),jpk) ) 645 ELSE 646 jfld = jfld + 1 647 dta_bdy(ib_bdy)%tem => bf(jfld)%fnow(:,1,:) 648 jfld = jfld + 1 649 dta_bdy(ib_bdy)%sal => bf(jfld)%fnow(:,1,:) 650 ENDIF 651 ENDIF 652 653 #if defined key_lim2 654 IF (nn_ice_lim2(ib_bdy) .gt. 0) THEN 655 IF( nn_ice_lim2_dta(ib_bdy) .eq. 0 ) THEN 656 IF( nn_ice_lim2(ib_bdy) .eq. jp_frs ) THEN 657 ilen0(1:3) = nblen(1:3) 658 ELSE 659 ilen0(1:3) = nblenrim(1:3) 660 ENDIF 661 ALLOCATE( dta_bdy(ib_bdy)%frld(ilen0(1)) ) 662 ALLOCATE( dta_bdy(ib_bdy)%hicif(ilen0(1)) ) 663 ALLOCATE( dta_bdy(ib_bdy)%hsnif(ilen0(1)) ) 664 ELSE 665 jfld = jfld + 1 666 dta_bdy(ib_bdy)%frld => bf(jfld)%fnow(:,1,1) 667 jfld = jfld + 1 668 dta_bdy(ib_bdy)%hicif => bf(jfld)%fnow(:,1,1) 669 jfld = jfld + 1 670 dta_bdy(ib_bdy)%hsnif => bf(jfld)%fnow(:,1,1) 671 ENDIF 672 ENDIF 673 #endif 674 675 ENDDO ! ib_bdy 676 677 IF( nn_timing == 1 ) CALL timing_stop('bdy_dta_init') 678 679 END SUBROUTINE bdy_dta_init 1217 680 1218 681 #else 1219 682 !!---------------------------------------------------------------------- 1220 !! Dummy module NO UnstructOpen Boundary Conditions683 !! Dummy module NO Open Boundary Conditions 1221 684 !!---------------------------------------------------------------------- 1222 685 CONTAINS 1223 SUBROUTINE bdy_dta_frs( kt ) ! Empty routine 1224 WRITE(*,*) 'bdy_dta_frs: You should not have seen this print! error?', kt 1225 END SUBROUTINE bdy_dta_frs 1226 SUBROUTINE bdy_dta_fla( kt, kit, icycle ) ! Empty routine 1227 WRITE(*,*) 'bdy_dta_frs: You should not have seen this print! error?', kt, kit 1228 END SUBROUTINE bdy_dta_fla 686 SUBROUTINE bdy_dta( kt, jit, time_offset ) ! Empty routine 687 INTEGER, INTENT( in ) :: kt 688 INTEGER, INTENT( in ), OPTIONAL :: jit 689 INTEGER, INTENT( in ), OPTIONAL :: time_offset 690 WRITE(*,*) 'bdy_dta: You should not have seen this print! error?', kt 691 END SUBROUTINE bdy_dta 692 SUBROUTINE bdy_dta_init() ! Empty routine 693 WRITE(*,*) 'bdy_dta_init: You should not have seen this print! error?' 694 END SUBROUTINE bdy_dta_init 1229 695 #endif 1230 696
Note: See TracChangeset
for help on using the changeset viewer.