Changeset 2800
- Timestamp:
- 2011-07-13T15:31:05+02:00 (13 years ago)
- Location:
- branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 1 added
- 9 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
r2715 r2800 150 150 LOGICAL, PUBLIC, PARAMETER :: lk_vvl = .FALSE. !: fixed grid flag 151 151 #endif 152 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hur , hvr !: inverse of u and v-points ocean depth (1/m)153 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu , hv !: depth at u- and v-points (meters)154 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_0 , hv_0 !: refernce depth at u- and v-points (meters)152 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: hur , hvr !: inverse of u and v-points ocean depth (1/m) 153 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu , hv !: depth at u- and v-points (meters) 154 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_0 , hv_0 !: refernce depth at u- and v-points (meters) 155 155 156 156 INTEGER, PUBLIC :: nla10 !: deepest W level Above ~10m (nlb10 - 1) -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r2797 r2800 30 30 USE domvvl ! variable volume 31 31 USE obc_oce ! ocean open boundary conditions 32 USE obcdyn3d ! open boundary condition for baroclinic velocities 33 USE obcdyn2d ! open boundary condition for barotropic variables 32 USE obcdyn ! open boundary condition for baroclinic velocities 34 33 USE obcvol ! ocean open boundary condition (obc_vol routines) 35 34 USE in_out_manager ! I/O manager … … 155 154 IF( .NOT. lk_dynspg_flt ) THEN 156 155 157 CALL obc_dyn3d( kt ) 158 ! 159 !!!! ENDA'S FIX: NEED TO THINK ABOUT THIS !!!! 160 CALL obc_dta( kt+1, jit=0 ) 161 CALL obc_dyn2d( kt, sshn_b ) 162 ! 163 !!gm ERROR - potential BUG: sshn should not be modified at this stage !! ssh_nxt not alrady called 164 CALL lbc_lnk( sshn, 'T', 1. ) ! Boundary conditions on sshn 165 ! 166 IF( ln_vol_cst ) CALL obc_vol( kt ) 167 ! 168 IF(ln_ctl) CALL prt_ctl( tab2d_1=sshn, clinfo1=' ssh : ', mask1=tmask ) 156 CALL obc_dyn( kt ) 157 158 !!$!!gm ERROR - potential BUG: sshn should not be modified at this stage !! ssh_nxt not alrady called 159 !!$ CALL lbc_lnk( sshn, 'T', 1. ) ! Boundary conditions on sshn 160 !!$ ! 161 !!$ IF( ln_vol_cst ) CALL obc_vol( kt ) 162 !!$ ! 163 !!$ IF(ln_ctl) CALL prt_ctl( tab2d_1=sshn, clinfo1=' ssh : ', mask1=tmask ) 164 169 165 ENDIF 170 166 ! -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_flt.F90
r2797 r2800 33 33 USE solpcg ! preconditionned conjugate gradient solver 34 34 USE solsor ! Successive Over-relaxation solver 35 USE obcdyn 3d ! ocean open boundary condition (obc_dyn3d routines)35 USE obcdyn ! ocean open boundary condition on dynamics 36 36 USE obcvol ! ocean open boundary condition (obc_vol routines) 37 37 USE cla ! cross land advection … … 180 180 181 181 #if defined key_obc 182 CALL obc_dyn 3d( kt )! Update velocities on each open boundary183 CALL obc_vol( kt ) ! Correction of the barotropic compon ant velocity to control the volume of the system182 CALL obc_dyn( kt ) ! Update velocities on each open boundary 183 CALL obc_vol( kt ) ! Correction of the barotropic component velocity to control the volume of the system 184 184 #endif 185 185 #if defined key_agrif -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_oce.F90
r2715 r2800 34 34 35 35 ! !!! Time splitting scheme (key_dynspg_ts) 36 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshn_e, ssha_e ! sea surface heigth (now, after, average)37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ua_e , va_e ! barotropic velocities (after)38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_e , hv_e ! now ocean depth ( = Ho+sshn_e )39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hur_e , hvr_e ! inverse of hu_e and hv_e40 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sshn_b! before field without time-filter36 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: sshn_e, ssha_e ! sea surface heigth (now, after, average) 37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: ua_e , va_e ! barotropic velocities (after) 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_e , hv_e ! now ocean depth ( = Ho+sshn_e ) 39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: hur_e , hvr_e ! inverse of hu_e and hv_e 40 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: sshn_b ! before field without time-filter 41 41 42 42 !!---------------------------------------------------------------------- -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r2797 r2800 466 466 ! ! ----------------------- 467 467 ! OBC open boundaries 468 IF( lk_obc .OR. ln_tides ) CALL obc_dyn2d( kt, sshn_e ) 468 #if defined key_obc 469 pssh => sshn_e 470 phur => hur_e 471 phvr => hvr_e 472 pu2d => ua_e 473 pv2d => va_e 474 475 IF( lk_obc ) CALL obc_dyn2d( kt ) 476 #endif 477 469 478 ! 470 479 CALL lbc_lnk( ua_e , 'U', -1. ) ! local domain boundaries -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obc_oce.F90
r2797 r2800 88 88 REAL(wp) :: obcsurftot !: Lateral surface of unstructured open boundary 89 89 90 REAL(wp), POINTER, DIMENSION(:,:) :: pssh !: 91 REAL(wp), POINTER, DIMENSION(:,:) :: phur !: 92 REAL(wp), POINTER, DIMENSION(:,:) :: phvr !: Pointers for barotropic fields 93 REAL(wp), POINTER, DIMENSION(:,:) :: pu2d !: 94 REAL(wp), POINTER, DIMENSION(:,:) :: pv2d !: 95 90 96 !!---------------------------------------------------------------------- 91 97 !! open boundary data variables -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcdta.F90
r2797 r2800 22 22 USE dom_oce ! ocean space and time domain 23 23 USE phycst ! physical constants 24 USE obc_oce ! ocean open boundary conditions 24 USE obc_oce ! ocean open boundary conditions 25 25 USE obctides ! tidal forcing at boundaries 26 26 USE fldread ! read input fields … … 40 40 INTEGER :: nb_obc_fld_sum ! Total number of fields to update for all boundary sets. 41 41 42 LOGICAL, DIMENSION(jp_obc) :: ln_full_vel_array ! =T => full velocities in 3D boundary conditions 43 ! =F => baroclinic velocities in 3D boundary conditions 44 42 45 TYPE(FLD), PUBLIC, ALLOCATABLE, DIMENSION(:), TARGET :: bf ! structure of input fields (file informations, fields read) 43 46 44 47 TYPE(MAP_POINTER), ALLOCATABLE, DIMENSION(:) :: nbmap_ptr ! array of pointers to nbmap 45 48 49 # include "domzgr_substitute.h90" 46 50 !!---------------------------------------------------------------------- 47 51 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 60 64 !! 61 65 !!---------------------------------------------------------------------- 66 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 67 USE wrk_nemo, ONLY: wrk_2d_9, wrk_2d_10 ! 2D workspace 68 !! 62 69 INTEGER, INTENT( in ) :: kt ! ocean time-step index 63 70 INTEGER, INTENT( in ), OPTIONAL :: jit ! subcycle time-step index (for timesplitting option) 64 71 !! 65 INTEGER :: ib_obc, jfld, jstart, jend ! local indices 66 INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts 72 INTEGER :: ib_obc, jfld, jstart, jend, ib, ii, ij, ik, igrd ! local indices 73 INTEGER, DIMENSION(jpbgrd) :: ilen1 74 INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts 67 75 !! 68 76 !!--------------------------------------------------------------------------- 77 78 IF(wrk_in_use(2, 9,10) ) THEN 79 CALL ctl_stop('obc_dta: ERROR: requested workspace arrays are unavailable.') ; RETURN 80 END IF 69 81 70 82 ! for nn_dtactl = 0, initialise data arrays once for all 71 83 ! from initial conditions 72 84 !------------------------------------------------------- 73 IF( kt .eq. 1 .and. .not. PRESENT(jit) ) THEN 74 85 IF( kt .eq. nit000 .and. .not. PRESENT(jit) ) THEN 86 87 ! Calculate depth-mean currents 88 !----------------------------- 89 pu2d => wrk_2d_9 90 pu2d => wrk_2d_10 91 92 pu2d(:,:) = 0.e0 93 pv2d(:,:) = 0.e0 94 95 DO ik = 1, jpkm1 !! Vertically integrated momentum trends 96 pu2d(:,:) = pu2d(:,:) + fse3u(:,:,ik) * umask(:,:,ik) * un(:,:,ik) 97 pv2d(:,:) = pv2d(:,:) + fse3v(:,:,ik) * vmask(:,:,ik) * vn(:,:,ik) 98 END DO 99 pu2d(:,:) = pu2d(:,:) * hur(:,:) 100 pv2d(:,:) = pv2d(:,:) * hvr(:,:) 101 75 102 DO ib_obc = 1, nb_obc 76 103 IF( nn_dtactl(ib_obc) .eq. 0 ) THEN 77 104 78 !!! TO BE DONE !!! 105 nblen => idx_obc(ib_obc)%nblen 106 nblenrim => idx_obc(ib_obc)%nblenrim 107 108 IF( nn_dyn2d(ib_obc) .gt. 0 ) THEN 109 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 110 ilen1(:) = nblen(:) 111 ELSE 112 ilen1(:) = nblenrim(:) 113 ENDIF 114 igrd = 1 115 DO ib = 1, ilen1(igrd) 116 ii = idx_obc(ib_obc)%nbi(ib,igrd) 117 ij = idx_obc(ib_obc)%nbj(ib,igrd) 118 dta_obc(ib_obc)%ssh(ib) = sshn(ii,ij) * tmask(ii,ij,1) 119 END DO 120 igrd = 2 121 DO ib = 1, ilen1(igrd) 122 ii = idx_obc(ib_obc)%nbi(ib,igrd) 123 ij = idx_obc(ib_obc)%nbj(ib,igrd) 124 dta_obc(ib_obc)%u2d(ib) = pu2d(ii,ij) * umask(ii,ij,1) 125 END DO 126 igrd = 3 127 DO ib = 1, ilen1(igrd) 128 ii = idx_obc(ib_obc)%nbi(ib,igrd) 129 ij = idx_obc(ib_obc)%nbj(ib,igrd) 130 dta_obc(ib_obc)%v2d(ib) = pv2d(ii,ij) * vmask(ii,ij,1) 131 END DO 132 ENDIF 133 134 IF( nn_dyn3d(ib_obc) .gt. 0 ) THEN 135 IF( nn_dyn3d(ib_obc) .eq. jp_frs ) THEN 136 ilen1(:) = nblen(:) 137 ELSE 138 ilen1(:) = nblenrim(:) 139 ENDIF 140 igrd = 2 141 DO ib = 1, ilen1(igrd) 142 DO ik = 1, jpkm1 143 ii = idx_obc(ib_obc)%nbi(ib,igrd) 144 ij = idx_obc(ib_obc)%nbj(ib,igrd) 145 dta_obc(ib_obc)%u3d(ib,ik) = ( un(ii,ij,ik) - pu2d(ii,ij) ) * umask(ii,ij,ik) 146 END DO 147 END DO 148 igrd = 3 149 DO ib = 1, ilen1(igrd) 150 DO ik = 1, jpkm1 151 ii = idx_obc(ib_obc)%nbi(ib,igrd) 152 ij = idx_obc(ib_obc)%nbj(ib,igrd) 153 dta_obc(ib_obc)%v3d(ib,ik) = ( vn(ii,ij,ik) - pv2d(ii,ij) ) * vmask(ii,ij,ik) 154 END DO 155 END DO 156 ENDIF 157 158 IF( nn_tra(ib_obc) .gt. 0 ) THEN 159 IF( nn_tra(ib_obc) .eq. jp_frs ) THEN 160 ilen1(:) = nblen(:) 161 ELSE 162 ilen1(:) = nblenrim(:) 163 ENDIF 164 igrd = 1 ! Everything is at T-points here 165 DO ib = 1, ilen1(igrd) 166 DO ik = 1, jpkm1 167 ii = idx_obc(ib_obc)%nbi(ib,igrd) 168 ij = idx_obc(ib_obc)%nbj(ib,igrd) 169 dta_obc(ib_obc)%tem(ib,ik) = tn(ii,ij,ik) * tmask(ii,ij,ik) 170 dta_obc(ib_obc)%sal(ib,ik) = sn(ii,ij,ik) * tmask(ii,ij,ik) 171 END DO 172 END DO 173 ENDIF 174 175 #if defined key_lim2 176 IF( nn_ice_lim2(ib_obc) .gt. 0 ) THEN 177 IF( nn_ice_lim2(ib_obc) .eq. jp_frs ) THEN 178 ilen1(:) = nblen(:) 179 ELSE 180 ilen1(:) = nblenrim(:) 181 ENDIF 182 igrd = 1 ! Everything is at T-points here 183 DO ib = 1, ilen1(igrd) 184 ii = idx_obc(ib_obc)%nbi(ib,igrd) 185 ij = idx_obc(ib_obc)%nbj(ib,igrd) 186 dta_obc(ib_obc)%frld(ib) = frld(ii,ij) * tmask(ii,ij,1) 187 dta_obc(ib_obc)%hicif(ib) = hicif(ii,ij) * tmask(ii,ij,1) 188 dta_obc(ib_obc)%hsnif(ib) = hsnif(ii,ij) * tmask(ii,ij,1) 189 END DO 190 ENDIF 191 #endif 79 192 80 193 ENDIF … … 103 216 jstart = jend+1 104 217 218 ! If full velocities in boundary data then split into barotropic and baroclinic data 219 ! (Note that we have already made sure that you can't use ln_full_vel = .true. at the same 220 ! time as the dynspg_ts option). 221 222 IF( ln_full_vel_array(ib_obc) ) THEN 223 224 igrd = 2 ! zonal velocity 225 dta_obc(ib_obc)%u2d(:) = 0.0 226 DO ib = 1, idx_obc(ib_obc)%nblen(igrd) 227 ii = idx_obc(ib_obc)%nbi(ib,igrd) 228 ij = idx_obc(ib_obc)%nbj(ib,igrd) 229 DO ik = 1, jpkm1 230 dta_obc(ib_obc)%u2d(ib) = dta_obc(ib_obc)%u2d(ib) & 231 & + fse3u(ii,ij,ik) * umask(ii,ij,ik) * dta_obc(ib_obc)%u3d(ib,ik) 232 END DO 233 dta_obc(ib_obc)%u2d(ib) = dta_obc(ib_obc)%u2d(ib) * hur(ii,ij) 234 DO ik = 1, jpkm1 235 dta_obc(ib_obc)%u3d(ib,ik) = dta_obc(ib_obc)%u3d(ib,ik) - dta_obc(ib_obc)%u2d(ib) 236 END DO 237 END DO 238 239 igrd = 3 ! meridional velocity 240 dta_obc(ib_obc)%v2d(:) = 0.0 241 DO ib = 1, idx_obc(ib_obc)%nblen(igrd) 242 ii = idx_obc(ib_obc)%nbi(ib,igrd) 243 ij = idx_obc(ib_obc)%nbj(ib,igrd) 244 DO ik = 1, jpkm1 245 dta_obc(ib_obc)%v2d(ib) = dta_obc(ib_obc)%v2d(ib) & 246 & + fse3v(ii,ij,ik) * vmask(ii,ij,ik) * dta_obc(ib_obc)%v3d(ib,ik) 247 END DO 248 dta_obc(ib_obc)%v2d(ib) = dta_obc(ib_obc)%v2d(ib) * hvr(ii,ij) 249 DO ik = 1, jpkm1 250 dta_obc(ib_obc)%v3d(ib,ik) = dta_obc(ib_obc)%v3d(ib,ik) - dta_obc(ib_obc)%v2d(ib) 251 END DO 252 END DO 253 254 ENDIF 255 105 256 END IF ! nn_dtactl(ib_obc) = 1 106 257 END DO ! ib_obc 258 259 IF(wrk_not_released(2, 9,10) ) CALL ctl_stop('obc_dta: ERROR: failed to release workspace arrays.') 107 260 108 261 END SUBROUTINE obc_dta … … 119 272 !! 120 273 !!---------------------------------------------------------------------- 274 USE dynspg_oce, ONLY: lk_dynspg_ts 275 !! 121 276 INTEGER :: ib_obc, jfld, jstart, jend, ierror ! local indices 122 277 !! 123 278 CHARACTER(len=100) :: cn_dir ! Root directory for location of data files 124 279 CHARACTER(len=100), DIMENSION(nb_obc) :: cn_dir_array ! Root directory for location of data files 280 LOGICAL :: ln_full_vel ! =T => full velocities in 3D boundary data 281 ! =F => baroclinic velocities in 3D boundary data 125 282 INTEGER :: ilen_global ! Max length required for global obc dta arrays 283 INTEGER, DIMENSION(jpbgrd) :: ilen0 ! size of local arrays 126 284 INTEGER, ALLOCATABLE, DIMENSION(:) :: ilen1, ilen3 ! size of 1st and 3rd dimensions of local arrays 127 285 INTEGER, ALLOCATABLE, DIMENSION(:) :: iobc ! obc set for a particular jfld … … 138 296 NAMELIST/namobc_dta/ bn_frld, bn_hicif, bn_hsnif 139 297 #endif 298 NAMELIST/namobc_dta/ ln_full_vel 140 299 !!--------------------------------------------------------------------------- 141 300 142 ! Work out how many fields there are to read in and allocate arrays143 ! ----------------------------------------------------------------- 301 ! Work out upper bound of how many fields there are to read in and allocate arrays 302 ! --------------------------------------------------------------------------- 144 303 ALLOCATE( nb_obc_fld(nb_obc) ) 145 304 nb_obc_fld(:) = 0 … … 189 348 ! set file information 190 349 cn_dir = './' ! directory in which the model is executed 350 ln_full_vel = .false. 191 351 ! ... default values (NB: frequency positive => hours, negative => months) 192 352 ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation ! … … 209 369 210 370 cn_dir_array(ib_obc) = cn_dir 371 ln_full_vel_array(ib_obc) = ln_full_vel 372 373 IF( ln_full_vel_array(ib_obc) .and. lk_dynspg_ts ) THEN 374 CALL ctl_stop( 'obc_dta_init: ERROR, cannot specify full velocities in boundary data',& 375 & 'with dynspg_ts option' ) ; RETURN 376 ENDIF 211 377 212 378 nblen => idx_obc(ib_obc)%nblen … … 217 383 IF( nn_dyn2d(ib_obc) .gt. 0 ) THEN 218 384 219 jfld = jfld + 1 220 blf_i(jfld) = bn_ssh 221 iobc(jfld) = ib_obc 222 igrid(jfld) = 1 223 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 224 ilen1(jfld) = nblen(igrid(jfld)) 225 ELSE 226 ilen1(jfld) = nblenrim(igrid(jfld)) 227 ENDIF 228 ilen3(jfld) = 1 229 230 jfld = jfld + 1 231 blf_i(jfld) = bn_u2d 232 iobc(jfld) = ib_obc 233 igrid(jfld) = 2 234 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 235 ilen1(jfld) = nblen(igrid(jfld)) 236 ELSE 237 ilen1(jfld) = nblenrim(igrid(jfld)) 238 ENDIF 239 ilen3(jfld) = 1 240 241 jfld = jfld + 1 242 blf_i(jfld) = bn_v2d 243 iobc(jfld) = ib_obc 244 igrid(jfld) = 3 245 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 246 ilen1(jfld) = nblen(igrid(jfld)) 247 ELSE 248 ilen1(jfld) = nblenrim(igrid(jfld)) 249 ENDIF 250 ilen3(jfld) = 1 385 IF( nn_dyn2d(ib_obc) .ne. jp_frs ) THEN 386 jfld = jfld + 1 387 blf_i(jfld) = bn_ssh 388 iobc(jfld) = ib_obc 389 igrid(jfld) = 1 390 ilen1(jfld) = nblenrim(igrid(jfld)) 391 ilen3(jfld) = 1 392 ENDIF 393 394 IF( .not. ln_full_vel_array(ib_obc) ) THEN 395 396 jfld = jfld + 1 397 blf_i(jfld) = bn_u2d 398 iobc(jfld) = ib_obc 399 igrid(jfld) = 2 400 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 401 ilen1(jfld) = nblen(igrid(jfld)) 402 ELSE 403 ilen1(jfld) = nblenrim(igrid(jfld)) 404 ENDIF 405 ilen3(jfld) = 1 406 407 jfld = jfld + 1 408 blf_i(jfld) = bn_v2d 409 iobc(jfld) = ib_obc 410 igrid(jfld) = 3 411 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 412 ilen1(jfld) = nblen(igrid(jfld)) 413 ELSE 414 ilen1(jfld) = nblenrim(igrid(jfld)) 415 ENDIF 416 ilen3(jfld) = 1 417 418 ENDIF 251 419 252 420 ENDIF 253 421 254 422 ! baroclinic velocities 255 IF( nn_dyn3d(ib_obc) .gt. 0 ) THEN 423 IF( nn_dyn3d(ib_obc) .gt. 0 .or. & 424 ( ln_full_vel_array(ib_obc) .and. nn_dyn2d(ib_obc) .gt. 0 ) ) THEN 256 425 257 426 jfld = jfld + 1 … … 345 514 ENDIF 346 515 #endif 516 ! Recalculate field counts 517 !------------------------- 518 nb_obc_fld_sum = 0 519 IF( ib_obc .eq. 1 ) THEN 520 nb_obc_fld(ib_obc) = jfld 521 nb_obc_fld_sum = jfld 522 ELSE 523 nb_obc_fld(ib_obc) = jfld - nb_obc_fld_sum 524 nb_obc_fld_sum = nb_obc_fld_sum + nb_obc_fld(ib_obc) 525 ENDIF 526 347 527 ENDIF ! nn_dtactl .eq. 1 348 528 ENDDO ! ib_obc 349 529 350 IF( jfld .ne. nb_obc_fld_sum ) THEN351 CALL ctl_stop( 'obc_dta: error in initialisation: jpfld .ne. nb_obc_fld_sum' ) ; RETURN352 ENDIF353 530 354 531 DO jfld = 1, nb_obc_fld_sum … … 385 562 IF (nn_dyn2d(ib_obc) .gt. 0) THEN 386 563 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 387 ilen1(1) = nblen(1) 388 ilen1(2) = nblen(2) 389 ilen1(3) = nblen(3) 390 ELSE 391 ilen1(1) = nblenrim(1) 392 ilen1(2) = nblenrim(2) 393 ilen1(3) = nblenrim(3) 394 ENDIF 395 ALLOCATE( dta_obc(ib_obc)%ssh(ilen1(1)) ) 396 ALLOCATE( dta_obc(ib_obc)%u2d(ilen1(2)) ) 397 ALLOCATE( dta_obc(ib_obc)%v2d(ilen1(3)) ) 564 ilen0(1:3) = nblen(1:3) 565 ELSE 566 ilen0(1:3) = nblenrim(1:3) 567 ENDIF 568 ALLOCATE( dta_obc(ib_obc)%ssh(ilen0(1)) ) 569 ALLOCATE( dta_obc(ib_obc)%u2d(ilen0(2)) ) 570 ALLOCATE( dta_obc(ib_obc)%v2d(ilen0(3)) ) 398 571 ENDIF 399 572 IF (nn_dyn3d(ib_obc) .gt. 0) THEN 400 573 IF( nn_dyn3d(ib_obc) .eq. jp_frs ) THEN 401 ilen1(2) = nblen(2) 402 ilen1(3) = nblen(3) 403 ELSE 404 ilen1(2) = nblenrim(2) 405 ilen1(3) = nblenrim(3) 406 ENDIF 407 ALLOCATE( dta_obc(ib_obc)%u3d(ilen1(2),jpk) ) 408 ALLOCATE( dta_obc(ib_obc)%v3d(ilen1(3),jpk) ) 574 ilen0(1:3) = nblen(1:3) 575 ELSE 576 ilen0(1:3) = nblenrim(1:3) 577 ENDIF 578 ALLOCATE( dta_obc(ib_obc)%u3d(ilen0(2),jpk) ) 579 ALLOCATE( dta_obc(ib_obc)%v3d(ilen0(3),jpk) ) 409 580 ENDIF 410 581 IF (nn_tra(ib_obc) .gt. 0) THEN 411 582 IF( nn_tra(ib_obc) .eq. jp_frs ) THEN 412 ilen 1(1) = nblen(1)413 ELSE 414 ilen 1(1) = nblenrim(1)415 ENDIF 416 ALLOCATE( dta_obc(ib_obc)%tem(ilen 1(1),jpk) )417 ALLOCATE( dta_obc(ib_obc)%sal(ilen 1(1),jpk) )583 ilen0(1:3) = nblen(1:3) 584 ELSE 585 ilen0(1:3) = nblenrim(1:3) 586 ENDIF 587 ALLOCATE( dta_obc(ib_obc)%tem(ilen0(1),jpk) ) 588 ALLOCATE( dta_obc(ib_obc)%sal(ilen0(1),jpk) ) 418 589 ENDIF 419 590 #if defined key_lim2 420 591 IF (nn_ice_lim2(ib_obc) .gt. 0) THEN 421 592 IF( nn_ice_lim2(ib_obc) .eq. jp_frs ) THEN 422 ilen 1(1) = nblen(igrid(jfld))423 ELSE 424 ilen 1(1) = nblenrim(igrid(jfld))425 ENDIF 426 ALLOCATE( dta_obc(ib_obc)%ssh(ilen 1(1)) )427 ALLOCATE( dta_obc(ib_obc)%u2d(ilen 1(1)) )428 ALLOCATE( dta_obc(ib_obc)%v2d(ilen 1(1)) )593 ilen0(1:3) = nblen(1:3) 594 ELSE 595 ilen0(1:3) = nblenrim(1:3) 596 ENDIF 597 ALLOCATE( dta_obc(ib_obc)%ssh(ilen0(1)) ) 598 ALLOCATE( dta_obc(ib_obc)%u2d(ilen0(1)) ) 599 ALLOCATE( dta_obc(ib_obc)%v2d(ilen0(1)) ) 429 600 ENDIF 430 601 #endif … … 436 607 !----------------------------------------------------------- 437 608 IF (nn_dyn2d(ib_obc) .gt. 0) THEN 438 jfld = jfld + 1 439 dta_obc(ib_obc)%ssh => bf(jfld)%fnow(:,1,1) 440 jfld = jfld + 1 441 dta_obc(ib_obc)%u2d => bf(jfld)%fnow(:,1,1) 442 jfld = jfld + 1 443 dta_obc(ib_obc)%v2d => bf(jfld)%fnow(:,1,1) 444 ENDIF 445 IF (nn_dyn3d(ib_obc) .gt. 0) THEN 609 IF( nn_dyn2d(ib_obc) .ne. jp_frs ) THEN 610 jfld = jfld + 1 611 dta_obc(ib_obc)%ssh => bf(jfld)%fnow(:,1,1) 612 ENDIF 613 IF( ln_full_vel_array(ib_obc) ) THEN 614 ! In this case we need space but we aren't reading it 615 ! directly from the external file. 616 IF( nn_dyn2d(ib_obc) .eq. jp_frs ) THEN 617 ilen0(2) = nblen(2) 618 ilen0(3) = nblen(3) 619 ELSE 620 ilen0(2) = nblenrim(2) 621 ilen0(3) = nblenrim(3) 622 ENDIF 623 ALLOCATE( dta_obc(ib_obc)%u2d(ilen0(2)) ) 624 ALLOCATE( dta_obc(ib_obc)%v2d(ilen0(3)) ) 625 ELSE 626 jfld = jfld + 1 627 dta_obc(ib_obc)%u2d => bf(jfld)%fnow(:,1,1) 628 jfld = jfld + 1 629 dta_obc(ib_obc)%v2d => bf(jfld)%fnow(:,1,1) 630 ENDIF 631 ENDIF 632 IF (nn_dyn3d(ib_obc) .gt. 0 .or. & 633 & ln_full_vel_array(ib_obc) .and. nn_dyn2d(ib_obc) .gt. 0 ) THEN 446 634 jfld = jfld + 1 447 635 dta_obc(ib_obc)%u3d => bf(jfld)%fnow(:,1,:) -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/OBC/obcdyn2d.F90
r2797 r2800 39 39 CONTAINS 40 40 41 SUBROUTINE obc_dyn2d( kt , pssh)41 SUBROUTINE obc_dyn2d( kt ) 42 42 !!---------------------------------------------------------------------- 43 43 !! *** SUBROUTINE obc_dyn2d *** … … 47 47 !!---------------------------------------------------------------------- 48 48 INTEGER, INTENT(in) :: kt ! Main time step counter 49 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pssh50 49 !! 51 50 INTEGER :: ib_obc ! Loop counter … … 56 55 CASE(jp_none) 57 56 CYCLE 57 CASE(jp_frs) 58 CALL obc_dyn2d_frs( idx_obc(ib_obc), dta_obc(ib_obc), kt ) 58 59 CASE(jp_flather) 59 CALL obc_dyn2d_fla( idx_obc(ib_obc), dta_obc(ib_obc) , pssh)60 CALL obc_dyn2d_fla( idx_obc(ib_obc), dta_obc(ib_obc) ) 60 61 CASE DEFAULT 61 62 CALL ctl_stop( 'obc_dyn3d : unrecognised option for open boundaries for barotropic variables' ) … … 65 66 END SUBROUTINE obc_dyn2d 66 67 67 SUBROUTINE obc_dyn2d_fla( idx, dta, pssh ) 68 !!---------------------------------------------------------------------- 69 !! *** SUBROUTINE obc_dyn_fla *** 68 SUBROUTINE obc_dyn2d_frs( idx, dta, kt ) 69 !!---------------------------------------------------------------------- 70 !! *** SUBROUTINE obc_dyn2d_frs *** 71 !! 72 !! ** Purpose : - Apply the Flow Relaxation Scheme for barotropic velocities 73 !! at open boundaries. 74 !! 75 !! References :- Engedahl H., 1995: Use of the flow relaxation scheme in 76 !! a three-dimensional baroclinic ocean model with realistic 77 !! topography. Tellus, 365-382. 78 !!---------------------------------------------------------------------- 79 INTEGER, INTENT(in) :: kt 80 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 81 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data 82 !! 83 INTEGER :: jb, jk ! dummy loop indices 84 INTEGER :: ii, ij, igrd ! local integers 85 REAL(wp) :: zwgt ! boundary weight 86 !!---------------------------------------------------------------------- 87 ! 88 ! 89 igrd = 2 ! Relaxation of zonal velocity 90 DO jb = 1, idx%nblen(igrd) 91 ii = idx%nbi(jb,igrd) 92 ij = idx%nbj(jb,igrd) 93 zwgt = idx%nbw(jb,igrd) 94 pu2d(ii,ij) = ( pu2d(ii,ij) + zwgt * ( dta%u2d(jb) - pu2d(ii,ij) ) ) * umask(ii,ij,1) 95 END DO 96 ! 97 igrd = 3 ! Relaxation of meridional velocity 98 DO jb = 1, idx%nblen(igrd) 99 ii = idx%nbi(jb,igrd) 100 ij = idx%nbj(jb,igrd) 101 zwgt = idx%nbw(jb,igrd) 102 pv2d(ii,ij) = ( pv2d(ii,ij) + zwgt * ( dta%v2d(jb) - pv2d(ii,ij) ) ) * vmask(ii,ij,1) 103 END DO 104 CALL lbc_lnk( pu2d, 'U', -1. ) 105 CALL lbc_lnk( pv2d, 'V', -1. ) ! Boundary points should be updated 106 ! 107 108 END SUBROUTINE obc_dyn2d_frs 109 110 111 SUBROUTINE obc_dyn2d_fla( idx, dta ) 112 !!---------------------------------------------------------------------- 113 !! *** SUBROUTINE obc_dyn2d_fla *** 70 114 !! 71 115 !! - Apply Flather boundary conditions on normal barotropic velocities … … 100 144 101 145 !!!!!!!!!!!! SOME WORK TO DO ON THE TIDES HERE !!!!!!!!!!!!! 146 147 !!! REPLACE spgu with nemo_wrk work space 102 148 103 149 ! Fill temporary array with ssh data (here spgu): … … 120 166 iip1 = ii - MIN( 0, INT( idx%flagu(jb) ) ) ! T pts i-indice outside the boundary 121 167 ! 122 zcorr = - idx%flagu(jb) * SQRT( grav * hur_e(ii, ij) ) * ( pssh(iim1, ij) - spgu(iip1,ij) )168 zcorr = - idx%flagu(jb) * SQRT( grav * phur(ii, ij) ) * ( pssh(iim1, ij) - spgu(iip1,ij) ) 123 169 !!$ zforc = dta%u2d(jb) + utide(jb) 124 170 zforc = dta%u2d(jb) 125 ua_e(ii,ij) = zforc + zcorr * umask(ii,ij,1)171 pu2d(ii,ij) = zforc + zcorr * umask(ii,ij,1) 126 172 END DO 127 173 ! … … 134 180 ijp1 = ij - MIN( 0, INT( idx%flagv(jb) ) ) ! T pts j-indice outside the boundary 135 181 ! 136 zcorr = - idx%flagv(jb) * SQRT( grav * hvr_e(ii, ij) ) * ( pssh(ii, ijm1) - spgu(ii,ijp1) )182 zcorr = - idx%flagv(jb) * SQRT( grav * phvr(ii, ij) ) * ( pssh(ii, ijm1) - spgu(ii,ijp1) ) 137 183 !!$ zforc = dta%v2d(jb) + vtide(jb) 138 184 zforc = dta%v2d(jb) 139 va_e(ii,ij) = zforc + zcorr * vmask(ii,ij,1)140 END DO 141 CALL lbc_lnk( ua_e, 'U', -1. ) ! Boundary points should be updated142 CALL lbc_lnk( va_e, 'V', -1. ) !185 pv2d(ii,ij) = zforc + zcorr * vmask(ii,ij,1) 186 END DO 187 CALL lbc_lnk( pu2d, 'U', -1. ) ! Boundary points should be updated 188 CALL lbc_lnk( pv2d, 'V', -1. ) ! 143 189 ! 144 190 ! -
branches/2011/UKMO_MERCATOR_obc_bdy_merge/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90
r2797 r2800 134 134 !!--------------------------------------------------------------------- 135 135 ! Note that shifting time to be centrered in the middle of sbc time step impacts only nsec_* variables of the calendar 136 IF( present(timeshift) ) THEN 136 IF( present(jit) ) THEN 137 ! ignore kn_fsbc in this case 138 isecsbc = nsec_year + nsec1jan000 + jit*rdt/REAL(nn_baro,wp) 139 ELSE IF( present(timeshift) ) THEN 137 140 isecsbc = nsec_year + nsec1jan000 + NINT(0.5 * REAL(kn_fsbc - 1,wp) * rdttra(1)) + timeshift * rdttra(1) ! middle of sbc time step 138 141 ELSE
Note: See TracChangeset
for help on using the changeset viewer.