Changeset 2865
- Timestamp:
- 2011-09-27T14:33:01+02:00 (13 years ago)
- Location:
- branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 13 edited
-
DYN/dynnxt.F90 (modified) (3 diffs)
-
DYN/dynspg_ts.F90 (modified) (1 diff)
-
OBC/obc_oce.F90 (modified) (2 diffs)
-
OBC/obcdta.F90 (modified) (12 diffs)
-
OBC/obcdyn.F90 (modified) (4 diffs)
-
OBC/obcdyn2d.F90 (modified) (3 diffs)
-
OBC/obcini.F90 (modified) (10 diffs)
-
OBC/obctides.F90 (modified) (5 diffs)
-
OBC/obctra.F90 (modified) (1 diff)
-
OBC/obcvol.F90 (modified) (1 diff)
-
SBC/fldread.F90 (modified) (5 diffs)
-
oce.F90 (modified) (1 diff)
-
step.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r2818 r2865 31 31 USE obc_oce ! ocean open boundary conditions 32 32 USE obcdta ! ocean open boundary conditions 33 USE obcdyn 3d! ocean open boundary conditions33 USE obcdyn ! ocean open boundary conditions 34 34 USE obcvol ! ocean open boundary condition (obc_vol routines) 35 35 USE in_out_manager ! I/O manager … … 153 153 # if defined key_obc 154 154 ! !* OBC open boundaries 155 IF( .NOT. lk_dynspg_flt ) THEN 156 157 CALL obc_dyn3d( kt ) 155 IF( lk_dynspg_exp ) CALL obc_dyn( kt ) 156 IF( lk_dynspg_ts ) CALL obc_dyn( kt, dyn3d_only=.true. ) 158 157 159 158 !!$!!gm ERROR - potential BUG: sshn should not be modified at this stage !! ssh_nxt not alrady called … … 163 162 !!$ ! 164 163 !!$ IF(ln_ctl) CALL prt_ctl( tab2d_1=sshn, clinfo1=' ssh : ', mask1=tmask ) 165 166 ENDIF167 164 ! 168 165 # endif -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r2814 r2865 355 355 ! !* Update the forcing (OBC and tides) 356 356 ! ! ------------------ 357 IF( lk_obc ) CALL obc_dta ( kt, jit=jn 357 IF( lk_obc ) CALL obc_dta ( kt, jit=jn, time_offset=+1 ) 358 358 359 359 ! !* after ssh_e -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obc_oce.F90
r2831 r2865 57 57 LOGICAL :: ln_mask_file !: =T read obcmask from file 58 58 LOGICAL :: ln_vol !: =T volume correction 59 LOGICAL, DIMENSION(jp_obc) :: ln_clim !: =T obc data files contain climatological data (time-cyclic)60 59 ! 61 60 INTEGER :: nb_obc !: number of open boundary sets 62 INTEGER, DIMENSION(jp_obc) :: nn_rimwidth !: boundary rim width 63 INTEGER, DIMENSION(jp_obc) :: nn_dtactl !: = 0 use the initial state as obc dta ; 64 !: = 1 read it in a NetCDF file 65 INTEGER, DIMENSION(jp_obc) :: nn_tides !: = 0 no tidal harmonic forcing 66 !: = 1 apply ONLY tidal harmonic forcing for barotropic solution 67 !: = 2 ADD tidal harmonic forcing to other barotropic boundary data 61 INTEGER, DIMENSION(jp_obc) :: nn_rimwidth !: boundary rim width for Flow Relaxation Scheme 68 62 INTEGER :: nn_volctl !: = 0 the total volume will have the variability of the surface Flux E-P 69 63 ! ! = 1 the volume will be constant during all the integration. … … 71 65 INTEGER, DIMENSION(jp_obc) :: nn_dyn2d_dta !: = 0 use the initial state as obc dta ; 72 66 !: = 1 read it in a NetCDF file 67 !: = 2 read tidal harmonic forcing from a NetCDF file 68 !: = 3 read external data AND tidal harmonic forcing from NetCDF files 73 69 INTEGER, DIMENSION(jp_obc) :: nn_dyn3d ! Choice of boundary condition for baroclinic velocities 74 70 INTEGER, DIMENSION(jp_obc) :: nn_dyn3d_dta !: = 0 use the initial state as obc dta ; -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcdta.F90
r2831 r2865 55 55 CONTAINS 56 56 57 SUBROUTINE obc_dta( kt, jit )57 SUBROUTINE obc_dta( kt, jit, time_offset ) 58 58 !!---------------------------------------------------------------------- 59 59 !! *** SUBROUTINE obc_dta *** … … 69 69 INTEGER, INTENT( in ) :: kt ! ocean time-step index 70 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. 71 77 !! 72 78 INTEGER :: ib_obc, jfld, jstart, jend, ib, ii, ij, ik, igrd ! local indices … … 202 208 IF( PRESENT(jit) ) THEN 203 209 ! Update barotropic boundary conditions only 204 ! jit is optional argument for fld_read 205 IF( nn_dyn2d(ib_obc) .gt. 0 .and. nn_dyn2d_dta(ib_obc) .eq. 1) THEN206 IF( nn_ tides(ib_obc) .eq. 1 ) THEN210 ! jit is optional argument for fld_read and tide_update 211 IF( nn_dyn2d(ib_obc) .gt. 0 ) THEN 212 IF( nn_dyn2d_dta(ib_obc) .eq. 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 207 213 dta_obc(ib_obc)%ssh(:) = 0.0 208 214 dta_obc(ib_obc)%u2d(:) = 0.0 209 215 dta_obc(ib_obc)%v2d(:) = 0.0 210 ELSE 216 ENDIF 217 IF( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn2d_dta(ib_obc) .eq. 3 ) THEN ! update external data 211 218 jend = jstart + 2 212 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), jit=jit )219 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), jit=jit, time_offset=time_offset ) 213 220 ENDIF 214 ENDIF215 IF( nn_tides(ib_obc) .gt. 0 ) THEN216 CALL tide_update( kt=kt, jit=jit, idx=idx_obc(ib_obc), dta=dta_obc(ib_obc), td=tides(ib_obc) )221 IF( nn_dyn2d_dta(ib_obc) .ge. 2 ) THEN ! update tidal harmonic forcing 222 CALL tide_update( kt=kt, idx=idx_obc(ib_obc), dta=dta_obc(ib_obc), td=tides(ib_obc), jit=jit, time_offset=time_offset ) 223 ENDIF 217 224 ENDIF 218 225 ELSE 219 IF( nb_obc_fld(ib_obc) .gt. 0 ) THEN 220 jend = jstart + nb_obc_fld(ib_obc) - 1 221 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), timeshift=1 ) 222 ENDIF 223 IF( nn_tides(ib_obc) .eq. 1 ) THEN 226 IF( nn_dyn2d(ib_obc) .gt. 0 .and. nn_dyn2d_dta(ib_obc) .eq. 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 224 227 dta_obc(ib_obc)%ssh(:) = 0.0 225 228 dta_obc(ib_obc)%u2d(:) = 0.0 226 229 dta_obc(ib_obc)%v2d(:) = 0.0 227 230 ENDIF 228 IF( nn_tides(ib_obc) .gt. 0 ) THEN 229 !!! THINK ABOUT kt, jit VALUES !!! 230 CALL tide_update( kt=kt, jit=0, idx=idx_obc(ib_obc), dta=dta_obc(ib_obc), td=tides(ib_obc) ) 231 IF( nb_obc_fld(ib_obc) .gt. 0 ) THEN ! update external data 232 jend = jstart + nb_obc_fld(ib_obc) - 1 233 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), time_offset=time_offset ) 234 ENDIF 235 IF( nn_dyn2d(ib_obc) .gt. 0 .and. nn_dyn2d_dta(ib_obc) .ge. 2 ) THEN ! update tidal harmonic forcing 236 CALL tide_update( kt=kt, idx=idx_obc(ib_obc), dta=dta_obc(ib_obc), td=tides(ib_obc), time_offset=time_offset ) 231 237 ENDIF 232 238 ENDIF … … 238 244 239 245 IF( ln_full_vel_array(ib_obc) .and. & 240 & ( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn 3d_dta(ib_obc) .eq. 1 ) ) THEN246 & ( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn2d_dta(ib_obc) .eq. 3 .or. nn_dyn3d_dta(ib_obc) .eq. 1 ) ) THEN 241 247 242 248 igrd = 2 ! zonal velocity … … 326 332 #endif 327 333 ) 334 IF(nn_dta(ib_obc) .gt. 1) nn_dta(ib_obc) = 1 328 335 END DO 329 336 … … 333 340 nb_obc_fld(:) = 0 334 341 DO ib_obc = 1, nb_obc 335 IF( nn_dyn2d(ib_obc) .gt. 0 .and. nn_dyn2d_dta(ib_obc) .eq. 1) THEN342 IF( nn_dyn2d(ib_obc) .gt. 0 .and. ( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn2d_dta(ib_obc) .eq. 3 ) ) THEN 336 343 nb_obc_fld(ib_obc) = nb_obc_fld(ib_obc) + 3 337 344 ENDIF … … 408 415 ! Only read in necessary fields for this set. 409 416 ! Important that barotropic variables come first. 410 IF( nn_dyn2d(ib_obc) .gt. 0 .and. nn_dyn2d_dta(ib_obc) .eq. 1) THEN411 412 IF( nn_dyn2d(ib_obc) .ne. jp_frs .and. nn_tides(ib_obc) .ne. 1) THEN417 IF( nn_dyn2d(ib_obc) .gt. 0 .and. ( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn2d_dta(ib_obc) .eq. 3 ) ) THEN 418 419 IF( nn_dyn2d(ib_obc) .ne. jp_frs ) THEN 413 420 jfld = jfld + 1 414 421 blf_i(jfld) = bn_ssh … … 419 426 ENDIF 420 427 421 IF( .not. ln_full_vel_array(ib_obc) .and. nn_tides(ib_obc) .ne. 1) THEN428 IF( .not. ln_full_vel_array(ib_obc) ) THEN 422 429 423 430 jfld = jfld + 1 … … 449 456 ! baroclinic velocities 450 457 IF( ( nn_dyn3d(ib_obc) .gt. 0 .and. nn_dyn3d_dta(ib_obc) .eq. 1 ) .or. & 451 ( ln_full_vel_array(ib_obc) .and. nn_dyn2d(ib_obc) .gt. 0 .and. nn_dyn2d_dta(ib_obc) .eq. 1 ) ) THEN 458 & ( ln_full_vel_array(ib_obc) .and. nn_dyn2d(ib_obc) .gt. 0 .and. & 459 & ( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn2d_dta(ib_obc) .eq. 3 ) ) ) THEN 452 460 453 461 jfld = jfld + 1 … … 583 591 584 592 IF (nn_dyn2d(ib_obc) .gt. 0) THEN 585 IF( nn_dyn2d_dta(ib_obc) .eq. 0 .or. ln_full_vel_array(ib_obc) .or. nn_tides(ib_obc) .eq. 1) THEN593 IF( nn_dyn2d_dta(ib_obc) .eq. 0 .or. nn_dyn2d_dta(ib_obc) .eq. 2 .or. ln_full_vel_array(ib_obc) ) THEN 586 594 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 587 595 ilen0(1:3) = nblen(1:3) … … 614 622 ENDIF 615 623 IF ( ( nn_dyn3d(ib_obc) .gt. 0 .and. nn_dyn3d_dta(ib_obc) .eq. 1 ).or. & 616 & ( ln_full_vel_array(ib_obc) .and. nn_dyn2d(ib_obc) .gt. 0 ) ) THEN 624 & ( ln_full_vel_array(ib_obc) .and. nn_dyn2d(ib_obc) .gt. 0 .and. & 625 & ( nn_dyn2d_dta(ib_obc) .eq. 1 .or. nn_dyn2d_dta(ib_obc) .eq. 3 ) ) ) THEN 617 626 jfld = jfld + 1 618 627 dta_obc(ib_obc)%u3d => bf(jfld)%fnow(:,1,:) … … 646 655 ilen0(1:3) = nblenrim(1:3) 647 656 ENDIF 648 ALLOCATE( dta_obc(ib_obc)% ssh(ilen0(1)) )649 ALLOCATE( dta_obc(ib_obc)% u2d(ilen0(1)) )650 ALLOCATE( dta_obc(ib_obc)% v2d(ilen0(1)) )657 ALLOCATE( dta_obc(ib_obc)%frld(ilen0(1)) ) 658 ALLOCATE( dta_obc(ib_obc)%hicif(ilen0(1)) ) 659 ALLOCATE( dta_obc(ib_obc)%hsnif(ilen0(1)) ) 651 660 ELSE 652 661 jfld = jfld + 1 -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcdyn.F90
r2800 r2865 41 41 CONTAINS 42 42 43 SUBROUTINE obc_dyn( kt )43 SUBROUTINE obc_dyn( kt, dyn3d_only ) 44 44 !!---------------------------------------------------------------------- 45 45 !! *** SUBROUTINE obc_dyn *** … … 51 51 USE wrk_nemo, ONLY: wrk_2d_7, wrk_2d_8 ! 2D workspace 52 52 !! 53 INTEGER, INTENT( in ) :: kt ! Main time step counter 53 INTEGER, INTENT( in ) :: kt ! Main time step counter 54 LOGICAL, INTENT( in ), OPTIONAL :: dyn3d_only ! T => only update baroclinic velocities 54 55 !! 55 56 INTEGER :: jk,ii,ij,ib,igrd ! Loop counter 57 LOGICAL :: ll_dyn2d, ll_dyn3d 56 58 !! 57 59 … … 60 62 END IF 61 63 64 ll_dyn2d = .true. 65 ll_dyn3d = .true. 66 67 IF( PRESENT(dyn3d_only) ) THEN 68 IF( dyn3d_only ) ll_dyn2d = .false. 69 ENDIF 70 62 71 !------------------------------------------------------- 63 72 ! Set pointers 64 73 !------------------------------------------------------- 65 74 66 pssh => sshn _b75 pssh => sshn 67 76 phur => hur 68 77 phvr => hvr … … 92 101 !------------------------------------------------------- 93 102 94 CALL obc_dyn2d( kt )103 IF( ll_dyn2d ) CALL obc_dyn2d( kt ) 95 104 96 CALL obc_dyn3d( kt )105 IF( ll_dyn3d ) CALL obc_dyn3d( kt ) 97 106 98 107 !------------------------------------------------------- -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcdyn2d.F90
r2818 r2865 55 55 CYCLE 56 56 CASE(jp_frs) 57 CALL obc_dyn2d_frs( idx_obc(ib_obc), dta_obc(ib_obc) , kt)57 CALL obc_dyn2d_frs( idx_obc(ib_obc), dta_obc(ib_obc) ) 58 58 CASE(jp_flather) 59 59 CALL obc_dyn2d_fla( idx_obc(ib_obc), dta_obc(ib_obc) ) … … 65 65 END SUBROUTINE obc_dyn2d 66 66 67 SUBROUTINE obc_dyn2d_frs( idx, dta , kt)67 SUBROUTINE obc_dyn2d_frs( idx, dta ) 68 68 !!---------------------------------------------------------------------- 69 69 !! *** SUBROUTINE obc_dyn2d_frs *** … … 76 76 !! topography. Tellus, 365-382. 77 77 !!---------------------------------------------------------------------- 78 INTEGER, INTENT(in) :: kt79 78 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 80 79 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcini.F90
r2831 r2865 51 51 !! ** Input : obc_init.nc, input file for unstructured open boundaries 52 52 !!---------------------------------------------------------------------- 53 INTEGER :: ib_obc, ii, ij, ik, igrd, ib, ir ! dummy loop indices 53 ! namelist variables 54 !------------------- 55 INTEGER, PARAMETER :: jp_nseg = 100 56 INTEGER :: nobcsege, nobcsegw, nobcsegn, nobcsegs 57 INTEGER, DIMENSION(jp_nseg) :: jpieob, jpjedt, jpjeft 58 INTEGER, DIMENSION(jp_nseg) :: jpiwob, jpjwdt, jpjwft 59 INTEGER, DIMENSION(jp_nseg) :: jpjnob, jpindt, jpinft 60 INTEGER, DIMENSION(jp_nseg) :: jpjsob, jpisdt, jpisft 61 62 ! local variables 63 !------------------- 64 INTEGER :: ib_obc, ii, ij, ik, igrd, ib, ir, iseg ! dummy loop indices 54 65 INTEGER :: icount, icountr, ibr_max, ilen1, ibm1 ! local integers 55 66 INTEGER :: iw, ie, is, in, inum, id_dummy ! - - … … 72 83 & nn_ice_lim2, nn_ice_lim2_dta, & 73 84 #endif 74 & nn_tides, ln_vol, ln_clim, nn_volctl,&85 & ln_vol, nn_volctl, & 75 86 & nn_rimwidth, nn_dmp2d_in, nn_dmp2d_out, & 76 87 & nn_dmp3d_in, nn_dmp3d_out 88 !! 89 NAMELIST/namobc_index/ nobcsege, jpieob, jpjedt, jpjeft, & 90 nobcsegw, jpiwob, jpjwdt, jpjwft, & 91 nobcsegn, jpjnob, jpindt, jpinft, & 92 nobcsegs, jpjsob, jpisdt, jpisft 93 77 94 !!---------------------------------------------------------------------- 78 95 … … 109 126 #endif 110 127 ln_vol = .false. 111 ln_clim(:) = .false.112 nn_tides(:) = 0 ! default to no tidal forcing113 128 nn_volctl = -1 ! uninitialised flag 114 129 nn_rimwidth(:) = -1 ! uninitialised flag … … 137 152 IF(lwp) WRITE(numout,*) ' ' 138 153 IF(lwp) WRITE(numout,*) '------ Open boundary data set ',ib_obc,'------' 154 155 IF( ln_coords_file(ib_obc) ) THEN 156 IF(lwp) WRITE(numout,*) 'Boundary definition read from file '//TRIM(cn_coords_file(ib_obc)) 157 ELSE 158 IF(lwp) WRITE(numout,*) 'Boundary defined in namelist.' 159 ENDIF 160 IF(lwp) WRITE(numout,*) 139 161 140 162 IF(lwp) WRITE(numout,*) 'Boundary conditions for barotropic solution: ' … … 149 171 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' initial state used for obc data' 150 172 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' boundary data taken from file' 151 CASE DEFAULT ; CALL ctl_stop( 'nn_dyn2d_dta must be 0 or 1' ) 173 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' tidal harmonic forcing taken from file' 174 CASE( 3 ) ; IF(lwp) WRITE(numout,*) ' boundary data AND tidal harmonic forcing taken from files' 175 CASE DEFAULT ; CALL ctl_stop( 'nn_dyn2d_dta must be between 0 and 3' ) 152 176 END SELECT 153 177 ENDIF … … 201 225 #endif 202 226 203 IF(lwp) WRITE(numout,*) 'Boundary rim width for the FRS nn_rimwidth = ', nn_rimwidth227 IF(lwp) WRITE(numout,*) 'Boundary rim width for the FRS scheme = ', nn_rimwidth(ib_obc) 204 228 IF(lwp) WRITE(numout,*) 205 206 SELECT CASE( nn_tides(ib_obc) )207 CASE(0)208 IF(lwp) WRITE(numout,*) 'No tidal harmonic forcing'209 IF(lwp) WRITE(numout,*)210 CASE(1)211 IF(lwp) WRITE(numout,*) 'Tidal harmonic forcing ONLY for barotropic solution'212 IF(lwp) WRITE(numout,*)213 CASE(2)214 IF(lwp) WRITE(numout,*) 'Tidal harmonic forcing ADDED to other barotropic boundary conditions'215 IF(lwp) WRITE(numout,*)216 CASE DEFAULT217 CALL ctl_stop( 'obc_ini: ERROR: incorrect value for nn_tides ' )218 END SELECT219 229 220 230 ENDDO … … 240 250 ! Work out global dimensions of boundary data 241 251 ! --------------------------------------------- 252 REWIND( numnam ) 242 253 DO ib_obc = 1, nb_obc 243 254 … … 245 256 IF( .NOT. ln_coords_file(ib_obc) ) THEN ! Work out size of global arrays from namelist parameters 246 257 247 248 !! 1. Read parameters from namelist 249 !! 2. Work out global size of boundary data arrays nblendta and jpbdta 258 ! No REWIND here because may need to read more than one namobc_index namelist. 259 READ ( numnam, namobc_index ) 260 261 ! Automatic boundary definition: if nobcsegX = -1 262 ! set boundary to whole side of model domain. 263 IF( nobcsege == -1 ) THEN 264 nobcsege = 1 265 jpieob(1) = jpiglo - 1 266 jpjedt(1) = 2 267 jpjeft(1) = jpjglo - 1 268 ENDIF 269 IF( nobcsegw == -1 ) THEN 270 nobcsegw = 1 271 jpiwob(1) = 2 272 jpjwdt(1) = 2 273 jpjwft(1) = jpjglo - 1 274 ENDIF 275 IF( nobcsegn == -1 ) THEN 276 nobcsegn = 1 277 jpjnob(1) = jpjglo - 1 278 jpindt(1) = 2 279 jpinft(1) = jpiglo - 1 280 ENDIF 281 IF( nobcsegs == -1 ) THEN 282 nobcsegs = 1 283 jpjsob(1) = 2 284 jpisdt(1) = 2 285 jpisft(1) = jpiglo - 1 286 ENDIF 287 288 nblendta(:,ib_obc) = 0 289 DO iseg = 1, nobcsege 290 igrd = 1 291 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpjeft(iseg) - jpjedt(iseg) + 1 292 igrd = 2 293 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpjeft(iseg) - jpjedt(iseg) + 1 294 igrd = 3 295 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpjeft(iseg) - jpjedt(iseg) 296 ENDDO 297 DO iseg = 1, nobcsegw 298 igrd = 1 299 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpjwft(iseg) - jpjwdt(iseg) + 1 300 igrd = 2 301 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpjwft(iseg) - jpjwdt(iseg) + 1 302 igrd = 3 303 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpjwft(iseg) - jpjwdt(iseg) 304 ENDDO 305 DO iseg = 1, nobcsegn 306 igrd = 1 307 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpinft(iseg) - jpindt(iseg) + 1 308 igrd = 2 309 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpinft(iseg) - jpindt(iseg) 310 igrd = 3 311 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpinft(iseg) - jpindt(iseg) + 1 312 ENDDO 313 DO iseg = 1, nobcsegs 314 igrd = 1 315 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpisft(iseg) - jpisdt(iseg) + 1 316 igrd = 2 317 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpisft(iseg) - jpisdt(iseg) 318 igrd = 3 319 nblendta(igrd,ib_obc) = nblendta(igrd,ib_obc) + jpisft(iseg) - jpisdt(iseg) + 1 320 ENDDO 321 322 nblendta(:,ib_obc) = nblendta(:,ib_obc) * nn_rimwidth(ib_obc) 323 jpbdta = MAXVAL(nblendta(:,ib_obc)) 250 324 251 325 … … 263 337 ENDIF 264 338 265 ENDDO 339 ENDDO ! ib_obc 266 340 267 341 ! Allocate arrays … … 274 348 ! Calculate global boundary index arrays or read in from file 275 349 !------------------------------------------------------------ 350 REWIND( numnam ) 276 351 DO ib_obc = 1, nb_obc 277 352 278 353 IF( .NOT. ln_coords_file(ib_obc) ) THEN ! Calculate global index arrays from namelist parameters 279 354 280 !! Calculate global index arrays from namelist parameters 355 ! No REWIND here because may need to read more than one namobc_index namelist. 356 READ ( numnam, namobc_index ) 357 358 ! Automatic boundary definition: if nobcsegX = -1 359 ! set boundary to whole side of model domain. 360 IF( nobcsege == -1 ) THEN 361 nobcsege = 1 362 jpieob(1) = jpiglo - 1 363 jpjedt(1) = 2 364 jpjeft(1) = jpjglo - 1 365 ENDIF 366 IF( nobcsegw == -1 ) THEN 367 nobcsegw = 1 368 jpiwob(1) = 2 369 jpjwdt(1) = 2 370 jpjwft(1) = jpjglo - 1 371 ENDIF 372 IF( nobcsegn == -1 ) THEN 373 nobcsegn = 1 374 jpjnob(1) = jpjglo - 1 375 jpindt(1) = 2 376 jpinft(1) = jpiglo - 1 377 ENDIF 378 IF( nobcsegs == -1 ) THEN 379 nobcsegs = 1 380 jpjsob(1) = 2 381 jpisdt(1) = 2 382 jpisft(1) = jpiglo - 1 383 ENDIF 384 385 ! ------------ T points ------------- 386 igrd = 1 387 icount = 0 388 DO ir = 1, nn_rimwidth(ib_obc) 389 ! east 390 DO iseg = 1, nobcsege 391 DO ij = jpjedt(iseg), jpjeft(iseg) 392 icount = icount + 1 393 nbidta(icount, igrd, ib_obc) = jpieob(iseg) - ir + 1 394 nbjdta(icount, igrd, ib_obc) = ij 395 nbrdta(icount, igrd, ib_obc) = ir 396 ENDDO 397 ENDDO 398 ! west 399 DO iseg = 1, nobcsegw 400 DO ij = jpjwdt(iseg), jpjwft(iseg) 401 icount = icount + 1 402 nbidta(icount, igrd, ib_obc) = jpiwob(iseg) + ir - 1 403 nbjdta(icount, igrd, ib_obc) = ij 404 nbrdta(icount, igrd, ib_obc) = ir 405 ENDDO 406 ENDDO 407 ! north 408 DO iseg = 1, nobcsegn 409 DO ii = jpindt(iseg), jpinft(iseg) 410 icount = icount + 1 411 nbidta(icount, igrd, ib_obc) = ii 412 nbjdta(icount, igrd, ib_obc) = jpjnob(iseg) - ir + 1 413 nbrdta(icount, igrd, ib_obc) = ir 414 ENDDO 415 ENDDO 416 ! south 417 DO iseg = 1, nobcsegs 418 DO ii = jpisdt(iseg), jpisft(iseg) 419 icount = icount + 1 420 nbidta(icount, igrd, ib_obc) = ii 421 nbjdta(icount, igrd, ib_obc) = jpjsob(iseg) + ir + 1 422 nbrdta(icount, igrd, ib_obc) = ir 423 ENDDO 424 ENDDO 425 ENDDO 426 427 ! ------------ U points ------------- 428 igrd = 2 429 icount = 0 430 DO ir = 1, nn_rimwidth(ib_obc) 431 ! east 432 DO iseg = 1, nobcsege 433 DO ij = jpjedt(iseg), jpjeft(iseg) 434 icount = icount + 1 435 nbidta(icount, igrd, ib_obc) = jpieob(iseg) - ir 436 nbjdta(icount, igrd, ib_obc) = ij 437 nbrdta(icount, igrd, ib_obc) = ir 438 ENDDO 439 ENDDO 440 ! west 441 DO iseg = 1, nobcsegw 442 DO ij = jpjwdt(iseg), jpjwft(iseg) 443 icount = icount + 1 444 nbidta(icount, igrd, ib_obc) = jpiwob(iseg) + ir - 1 445 nbjdta(icount, igrd, ib_obc) = ij 446 nbrdta(icount, igrd, ib_obc) = ir 447 ENDDO 448 ENDDO 449 ! north 450 DO iseg = 1, nobcsegn 451 DO ii = jpindt(iseg), jpinft(iseg) - 1 452 icount = icount + 1 453 nbidta(icount, igrd, ib_obc) = ii 454 nbjdta(icount, igrd, ib_obc) = jpjnob(iseg) - ir + 1 455 nbrdta(icount, igrd, ib_obc) = ir 456 ENDDO 457 ENDDO 458 ! south 459 DO iseg = 1, nobcsegs 460 DO ii = jpisdt(iseg), jpisft(iseg) - 1 461 icount = icount + 1 462 nbidta(icount, igrd, ib_obc) = ii 463 nbjdta(icount, igrd, ib_obc) = jpjsob(iseg) + ir + 1 464 nbrdta(icount, igrd, ib_obc) = ir 465 ENDDO 466 ENDDO 467 ENDDO 468 469 ! ------------ V points ------------- 470 igrd = 3 471 icount = 0 472 DO ir = 1, nn_rimwidth(ib_obc) 473 ! east 474 DO iseg = 1, nobcsege 475 DO ij = jpjedt(iseg), jpjeft(iseg) - 1 476 icount = icount + 1 477 nbidta(icount, igrd, ib_obc) = jpieob(iseg) - ir + 1 478 nbjdta(icount, igrd, ib_obc) = ij 479 nbrdta(icount, igrd, ib_obc) = ir 480 ENDDO 481 ENDDO 482 ! west 483 DO iseg = 1, nobcsegw 484 DO ij = jpjwdt(iseg), jpjwft(iseg) - 1 485 icount = icount + 1 486 nbidta(icount, igrd, ib_obc) = jpiwob(iseg) + ir - 1 487 nbjdta(icount, igrd, ib_obc) = ij 488 nbrdta(icount, igrd, ib_obc) = ir 489 ENDDO 490 ENDDO 491 ! north 492 DO iseg = 1, nobcsegn 493 DO ii = jpindt(iseg), jpinft(iseg) 494 icount = icount + 1 495 nbidta(icount, igrd, ib_obc) = ii 496 nbjdta(icount, igrd, ib_obc) = jpjnob(iseg) - ir 497 nbrdta(icount, igrd, ib_obc) = ir 498 ENDDO 499 ENDDO 500 ! south 501 DO iseg = 1, nobcsegs 502 DO ii = jpisdt(iseg), jpisft(iseg) 503 icount = icount + 1 504 nbidta(icount, igrd, ib_obc) = ii 505 nbjdta(icount, igrd, ib_obc) = jpjsob(iseg) + ir + 1 506 nbrdta(icount, igrd, ib_obc) = ir 507 ENDDO 508 ENDDO 509 ENDDO 281 510 282 511 ELSE ! Read global index arrays from boundary coordinates file. -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obctides.F90
r2814 r2865 94 94 REWIND(numnam) 95 95 DO ib_obc = 1, nb_obc 96 IF( nn_ tides(ib_obc) .gt. 0) THEN96 IF( nn_dyn2d_dta(ib_obc) .ge. 2 ) THEN 97 97 98 98 td => tides(ib_obc) … … 260 260 ENDIF ! date correction 261 261 ! 262 ENDIF ! nn_ tides .gt. 0262 ENDIF ! nn_dyn2d_dta(ib_obc) .ge. 2 263 263 ! 264 264 END DO ! loop on ib_obc … … 267 267 268 268 269 SUBROUTINE tide_update ( kt, jit, idx, dta, td)269 SUBROUTINE tide_update ( kt, idx, dta, td, jit, time_offset ) 270 270 !!---------------------------------------------------------------------- 271 271 !! *** SUBROUTINE tide_update *** … … 276 276 INTEGER, INTENT( in ) :: kt ! Main timestep counter 277 277 !!gm doctor jit ==> kit 278 INTEGER, INTENT( in ) :: jit ! Barotropic timestep counter (for timesplitting option)279 278 TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices 280 279 TYPE(OBC_DATA), INTENT(inout) :: dta ! OBC external data 281 280 TYPE(TIDES_DATA),INTENT( in ) :: td ! tidal harmonics data 281 INTEGER,INTENT(in),OPTIONAL :: jit ! Barotropic timestep counter (for timesplitting option) 282 INTEGER,INTENT( in ), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit 283 ! is present then units = subcycle timesteps. 284 ! time_offset = 0 => get data at "now" time level 285 ! time_offset = -1 => get data at "before" time level 286 ! time_offset = +1 => get data at "after" time level 287 ! etc. 282 288 !! 283 289 INTEGER :: itide, igrd, ib ! dummy loop indices 290 INTEGER :: time_add ! time offset in units of timesteps 284 291 REAL(wp) :: z_arg, z_sarg 285 292 REAL(wp), DIMENSION(jptides_max) :: z_sist, z_cost 286 293 !!---------------------------------------------------------------------- 287 294 295 time_add = 0 296 IF( PRESENT(time_offset) ) THEN 297 time_add = time_offset 298 ENDIF 299 288 300 ! Note tide phase speeds are in deg/hour, so we need to convert the 289 301 ! elapsed time in seconds to hours by dividing by 3600.0 290 IF( jit == 0 ) THEN 291 z_arg = kt * rdt * rad / 3600.0 292 ELSE ! we are in a barotropic subcycle (for timesplitting option) 293 ! z_arg = ( (kt-1) * rdt + jit * rdt / REAL(nn_baro,lwp) ) * rad / 3600.0 294 z_arg = ( (kt-1) * rdt + jit * rdt / REAL(nn_baro,wp) ) * rad / 3600.0 302 IF( PRESENT(jit) ) THEN 303 z_arg = ( (kt-1) * rdt + (jit+time_add) * rdt / REAL(nn_baro,wp) ) * rad / 3600.0 304 ELSE 305 z_arg = (kt+time_add) * rdt * rad / 3600.0 295 306 ENDIF 296 307 … … 301 312 END DO 302 313 303 !304 314 DO itide = 1, td%ncpt 305 315 igrd=1 ! SSH on tracer grid. -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obctra.F90
r2797 r2865 52 52 CALL obc_tra_frs( idx_obc(ib_obc), dta_obc(ib_obc), kt ) 53 53 CASE DEFAULT 54 CALL ctl_stop( 'obc_tra : unrecognised option for open boundaries for T an S' )54 CALL ctl_stop( 'obc_tra : unrecognised option for open boundaries for T and S' ) 55 55 END SELECT 56 56 ENDDO -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcvol.F90
r2797 r2865 10 10 !!---------------------------------------------------------------------- 11 11 #if defined key_obc && defined key_dynspg_flt 12 !!$ !!---------------------------------------------------------------------- 13 !!$ !! 'key_obc' AND unstructured open boundary conditions 14 !!$ !! 'key_dynspg_flt' filtered free surface 15 !!$ !!---------------------------------------------------------------------- 16 !!$ USE oce ! ocean dynamics and tracers 17 !!$ USE dom_oce ! ocean space and time domain 18 !!$ USE phycst ! physical constants 19 !!$ USE obc_oce ! ocean open boundary conditions 20 !!$ USE lib_mpp ! for mppsum 21 !!$ USE in_out_manager ! I/O manager 22 !!$ USE sbc_oce ! ocean surface boundary conditions 23 !!$ 24 !!$ IMPLICIT NONE 25 !!$ PRIVATE 26 !!$ 27 !!$ PUBLIC obc_vol ! routine called by dynspg_flt.h90 28 !!$ 29 !!$ !! * Substitutions 30 !!$# include "domzgr_substitute.h90" 31 !!$ !!---------------------------------------------------------------------- 32 !!$ !! NEMO/OPA 3.3 , NEMO Consortium (2010) 33 !!$ !! $Id$ 34 !!$ !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 35 !!$ !!---------------------------------------------------------------------- 36 !!$CONTAINS 37 !!$ 38 !!$ SUBROUTINE obc_vol( kt ) 39 !!$ !!---------------------------------------------------------------------- 40 !!$ !! *** ROUTINE obcvol *** 41 !!$ !! 42 !!$ !! ** Purpose : This routine is called in dynspg_flt to control 43 !!$ !! the volume of the system. A correction velocity is calculated 44 !!$ !! to correct the total transport through the unstructured OBC. 45 !!$ !! The total depth used is constant (H0) to be consistent with the 46 !!$ !! linear free surface coded in OPA 8.2 47 !!$ !! 48 !!$ !! ** Method : The correction velocity (zubtpecor here) is defined calculating 49 !!$ !! the total transport through all open boundaries (trans_obc) minus 50 !!$ !! the cumulate E-P flux (z_cflxemp) divided by the total lateral 51 !!$ !! surface (obcsurftot) of the unstructured boundary. 52 !!$ !! zubtpecor = [trans_obc - z_cflxemp ]*(1./obcsurftot) 53 !!$ !! with z_cflxemp => sum of (Evaporation minus Precipitation) 54 !!$ !! over all the domain in m3/s at each time step. 55 !!$ !! z_cflxemp < 0 when precipitation dominate 56 !!$ !! z_cflxemp > 0 when evaporation dominate 57 !!$ !! 58 !!$ !! There are 2 options (user's desiderata): 59 !!$ !! 1/ The volume changes according to E-P, this is the default 60 !!$ !! option. In this case the cumulate E-P flux are setting to 61 !!$ !! zero (z_cflxemp=0) to calculate the correction velocity. So 62 !!$ !! it will only balance the flux through open boundaries. 63 !!$ !! (set nn_volctl to 0 in tne namelist for this option) 64 !!$ !! 2/ The volume is constant even with E-P flux. In this case 65 !!$ !! the correction velocity must balance both the flux 66 !!$ !! through open boundaries and the ones through the free 67 !!$ !! surface. 68 !!$ !! (set nn_volctl to 1 in tne namelist for this option) 69 !!$ !!---------------------------------------------------------------------- 70 !!$ INTEGER, INTENT( in ) :: kt ! ocean time-step index 71 !!$ !! 72 !!$ INTEGER :: ji, jj, jk, jb, jgrd 73 !!$ INTEGER :: ii, ij 74 !!$ REAL(wp) :: zubtpecor, z_cflxemp, ztranst 75 !!$ !!----------------------------------------------------------------------------- 76 !!$ 77 !!$ IF( ln_vol ) THEN 78 !!$ 79 !!$ IF( kt == nit000 ) THEN 80 !!$ IF(lwp) WRITE(numout,*) 81 !!$ IF(lwp) WRITE(numout,*)'obc_vol : Correction of velocities along unstructured OBC' 82 !!$ IF(lwp) WRITE(numout,*)'~~~~~~~' 83 !!$ END IF 84 !!$ 85 !!$ ! Calculate the cumulate surface Flux z_cflxemp (m3/s) over all the domain 86 !!$ ! ----------------------------------------------------------------------- 87 !!$ z_cflxemp = SUM ( ( emp(:,:)-rnf(:,:) ) * obctmask(:,:) * e1t(:,:) * e2t(:,:) ) / rau0 88 !!$ IF( lk_mpp ) CALL mpp_sum( z_cflxemp ) ! sum over the global domain 89 !!$ 90 !!$ ! Transport through the unstructured open boundary 91 !!$ ! ------------------------------------------------ 92 !!$ zubtpecor = 0.e0 93 !!$ jgrd = 2 ! cumulate u component contribution first 94 !!$ DO jb = 1, nblenrim(jgrd) 95 !!$ DO jk = 1, jpkm1 96 !!$ ii = nbi(jb,jgrd) 97 !!$ ij = nbj(jb,jgrd) 98 !!$ zubtpecor = zubtpecor + flagu(jb) * ua(ii,ij, jk) * e2u(ii,ij) * fse3u(ii,ij,jk) 99 !!$ END DO 100 !!$ END DO 101 !!$ jgrd = 3 ! then add v component contribution 102 !!$ DO jb = 1, nblenrim(jgrd) 103 !!$ DO jk = 1, jpkm1 104 !!$ ii = nbi(jb,jgrd) 105 !!$ ij = nbj(jb,jgrd) 106 !!$ zubtpecor = zubtpecor + flagv(jb) * va(ii,ij, jk) * e1v(ii,ij) * fse3v(ii,ij,jk) 107 !!$ END DO 108 !!$ END DO 109 !!$ IF( lk_mpp ) CALL mpp_sum( zubtpecor ) ! sum over the global domain 110 !!$ 111 !!$ ! The normal velocity correction 112 !!$ ! ------------------------------ 113 !!$ IF( nn_volctl==1 ) THEN ; zubtpecor = ( zubtpecor - z_cflxemp) / obcsurftot 114 !!$ ELSE ; zubtpecor = zubtpecor / obcsurftot 115 !!$ END IF 116 !!$ 117 !!$ ! Correction of the total velocity on the unstructured boundary to respect the mass flux conservation 118 !!$ ! ------------------------------------------------------------- 119 !!$ ztranst = 0.e0 120 !!$ jgrd = 2 ! correct u component 121 !!$ DO jb = 1, nblenrim(jgrd) 122 !!$ DO jk = 1, jpkm1 123 !!$ ii = nbi(jb,jgrd) 124 !!$ ij = nbj(jb,jgrd) 125 !!$ ua(ii,ij,jk) = ua(ii,ij,jk) - flagu(jb) * zubtpecor * umask(ii,ij,jk) 126 !!$ ztranst = ztranst + flagu(jb) * ua(ii,ij,jk) * e2u(ii,ij) * fse3u(ii,ij,jk) 127 !!$ END DO 128 !!$ END DO 129 !!$ jgrd = 3 ! correct v component 130 !!$ DO jb = 1, nblenrim(jgrd) 131 !!$ DO jk = 1, jpkm1 132 !!$ ii = nbi(jb,jgrd) 133 !!$ ij = nbj(jb,jgrd) 134 !!$ va(ii,ij,jk) = va(ii,ij,jk) -flagv(jb) * zubtpecor * vmask(ii,ij,jk) 135 !!$ ztranst = ztranst + flagv(jb) * va(ii,ij,jk) * e1v(ii,ij) * fse3v(ii,ij,jk) 136 !!$ END DO 137 !!$ END DO 138 !!$ IF( lk_mpp ) CALL mpp_sum( ztranst ) ! sum over the global domain 139 !!$ 140 !!$ ! Check the cumulated transport through unstructured OBC once barotropic velocities corrected 141 !!$ ! ------------------------------------------------------ 142 !!$ IF( lwp .AND. MOD( kt, nwrite ) == 0) THEN 143 !!$ IF(lwp) WRITE(numout,*) 144 !!$ IF(lwp) WRITE(numout,*)'obc_vol : time step :', kt 145 !!$ IF(lwp) WRITE(numout,*)'~~~~~~~ ' 146 !!$ IF(lwp) WRITE(numout,*)' cumulate flux EMP =', z_cflxemp , ' (m3/s)' 147 !!$ IF(lwp) WRITE(numout,*)' total lateral surface of OBC =', obcsurftot, '(m2)' 148 !!$ IF(lwp) WRITE(numout,*)' correction velocity zubtpecor =', zubtpecor , '(m/s)' 149 !!$ IF(lwp) WRITE(numout,*)' cumulated transport ztranst =', ztranst , '(m3/s)' 150 !!$ END IF 151 !!$ ! 152 !!$ END IF ! ln_vol 153 !!$ 154 !!$ END SUBROUTINE obc_vol 155 !!$ 156 !!$#else 12 !!---------------------------------------------------------------------- 13 !! 'key_obc' AND unstructured open boundary conditions 14 !! 'key_dynspg_flt' filtered free surface 15 !!---------------------------------------------------------------------- 16 USE oce ! ocean dynamics and tracers 17 USE dom_oce ! ocean space and time domain 18 USE phycst ! physical constants 19 USE obc_oce ! ocean open boundary conditions 20 USE lib_mpp ! for mppsum 21 USE in_out_manager ! I/O manager 22 USE sbc_oce ! ocean surface boundary conditions 23 24 IMPLICIT NONE 25 PRIVATE 26 27 PUBLIC obc_vol ! routine called by dynspg_flt.h90 28 29 !! * Substitutions 30 # include "domzgr_substitute.h90" 31 !!---------------------------------------------------------------------- 32 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 33 !! $Id$ 34 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 35 !!---------------------------------------------------------------------- 36 CONTAINS 37 38 SUBROUTINE obc_vol( kt ) 39 !!---------------------------------------------------------------------- 40 !! *** ROUTINE obcvol *** 41 !! 42 !! ** Purpose : This routine is called in dynspg_flt to control 43 !! the volume of the system. A correction velocity is calculated 44 !! to correct the total transport through the unstructured OBC. 45 !! The total depth used is constant (H0) to be consistent with the 46 !! linear free surface coded in OPA 8.2 47 !! 48 !! ** Method : The correction velocity (zubtpecor here) is defined calculating 49 !! the total transport through all open boundaries (trans_obc) minus 50 !! the cumulate E-P flux (z_cflxemp) divided by the total lateral 51 !! surface (obcsurftot) of the unstructured boundary. 52 !! zubtpecor = [trans_obc - z_cflxemp ]*(1./obcsurftot) 53 !! with z_cflxemp => sum of (Evaporation minus Precipitation) 54 !! over all the domain in m3/s at each time step. 55 !! z_cflxemp < 0 when precipitation dominate 56 !! z_cflxemp > 0 when evaporation dominate 57 !! 58 !! There are 2 options (user's desiderata): 59 !! 1/ The volume changes according to E-P, this is the default 60 !! option. In this case the cumulate E-P flux are setting to 61 !! zero (z_cflxemp=0) to calculate the correction velocity. So 62 !! it will only balance the flux through open boundaries. 63 !! (set nn_volctl to 0 in tne namelist for this option) 64 !! 2/ The volume is constant even with E-P flux. In this case 65 !! the correction velocity must balance both the flux 66 !! through open boundaries and the ones through the free 67 !! surface. 68 !! (set nn_volctl to 1 in tne namelist for this option) 69 !!---------------------------------------------------------------------- 70 INTEGER, INTENT( in ) :: kt ! ocean time-step index 71 !! 72 INTEGER :: ji, jj, jk, jb, jgrd 73 INTEGER :: ib_obc, ii, ij 74 REAL(wp) :: zubtpecor, z_cflxemp, ztranst 75 TYPE(OBC_INDEX), POINTER :: idx 76 !!----------------------------------------------------------------------------- 77 78 IF( ln_vol ) THEN 79 80 IF( kt == nit000 ) THEN 81 IF(lwp) WRITE(numout,*) 82 IF(lwp) WRITE(numout,*)'obc_vol : Correction of velocities along unstructured OBC' 83 IF(lwp) WRITE(numout,*)'~~~~~~~' 84 END IF 85 86 ! Calculate the cumulate surface Flux z_cflxemp (m3/s) over all the domain 87 ! ----------------------------------------------------------------------- 88 z_cflxemp = SUM ( ( emp(:,:)-rnf(:,:) ) * obctmask(:,:) * e1t(:,:) * e2t(:,:) ) / rau0 89 IF( lk_mpp ) CALL mpp_sum( z_cflxemp ) ! sum over the global domain 90 91 ! Transport through the unstructured open boundary 92 ! ------------------------------------------------ 93 zubtpecor = 0.e0 94 DO ib_obc = 1, nb_obc 95 idx => idx_obc(ib_obc) 96 97 jgrd = 2 ! cumulate u component contribution first 98 DO jb = 1, idx%nblenrim(jgrd) 99 DO jk = 1, jpkm1 100 ii = idx%nbi(jb,jgrd) 101 ij = idx%nbj(jb,jgrd) 102 zubtpecor = zubtpecor + idx%flagu(jb) * ua(ii,ij, jk) * e2u(ii,ij) * fse3u(ii,ij,jk) 103 END DO 104 END DO 105 jgrd = 3 ! then add v component contribution 106 DO jb = 1, idx%nblenrim(jgrd) 107 DO jk = 1, jpkm1 108 ii = idx%nbi(jb,jgrd) 109 ij = idx%nbj(jb,jgrd) 110 zubtpecor = zubtpecor + idx%flagv(jb) * va(ii,ij, jk) * e1v(ii,ij) * fse3v(ii,ij,jk) 111 END DO 112 END DO 113 114 END DO 115 IF( lk_mpp ) CALL mpp_sum( zubtpecor ) ! sum over the global domain 116 117 ! The normal velocity correction 118 ! ------------------------------ 119 IF( nn_volctl==1 ) THEN ; zubtpecor = ( zubtpecor - z_cflxemp) / obcsurftot 120 ELSE ; zubtpecor = zubtpecor / obcsurftot 121 END IF 122 123 ! Correction of the total velocity on the unstructured boundary to respect the mass flux conservation 124 ! ------------------------------------------------------------- 125 ztranst = 0.e0 126 DO ib_obc = 1, nb_obc 127 idx => idx_obc(ib_obc) 128 129 jgrd = 2 ! correct u component 130 DO jb = 1, idx%nblenrim(jgrd) 131 DO jk = 1, jpkm1 132 ii = idx%nbi(jb,jgrd) 133 ij = idx%nbj(jb,jgrd) 134 ua(ii,ij,jk) = ua(ii,ij,jk) - idx%flagu(jb) * zubtpecor * umask(ii,ij,jk) 135 ztranst = ztranst + idx%flagu(jb) * ua(ii,ij,jk) * e2u(ii,ij) * fse3u(ii,ij,jk) 136 END DO 137 END DO 138 jgrd = 3 ! correct v component 139 DO jb = 1, idx%nblenrim(jgrd) 140 DO jk = 1, jpkm1 141 ii = idx%nbi(jb,jgrd) 142 ij = idx%nbj(jb,jgrd) 143 va(ii,ij,jk) = va(ii,ij,jk) -idx%flagv(jb) * zubtpecor * vmask(ii,ij,jk) 144 ztranst = ztranst + idx%flagv(jb) * va(ii,ij,jk) * e1v(ii,ij) * fse3v(ii,ij,jk) 145 END DO 146 END DO 147 148 END DO 149 IF( lk_mpp ) CALL mpp_sum( ztranst ) ! sum over the global domain 150 151 ! Check the cumulated transport through unstructured OBC once barotropic velocities corrected 152 ! ------------------------------------------------------ 153 IF( lwp .AND. MOD( kt, nwrite ) == 0) THEN 154 IF(lwp) WRITE(numout,*) 155 IF(lwp) WRITE(numout,*)'obc_vol : time step :', kt 156 IF(lwp) WRITE(numout,*)'~~~~~~~ ' 157 IF(lwp) WRITE(numout,*)' cumulate flux EMP =', z_cflxemp , ' (m3/s)' 158 IF(lwp) WRITE(numout,*)' total lateral surface of OBC =', obcsurftot, '(m2)' 159 IF(lwp) WRITE(numout,*)' correction velocity zubtpecor =', zubtpecor , '(m/s)' 160 IF(lwp) WRITE(numout,*)' cumulated transport ztranst =', ztranst , '(m3/s)' 161 END IF 162 ! 163 END IF ! ln_vol 164 165 END SUBROUTINE obc_vol 166 167 #else 157 168 !!---------------------------------------------------------------------- 158 169 !! Dummy module NO Unstruct Open Boundary Conditions -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90
r2818 r2865 104 104 CONTAINS 105 105 106 SUBROUTINE fld_read( kt, kn_fsbc, sd, map, jit, time shift )106 SUBROUTINE fld_read( kt, kn_fsbc, sd, map, jit, time_offset ) 107 107 !!--------------------------------------------------------------------- 108 108 !! *** ROUTINE fld_read *** … … 121 121 TYPE(MAP_POINTER),INTENT(in), OPTIONAL, DIMENSION(:) :: map ! global-to-local mapping index 122 122 INTEGER , INTENT(in ), OPTIONAL :: jit ! subcycle timestep for timesplitting option 123 INTEGER , INTENT(in ), OPTIONAL :: timeshift ! provide fields at time other than "now" 123 INTEGER , INTENT(in ), OPTIONAL :: time_offset ! provide fields at time other than "now" 124 ! time_offset = -1 => fields at "before" time level 125 ! time_offset = +1 => fields at "after" time levels 126 ! etc. 124 127 !! 125 128 INTEGER :: imf ! size of the structure sd … … 128 131 INTEGER :: isecend ! number of second since Jan. 1st 00h of nit000 year at nitend 129 132 INTEGER :: isecsbc ! number of seconds between Jan. 1st 00h of nit000 year and the middle of sbc time step 133 INTEGER :: time_add ! local time_offset variable 130 134 LOGICAL :: llnxtyr ! open next year file? 131 135 LOGICAL :: llnxtmth ! open next month file? … … 143 147 ENDIF 144 148 149 time_add = 0 150 IF( PRESENT(time_offset) ) THEN 151 time_add = time_offset 152 ENDIF 153 145 154 ! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar 146 155 IF( present(jit) ) THEN 147 156 ! ignore kn_fsbc in this case 148 isecsbc = nsec_year + nsec1jan000 + jit*rdt/REAL(nn_baro,wp) 149 ELSE IF( present(timeshift) ) THEN 150 isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) + timeshift * rdttra(1) ! middle of sbc time step 157 isecsbc = nsec_year + nsec1jan000 + (jit+time_add)*rdt/REAL(nn_baro,wp) 151 158 ELSE 152 isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) ! middle of sbc time step159 isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) + time_add * rdttra(1) ! middle of sbc time step 153 160 ENDIF 154 161 imf = SIZE( sd ) … … 263 270 clfmt = "('fld_read: var ', a, ' kt = ', i8, ' (', f7.2,' days), Y/M/D = ', i4.4,'/', i2.2,'/', i2.2," // & 264 271 & "', records b/a: ', i4.4, '/', i4.4, ' (days ', f7.2,'/', f7.2, ')')" 265 WRITE(numout, clfmt) TRIM( sd(jf)%clvar ), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 272 WRITE(numout, clfmt) TRIM( sd(jf)%clvar ), kt, REAL(isecsbc,wp)/rday, nyear, nmonth, nday, & 266 273 & sd(jf)%nrec_b(1), sd(jf)%nrec_a(1), REAL(sd(jf)%nrec_b(2),wp)/rday, REAL(sd(jf)%nrec_a(2),wp)/rday 274 WRITE(numout, *) 'time_add is : ',time_add 267 275 ENDIF 268 276 ! temporal interpolation weights -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/oce.F90
r2715 r2865 35 35 !! free surface ! before ! now ! after ! 36 36 !! ------------ ! fields ! fields ! trends ! 37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshb , sshn , ssha !: sea surface height at t-point [m]38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshu_b , sshu_n , sshu_a !: sea surface height at u-point [m]39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshv_b , sshv_n , sshv_a !: sea surface height at u-point [m]40 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshf_n !: sea surface height at f-point [m]37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: sshb , sshn , ssha !: sea surface height at t-point [m] 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshu_b , sshu_n , sshu_a !: sea surface height at u-point [m] 39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshv_b , sshv_n , sshv_a !: sea surface height at u-point [m] 40 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshf_n !: sea surface height at f-point [m] 41 41 ! 42 42 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: spgu, spgv !: horizontal surface pressure gradient -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/step.F90
r2797 r2865 97 97 IF( lk_dtasal ) CALL dta_sal( kstp ) ! update 3D salinity data 98 98 CALL sbc ( kstp ) ! Sea Boundary Condition (including sea-ice) 99 IF( lk_obc ) CALL obc_dta( kstp )! update dynamic and tracer data at open boundaries99 IF( lk_obc ) CALL obc_dta( kstp, time_offset=+1 ) ! update dynamic and tracer data at open boundaries 100 100 101 101 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
Note: See TracChangeset
for help on using the changeset viewer.
