- Timestamp:
- 2019-11-26T15:11:43+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/ENHANCE-02_ISF_nemo_TEST_MERGE/src/OCE/BDY/bdydta.F90
r11229 r11967 43 43 PUBLIC bdy_dta_init ! routine called by nemogcm.F90 44 44 45 INTEGER, ALLOCATABLE, DIMENSION(:) :: nb_bdy_fld ! Number of fields to update for each boundary set. 46 INTEGER :: nb_bdy_fld_sum ! Total number of fields to update for all boundary sets. 47 LOGICAL, DIMENSION(jp_bdy) :: ln_full_vel_array ! =T => full velocities in 3D boundary conditions 48 ! =F => baroclinic velocities in 3D boundary conditions 45 INTEGER , PARAMETER :: jpbdyfld = 16 ! maximum number of files to read 46 INTEGER , PARAMETER :: jp_bdyssh = 1 ! 47 INTEGER , PARAMETER :: jp_bdyu2d = 2 ! 48 INTEGER , PARAMETER :: jp_bdyv2d = 3 ! 49 INTEGER , PARAMETER :: jp_bdyu3d = 4 ! 50 INTEGER , PARAMETER :: jp_bdyv3d = 5 ! 51 INTEGER , PARAMETER :: jp_bdytem = 6 ! 52 INTEGER , PARAMETER :: jp_bdysal = 7 ! 53 INTEGER , PARAMETER :: jp_bdya_i = 8 ! 54 INTEGER , PARAMETER :: jp_bdyh_i = 9 ! 55 INTEGER , PARAMETER :: jp_bdyh_s = 10 ! 56 INTEGER , PARAMETER :: jp_bdyt_i = 11 ! 57 INTEGER , PARAMETER :: jp_bdyt_s = 12 ! 58 INTEGER , PARAMETER :: jp_bdytsu = 13 ! 59 INTEGER , PARAMETER :: jp_bdys_i = 14 ! 60 INTEGER , PARAMETER :: jp_bdyaip = 15 ! 61 INTEGER , PARAMETER :: jp_bdyhip = 16 ! 62 #if ! defined key_si3 63 INTEGER , PARAMETER :: jpl = 1 64 #endif 65 49 66 !$AGRIF_DO_NOT_TREAT 50 TYPE(FLD), PUBLIC, ALLOCATABLE, DIMENSION(: ), TARGET :: bf! structure of input fields (file informations, fields read)67 TYPE(FLD), PUBLIC, ALLOCATABLE, DIMENSION(:,:), TARGET :: bf ! structure of input fields (file informations, fields read) 51 68 !$AGRIF_END_DO_NOT_TREAT 52 TYPE(MAP_POINTER), ALLOCATABLE, DIMENSION(:) :: nbmap_ptr ! array of pointers to nbmap53 54 #if defined key_si355 INTEGER :: nice_cat ! number of categories in the input file56 INTEGER :: jfld_hti, jfld_hts, jfld_ai ! indices of ice thickness, snow thickness and concentration in bf structure57 INTEGER, DIMENSION(jp_bdy) :: jfld_htit, jfld_htst, jfld_ait58 #endif59 69 60 70 !!---------------------------------------------------------------------- … … 65 75 CONTAINS 66 76 67 SUBROUTINE bdy_dta( kt, jit, time_offset )77 SUBROUTINE bdy_dta( kt, kit, kt_offset ) 68 78 !!---------------------------------------------------------------------- 69 79 !! *** SUBROUTINE bdy_dta *** … … 75 85 !!---------------------------------------------------------------------- 76 86 INTEGER, INTENT(in) :: kt ! ocean time-step index 77 INTEGER, INTENT(in), OPTIONAL :: jit ! subcycle time-step index (for timesplitting option)78 INTEGER, INTENT(in), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit87 INTEGER, INTENT(in), OPTIONAL :: kit ! subcycle time-step index (for timesplitting option) 88 INTEGER, INTENT(in), OPTIONAL :: kt_offset ! time offset in units of timesteps. NB. if kit 79 89 ! ! is present then units = subcycle timesteps. 80 ! ! time_offset = 0 => get data at "now" time level81 ! ! time_offset = -1 => get data at "before" time level82 ! ! time_offset = +1 => get data at "after" time level90 ! ! kt_offset = 0 => get data at "now" time level 91 ! ! kt_offset = -1 => get data at "before" time level 92 ! ! kt_offset = +1 => get data at "after" time level 83 93 ! ! etc. 84 94 ! 85 INTEGER :: jbdy, jfld, jstart, jend, ib, jl ! dummy loop indices 86 INTEGER :: ii, ij, ik, igrd ! local integers 87 INTEGER, DIMENSION(jpbgrd) :: ilen1 88 INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts 89 TYPE(OBC_DATA), POINTER :: dta ! short cut 95 INTEGER :: jbdy, jfld, jstart, jend, ib, jl ! dummy loop indices 96 INTEGER :: ii, ij, ik, igrd, ipl ! local integers 97 INTEGER, DIMENSION(jpbgrd) :: ilen1 98 INTEGER, DIMENSION(:), POINTER :: nblen, nblenrim ! short cuts 99 TYPE(OBC_DATA) , POINTER :: dta_alias ! short cut 100 TYPE(FLD), DIMENSION(:), POINTER :: bf_alias 90 101 !!--------------------------------------------------------------------------- 91 102 ! … … 94 105 ! Initialise data arrays once for all from initial conditions where required 95 106 !--------------------------------------------------------------------------- 96 IF( kt == nit000 .AND. .NOT.PRESENT( jit) ) THEN107 IF( kt == nit000 .AND. .NOT.PRESENT(kit) ) THEN 97 108 98 109 ! Calculate depth-mean currents 99 110 !----------------------------- 100 111 101 112 DO jbdy = 1, nb_bdy 102 113 ! 103 114 nblen => idx_bdy(jbdy)%nblen 104 115 nblenrim => idx_bdy(jbdy)%nblenrim 105 dta => dta_bdy(jbdy)106 116 ! 107 117 IF( nn_dyn2d_dta(jbdy) == 0 ) THEN 108 118 ilen1(:) = nblen(:) 109 IF( dta %ll_ssh ) THEN119 IF( dta_bdy(jbdy)%lneed_ssh ) THEN 110 120 igrd = 1 111 121 DO ib = 1, ilen1(igrd) … … 113 123 ij = idx_bdy(jbdy)%nbj(ib,igrd) 114 124 dta_bdy(jbdy)%ssh(ib) = sshn(ii,ij) * tmask(ii,ij,1) 115 END DO 116 ENDIF 117 IF( dta %ll_u2d) THEN125 END DO 126 ENDIF 127 IF( dta_bdy(jbdy)%lneed_dyn2d) THEN 118 128 igrd = 2 119 129 DO ib = 1, ilen1(igrd) … … 121 131 ij = idx_bdy(jbdy)%nbj(ib,igrd) 122 132 dta_bdy(jbdy)%u2d(ib) = un_b(ii,ij) * umask(ii,ij,1) 123 END DO 124 ENDIF 125 IF( dta%ll_v2d ) THEN 133 END DO 126 134 igrd = 3 127 135 DO ib = 1, ilen1(igrd) … … 129 137 ij = idx_bdy(jbdy)%nbj(ib,igrd) 130 138 dta_bdy(jbdy)%v2d(ib) = vn_b(ii,ij) * vmask(ii,ij,1) 131 END DO 139 END DO 132 140 ENDIF 133 141 ENDIF … … 135 143 IF( nn_dyn3d_dta(jbdy) == 0 ) THEN 136 144 ilen1(:) = nblen(:) 137 IF( dta %ll_u3d ) THEN145 IF( dta_bdy(jbdy)%lneed_dyn3d ) THEN 138 146 igrd = 2 139 147 DO ib = 1, ilen1(igrd) … … 143 151 dta_bdy(jbdy)%u3d(ib,ik) = ( un(ii,ij,ik) - un_b(ii,ij) ) * umask(ii,ij,ik) 144 152 END DO 145 END DO 146 ENDIF 147 IF( dta%ll_v3d ) THEN 153 END DO 148 154 igrd = 3 149 155 DO ib = 1, ilen1(igrd) … … 152 158 ij = idx_bdy(jbdy)%nbj(ib,igrd) 153 159 dta_bdy(jbdy)%v3d(ib,ik) = ( vn(ii,ij,ik) - vn_b(ii,ij) ) * vmask(ii,ij,ik) 154 155 END DO 160 END DO 161 END DO 156 162 ENDIF 157 163 ENDIF … … 159 165 IF( nn_tra_dta(jbdy) == 0 ) THEN 160 166 ilen1(:) = nblen(:) 161 IF( dta %ll_tem) THEN167 IF( dta_bdy(jbdy)%lneed_tra ) THEN 162 168 igrd = 1 163 169 DO ib = 1, ilen1(igrd) … … 165 171 ii = idx_bdy(jbdy)%nbi(ib,igrd) 166 172 ij = idx_bdy(jbdy)%nbj(ib,igrd) 167 dta_bdy(jbdy)%tem(ib,ik) = tsn(ii,ij,ik,jp_tem) * tmask(ii,ij,ik) 173 dta_bdy(jbdy)%tem(ib,ik) = tsn(ii,ij,ik,jp_bdytem) * tmask(ii,ij,ik) 174 dta_bdy(jbdy)%sal(ib,ik) = tsn(ii,ij,ik,jp_bdysal) * tmask(ii,ij,ik) 168 175 END DO 169 END DO 170 ENDIF 171 IF( dta%ll_sal ) THEN 172 igrd = 1 173 DO ib = 1, ilen1(igrd) 174 DO ik = 1, jpkm1 175 ii = idx_bdy(jbdy)%nbi(ib,igrd) 176 ij = idx_bdy(jbdy)%nbj(ib,igrd) 177 dta_bdy(jbdy)%sal(ib,ik) = tsn(ii,ij,ik,jp_sal) * tmask(ii,ij,ik) 178 END DO 179 END DO 176 END DO 180 177 ENDIF 181 178 ENDIF … … 184 181 IF( nn_ice_dta(jbdy) == 0 ) THEN ! set ice to initial values 185 182 ilen1(:) = nblen(:) 186 IF( dta %ll_a_i) THEN183 IF( dta_bdy(jbdy)%lneed_ice ) THEN 187 184 igrd = 1 188 185 DO jl = 1, jpl … … 190 187 ii = idx_bdy(jbdy)%nbi(ib,igrd) 191 188 ij = idx_bdy(jbdy)%nbj(ib,igrd) 192 dta_bdy(jbdy)%a_i (ib,jl) = a_i(ii,ij,jl) * tmask(ii,ij,1) 193 END DO 194 END DO 195 ENDIF 196 IF( dta%ll_h_i ) THEN 197 igrd = 1 198 DO jl = 1, jpl 199 DO ib = 1, ilen1(igrd) 200 ii = idx_bdy(jbdy)%nbi(ib,igrd) 201 ij = idx_bdy(jbdy)%nbj(ib,igrd) 202 dta_bdy(jbdy)%h_i (ib,jl) = h_i(ii,ij,jl) * tmask(ii,ij,1) 203 END DO 204 END DO 205 ENDIF 206 IF( dta%ll_h_s ) THEN 207 igrd = 1 208 DO jl = 1, jpl 209 DO ib = 1, ilen1(igrd) 210 ii = idx_bdy(jbdy)%nbi(ib,igrd) 211 ij = idx_bdy(jbdy)%nbj(ib,igrd) 212 dta_bdy(jbdy)%h_s (ib,jl) = h_s(ii,ij,jl) * tmask(ii,ij,1) 189 dta_bdy(jbdy)%a_i(ib,jl) = a_i (ii,ij,jl) * tmask(ii,ij,1) 190 dta_bdy(jbdy)%h_i(ib,jl) = h_i (ii,ij,jl) * tmask(ii,ij,1) 191 dta_bdy(jbdy)%h_s(ib,jl) = h_s (ii,ij,jl) * tmask(ii,ij,1) 192 dta_bdy(jbdy)%t_i(ib,jl) = SUM(t_i (ii,ij,:,jl)) * r1_nlay_i * tmask(ii,ij,1) 193 dta_bdy(jbdy)%t_s(ib,jl) = SUM(t_s (ii,ij,:,jl)) * r1_nlay_s * tmask(ii,ij,1) 194 dta_bdy(jbdy)%tsu(ib,jl) = t_su(ii,ij,jl) * tmask(ii,ij,1) 195 dta_bdy(jbdy)%s_i(ib,jl) = s_i (ii,ij,jl) * tmask(ii,ij,1) 196 ! melt ponds 197 dta_bdy(jbdy)%aip(ib,jl) = a_ip(ii,ij,jl) * tmask(ii,ij,1) 198 dta_bdy(jbdy)%hip(ib,jl) = h_ip(ii,ij,jl) * tmask(ii,ij,1) 213 199 END DO 214 200 END DO … … 222 208 ! update external data from files 223 209 !-------------------------------- 224 225 jstart = 1 226 DO jbdy = 1, nb_bdy 227 dta => dta_bdy(jbdy) 228 IF( nn_dta(jbdy) == 1 ) THEN ! skip this bit if no external data required 229 230 IF( PRESENT(jit) ) THEN 231 ! Update barotropic boundary conditions only 232 ! jit is optional argument for fld_read and bdytide_update 233 IF( cn_dyn2d(jbdy) /= 'none' ) THEN 234 IF( nn_dyn2d_dta(jbdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 235 IF( dta%ll_ssh ) dta%ssh(:) = 0._wp 236 IF( dta%ll_u2d ) dta%u2d(:) = 0._wp 237 IF( dta%ll_u3d ) dta%v2d(:) = 0._wp 238 ENDIF 239 IF (cn_tra(jbdy) /= 'runoff') THEN 240 IF( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 ) THEN 241 242 jend = jstart + dta%nread(2) - 1 243 IF( ln_full_vel_array(jbdy) ) THEN 244 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), & 245 & kit=jit, kt_offset=time_offset , jpk_bdy=nb_jpk_bdy(jbdy), & 246 & fvl=ln_full_vel_array(jbdy) ) 247 ELSE 248 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), & 249 & kit=jit, kt_offset=time_offset ) 250 ENDIF 251 252 ! If full velocities in boundary data then extract barotropic velocities from 3D fields 253 IF( ln_full_vel_array(jbdy) .AND. & 254 & ( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 .OR. & 255 & nn_dyn3d_dta(jbdy) == 1 ) )THEN 256 257 igrd = 2 ! zonal velocity 258 dta%u2d(:) = 0._wp 259 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 260 ii = idx_bdy(jbdy)%nbi(ib,igrd) 261 ij = idx_bdy(jbdy)%nbj(ib,igrd) 262 DO ik = 1, jpkm1 263 dta%u2d(ib) = dta%u2d(ib) & 264 & + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik) 265 END DO 266 dta%u2d(ib) = dta%u2d(ib) * r1_hu_n(ii,ij) 267 END DO 268 igrd = 3 ! meridional velocity 269 dta%v2d(:) = 0._wp 270 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 271 ii = idx_bdy(jbdy)%nbi(ib,igrd) 272 ij = idx_bdy(jbdy)%nbj(ib,igrd) 273 DO ik = 1, jpkm1 274 dta%v2d(ib) = dta%v2d(ib) & 275 & + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik) 276 END DO 277 dta%v2d(ib) = dta%v2d(ib) * r1_hv_n(ii,ij) 278 END DO 279 ENDIF 280 ENDIF 281 IF( nn_dyn2d_dta(jbdy) .ge. 2 ) THEN ! update tidal harmonic forcing 282 CALL bdytide_update( kt=kt, idx=idx_bdy(jbdy), dta=dta, td=tides(jbdy), & 283 & jit=jit, time_offset=time_offset ) 284 ENDIF 285 ENDIF 286 ENDIF 287 ELSE 288 IF (cn_tra(jbdy) == 'runoff') then ! runoff condition 289 jend = nb_bdy_fld(jbdy) 290 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 291 & map=nbmap_ptr(jstart:jend), kt_offset=time_offset ) 292 ! 293 igrd = 2 ! zonal velocity 294 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 295 ii = idx_bdy(jbdy)%nbi(ib,igrd) 296 ij = idx_bdy(jbdy)%nbj(ib,igrd) 297 dta%u2d(ib) = dta%u2d(ib) / ( e2u(ii,ij) * hu_0(ii,ij) ) 298 END DO 299 ! 300 igrd = 3 ! meridional velocity 301 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 302 ii = idx_bdy(jbdy)%nbi(ib,igrd) 303 ij = idx_bdy(jbdy)%nbj(ib,igrd) 304 dta%v2d(ib) = dta%v2d(ib) / ( e1v(ii,ij) * hv_0(ii,ij) ) 305 END DO 306 ELSE 307 IF( nn_dyn2d_dta(jbdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 308 IF( dta%ll_ssh ) dta%ssh(:) = 0._wp 309 IF( dta%ll_u2d ) dta%u2d(:) = 0._wp 310 IF( dta%ll_v2d ) dta%v2d(:) = 0._wp 311 ENDIF 312 IF( dta%nread(1) .gt. 0 ) THEN ! update external data 313 jend = jstart + dta%nread(1) - 1 314 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 315 & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy(jbdy), & 316 & fvl=ln_full_vel_array(jbdy) ) 317 ENDIF 318 ! If full velocities in boundary data then split into barotropic and baroclinic data 319 IF( ln_full_vel_array(jbdy) .and. & 320 & ( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 .OR. & 321 & nn_dyn3d_dta(jbdy) == 1 ) ) THEN 322 igrd = 2 ! zonal velocity 323 dta%u2d(:) = 0._wp 324 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 325 ii = idx_bdy(jbdy)%nbi(ib,igrd) 326 ij = idx_bdy(jbdy)%nbj(ib,igrd) 327 DO ik = 1, jpkm1 328 dta%u2d(ib) = dta%u2d(ib) & 329 & + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik) 330 END DO 331 dta%u2d(ib) = dta%u2d(ib) * r1_hu_n(ii,ij) 332 DO ik = 1, jpkm1 333 dta%u3d(ib,ik) = dta%u3d(ib,ik) - dta%u2d(ib) 334 END DO 335 END DO 336 igrd = 3 ! meridional velocity 337 dta%v2d(:) = 0._wp 338 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 339 ii = idx_bdy(jbdy)%nbi(ib,igrd) 340 ij = idx_bdy(jbdy)%nbj(ib,igrd) 341 DO ik = 1, jpkm1 342 dta%v2d(ib) = dta%v2d(ib) & 343 & + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik) 344 END DO 345 dta%v2d(ib) = dta%v2d(ib) * r1_hv_n(ii,ij) 346 DO ik = 1, jpkm1 347 dta%v3d(ib,ik) = dta%v3d(ib,ik) - dta%v2d(ib) 348 END DO 349 END DO 350 ENDIF 351 352 ENDIF 210 211 DO jbdy = 1, nb_bdy 212 213 dta_alias => dta_bdy(jbdy) 214 bf_alias => bf(:,jbdy) 215 216 ! read/update all bdy data 217 ! ------------------------ 218 CALL fld_read( kt, 1, bf_alias, kit = kit, kt_offset = kt_offset ) 219 220 ! apply some corrections in some specific cases... 221 ! -------------------------------------------------- 222 ! 223 ! if runoff condition: change river flow we read (in m3/s) into barotropic velocity (m/s) 224 IF( cn_tra(jbdy) == 'runoff' .AND. TRIM(bf_alias(jp_bdyu2d)%clrootname) /= 'NOT USED' ) THEN ! runoff and we read u/v2d 225 ! 226 igrd = 2 ! zonal flow (m3/s) to barotropic zonal velocity (m/s) 227 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 228 ii = idx_bdy(jbdy)%nbi(ib,igrd) 229 ij = idx_bdy(jbdy)%nbj(ib,igrd) 230 dta_alias%u2d(ib) = dta_alias%u2d(ib) / ( e2u(ii,ij) * hu_0(ii,ij) ) 231 END DO 232 igrd = 3 ! meridional flow (m3/s) to barotropic meridional velocity (m/s) 233 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 234 ii = idx_bdy(jbdy)%nbi(ib,igrd) 235 ij = idx_bdy(jbdy)%nbj(ib,igrd) 236 dta_alias%v2d(ib) = dta_alias%v2d(ib) / ( e1v(ii,ij) * hv_0(ii,ij) ) 237 END DO 238 ENDIF 239 240 ! tidal harmonic forcing ONLY: initialise arrays 241 IF( nn_dyn2d_dta(jbdy) == 2 ) THEN ! we did not read ssh, u/v2d 242 IF( dta_alias%lneed_ssh ) dta_alias%ssh(:) = 0._wp 243 IF( dta_alias%lneed_dyn2d ) dta_alias%u2d(:) = 0._wp 244 IF( dta_alias%lneed_dyn2d ) dta_alias%v2d(:) = 0._wp 245 ENDIF 246 247 ! If full velocities in boundary data, then split it into barotropic and baroclinic component 248 IF( bf_alias(jp_bdyu3d)%ltotvel ) THEN ! if we read 3D total velocity (can be true only if u3d was read) 249 ! 250 igrd = 2 ! zonal velocity 251 dta_alias%u2d(:) = 0._wp ! compute barotrope zonal velocity and put it in u2d 252 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 253 ii = idx_bdy(jbdy)%nbi(ib,igrd) 254 ij = idx_bdy(jbdy)%nbj(ib,igrd) 255 DO ik = 1, jpkm1 256 dta_alias%u2d(ib) = dta_alias%u2d(ib) + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta_alias%u3d(ib,ik) 257 END DO 258 dta_alias%u2d(ib) = dta_alias%u2d(ib) * r1_hu_n(ii,ij) 259 DO ik = 1, jpkm1 ! compute barocline zonal velocity and put it in u3d 260 dta_alias%u3d(ib,ik) = dta_alias%u3d(ib,ik) - dta_alias%u2d(ib) 261 END DO 262 END DO 263 igrd = 3 ! meridional velocity 264 dta_alias%v2d(:) = 0._wp ! compute barotrope meridional velocity and put it in v2d 265 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 266 ii = idx_bdy(jbdy)%nbi(ib,igrd) 267 ij = idx_bdy(jbdy)%nbj(ib,igrd) 268 DO ik = 1, jpkm1 269 dta_alias%v2d(ib) = dta_alias%v2d(ib) + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta_alias%v3d(ib,ik) 270 END DO 271 dta_alias%v2d(ib) = dta_alias%v2d(ib) * r1_hv_n(ii,ij) 272 DO ik = 1, jpkm1 ! compute barocline meridional velocity and put it in v3d 273 dta_alias%v3d(ib,ik) = dta_alias%v3d(ib,ik) - dta_alias%v2d(ib) 274 END DO 275 END DO 276 ENDIF ! ltotvel 277 278 ! update tidal harmonic forcing 279 IF( PRESENT(kit) .AND. nn_dyn2d_dta(jbdy) .GE. 2 ) THEN 280 CALL bdytide_update( kt = kt, idx = idx_bdy(jbdy), dta = dta_alias, td = tides(jbdy), & 281 & kit = kit, kt_offset = kt_offset ) 282 ENDIF 283 284 ! atm surface pressure : add inverted barometer effect to ssh if it was read 285 IF ( ln_apr_obc .AND. TRIM(bf_alias(jp_bdyssh)%clrootname) /= 'NOT USED' ) THEN 286 igrd = 1 287 DO ib = 1, idx_bdy(jbdy)%nblenrim(igrd) ! ssh is used only on the rim 288 ii = idx_bdy(jbdy)%nbi(ib,igrd) 289 ij = idx_bdy(jbdy)%nbj(ib,igrd) 290 dta_alias%ssh(ib) = dta_alias%ssh(ib) + ssh_ib(ii,ij) 291 END DO 292 ENDIF 293 353 294 #if defined key_si3 354 ! convert N-cat fields (input) into jpl-cat (output) 355 IF( cn_ice(jbdy) /= 'none' .AND. nn_ice_dta(jbdy) == 1 ) THEN 356 jfld_hti = jfld_htit(jbdy) 357 jfld_hts = jfld_htst(jbdy) 358 jfld_ai = jfld_ait(jbdy) 359 CALL ice_var_itd( bf(jfld_hti)%fnow(:,1,:), bf(jfld_hts)%fnow(:,1,:), bf(jfld_ai)%fnow(:,1,:), & 360 & dta_bdy(jbdy)%h_i , dta_bdy(jbdy)%h_s , dta_bdy(jbdy)%a_i ) 361 ENDIF 295 IF( dta_alias%lneed_ice ) THEN 296 ! fill temperature and salinity arrays 297 IF( TRIM(bf_alias(jp_bdyt_i)%clrootname) == 'NOT USED' ) bf_alias(jp_bdyt_i)%fnow(:,1,:) = rice_tem (jbdy) 298 IF( TRIM(bf_alias(jp_bdyt_s)%clrootname) == 'NOT USED' ) bf_alias(jp_bdyt_s)%fnow(:,1,:) = rice_tem (jbdy) 299 IF( TRIM(bf_alias(jp_bdytsu)%clrootname) == 'NOT USED' ) bf_alias(jp_bdytsu)%fnow(:,1,:) = rice_tem (jbdy) 300 IF( TRIM(bf_alias(jp_bdys_i)%clrootname) == 'NOT USED' ) bf_alias(jp_bdys_i)%fnow(:,1,:) = rice_sal (jbdy) 301 IF( TRIM(bf_alias(jp_bdyaip)%clrootname) == 'NOT USED' ) bf_alias(jp_bdyaip)%fnow(:,1,:) = rice_apnd(jbdy) * & ! rice_apnd is the pond fraction 302 & bf_alias(jp_bdya_i)%fnow(:,1,:) ! ( a_ip = rice_apnd * a_i ) 303 IF( TRIM(bf_alias(jp_bdyhip)%clrootname) == 'NOT USED' ) bf_alias(jp_bdyhip)%fnow(:,1,:) = rice_hpnd(jbdy) 304 ! if T_su is read and not T_i, set T_i = (T_su + T_freeze)/2 305 IF( TRIM(bf_alias(jp_bdytsu)%clrootname) /= 'NOT USED' .AND. TRIM(bf_alias(jp_bdyt_i)%clrootname) == 'NOT USED' ) & 306 & bf_alias(jp_bdyt_i)%fnow(:,1,:) = 0.5_wp * ( bf_alias(jp_bdytsu)%fnow(:,1,:) + 271.15 ) 307 ! if T_su is read and not T_s, set T_s = T_su 308 IF( TRIM(bf_alias(jp_bdytsu)%clrootname) /= 'NOT USED' .AND. TRIM(bf_alias(jp_bdyt_s)%clrootname) == 'NOT USED' ) & 309 & bf_alias(jp_bdyt_s)%fnow(:,1,:) = bf_alias(jp_bdytsu)%fnow(:,1,:) 310 ! if T_s is read and not T_su, set T_su = T_s 311 IF( TRIM(bf_alias(jp_bdyt_s)%clrootname) /= 'NOT USED' .AND. TRIM(bf_alias(jp_bdytsu)%clrootname) == 'NOT USED' ) & 312 & bf_alias(jp_bdytsu)%fnow(:,1,:) = bf_alias(jp_bdyt_s)%fnow(:,1,:) 313 ! if T_s is read and not T_i, set T_i = (T_s + T_freeze)/2 314 IF( TRIM(bf_alias(jp_bdyt_s)%clrootname) /= 'NOT USED' .AND. TRIM(bf_alias(jp_bdyt_i)%clrootname) == 'NOT USED' ) & 315 & bf_alias(jp_bdyt_i)%fnow(:,1,:) = 0.5_wp * ( bf_alias(jp_bdyt_s)%fnow(:,1,:) + 271.15 ) 316 317 ! make sure ponds = 0 if no ponds scheme 318 IF ( .NOT.ln_pnd ) THEN 319 bf_alias(jp_bdyaip)%fnow(:,1,:) = 0._wp 320 bf_alias(jp_bdyhip)%fnow(:,1,:) = 0._wp 321 ENDIF 322 323 ! convert N-cat fields (input) into jpl-cat (output) 324 ipl = SIZE(bf_alias(jp_bdya_i)%fnow, 3) 325 IF( ipl /= jpl ) THEN ! ice: convert N-cat fields (input) into jpl-cat (output) 326 CALL ice_var_itd( bf_alias(jp_bdyh_i)%fnow(:,1,:), bf_alias(jp_bdyh_s)%fnow(:,1,:), bf_alias(jp_bdya_i)%fnow(:,1,:), & 327 & dta_alias%h_i , dta_alias%h_s , dta_alias%a_i , & 328 & bf_alias(jp_bdyt_i)%fnow(:,1,:), bf_alias(jp_bdyt_s)%fnow(:,1,:), & 329 & bf_alias(jp_bdytsu)%fnow(:,1,:), bf_alias(jp_bdys_i)%fnow(:,1,:), & 330 & bf_alias(jp_bdyaip)%fnow(:,1,:), bf_alias(jp_bdyhip)%fnow(:,1,:), & 331 & dta_alias%t_i , dta_alias%t_s , & 332 & dta_alias%tsu , dta_alias%s_i , & 333 & dta_alias%aip , dta_alias%hip ) 334 ENDIF 335 ENDIF 362 336 #endif 363 ENDIF364 jstart = jstart + dta%nread(1)365 ENDIF ! nn_dta(jbdy) = 1366 337 END DO ! jbdy 367 368 IF ( ln_apr_obc ) THEN369 DO jbdy = 1, nb_bdy370 IF (cn_tra(jbdy) /= 'runoff')THEN371 igrd = 1 ! meridional velocity372 DO ib = 1, idx_bdy(jbdy)%nblenrim(igrd)373 ii = idx_bdy(jbdy)%nbi(ib,igrd)374 ij = idx_bdy(jbdy)%nbj(ib,igrd)375 dta_bdy(jbdy)%ssh(ib) = dta_bdy(jbdy)%ssh(ib) + ssh_ib(ii,ij)376 END DO377 ENDIF378 END DO379 ENDIF380 338 381 339 IF ( ln_tide ) THEN 382 340 IF (ln_dynspg_ts) THEN ! Fill temporary arrays with slow-varying bdy data 383 DO jbdy = 1, nb_bdy ! Tidal component added in ts loop384 IF ( nn_dyn2d_dta(jbdy) . ge. 2 ) THEN341 DO jbdy = 1, nb_bdy ! Tidal component added in ts loop 342 IF ( nn_dyn2d_dta(jbdy) .GE. 2 ) THEN 385 343 nblen => idx_bdy(jbdy)%nblen 386 344 nblenrim => idx_bdy(jbdy)%nblenrim 387 IF( cn_dyn2d(jbdy) == 'frs' ) THEN; ilen1(:)=nblen(:) ; ELSE ; ilen1(:)=nblenrim(:) ; ENDIF 388 IF ( dta_bdy(jbdy)%ll_ssh ) dta_bdy_s(jbdy)%ssh(1:ilen1(1)) = dta_bdy(jbdy)%ssh(1:ilen1(1)) 389 IF ( dta_bdy(jbdy)%ll_u2d ) dta_bdy_s(jbdy)%u2d(1:ilen1(2)) = dta_bdy(jbdy)%u2d(1:ilen1(2)) 390 IF ( dta_bdy(jbdy)%ll_v2d ) dta_bdy_s(jbdy)%v2d(1:ilen1(3)) = dta_bdy(jbdy)%v2d(1:ilen1(3)) 391 ENDIF 392 END DO 393 ELSE ! Add tides if not split-explicit free surface else this is done in ts loop 394 ! 395 CALL bdy_dta_tides( kt=kt, time_offset=time_offset ) 396 ENDIF 397 ENDIF 398 399 ! 400 IF( ln_timing ) CALL timing_stop('bdy_dta') 401 ! 402 END SUBROUTINE bdy_dta 345 IF( cn_dyn2d(jbdy) == 'frs' ) THEN ; ilen1(:)=nblen(:) ; ELSE ; ilen1(:)=nblenrim(:) ; ENDIF 346 IF ( dta_bdy(jbdy)%lneed_ssh ) dta_bdy_s(jbdy)%ssh(1:ilen1(1)) = dta_bdy(jbdy)%ssh(1:ilen1(1)) 347 IF ( dta_bdy(jbdy)%lneed_dyn2d ) dta_bdy_s(jbdy)%u2d(1:ilen1(2)) = dta_bdy(jbdy)%u2d(1:ilen1(2)) 348 IF ( dta_bdy(jbdy)%lneed_dyn2d ) dta_bdy_s(jbdy)%v2d(1:ilen1(3)) = dta_bdy(jbdy)%v2d(1:ilen1(3)) 349 ENDIF 350 END DO 351 ELSE ! Add tides if not split-explicit free surface else this is done in ts loop 352 ! 353 CALL bdy_dta_tides( kt=kt, kt_offset=kt_offset ) 354 ENDIF 355 ENDIF 356 ! 357 IF( ln_timing ) CALL timing_stop('bdy_dta') 358 ! 359 END SUBROUTINE bdy_dta 403 360 404 361 … … 413 370 !! 414 371 !!---------------------------------------------------------------------- 415 INTEGER :: jbdy, jfld, jstart, jend, ierror, ios ! Local integers 372 INTEGER :: jbdy, jfld ! Local integers 373 INTEGER :: ierror, ios ! 416 374 ! 375 CHARACTER(len=3) :: cl3 ! 417 376 CHARACTER(len=100) :: cn_dir ! Root directory for location of data files 418 CHARACTER(len=100), DIMENSION(nb_bdy) :: cn_dir_array ! Root directory for location of data files419 CHARACTER(len = 256):: clname ! temporary file name420 377 LOGICAL :: ln_full_vel ! =T => full velocities in 3D boundary data 421 378 ! ! =F => baroclinic velocities in 3D boundary data 422 INTEGER :: ilen_global ! Max length required for global bdy dta arrays423 INTEGER, ALLOCATABLE, DIMENSION(:) :: ilen1, ilen3 ! size of 1st and 3rd dimensions of local arrays424 INTEGER , ALLOCATABLE, DIMENSION(:) :: ibdy ! bdy set for a particular jfld425 INTEGER , ALLOCATABLE, DIMENSION(:) :: igrid ! index for grid type (1,2,3 = T,U,V)426 INTEGER , POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts427 TYPE(OBC_DATA), POINTER :: dta ! short cut428 #if defined key_si3 429 INTEGER :: kndims ! number of dimensions in an array (i.e. 3 = wo ice cat; 4 = w ice cat)430 INTEGER, DIMENSION(4) :: kdimsz ! size of dimensions431 INTEGER :: inum,id1 ! local integer432 #endif 433 TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: blf_i ! array of namelist information structures434 TYPE(FLD_N) :: bn_tem, bn_sal, bn_u3d, bn_v3d !435 TYPE(FLD_N) :: bn_ssh, bn_u2d, bn_v2d ! informations about the fields to be read436 #if defined key_si3 437 TYPE(FLD _N) :: bn_a_i, bn_h_i, bn_h_s438 #endif 379 LOGICAL :: ln_zinterp ! =T => requires a vertical interpolation of the bdydta 380 REAL(wp) :: rn_ice_tem, rn_ice_sal, rn_ice_age, rn_ice_apnd, rn_ice_hpnd 381 INTEGER :: ipk,ipl ! 382 INTEGER :: idvar ! variable ID 383 INTEGER :: indims ! number of dimensions of the variable 384 INTEGER :: iszdim ! number of dimensions of the variable 385 INTEGER, DIMENSION(4) :: i4dimsz ! size of variable dimensions 386 INTEGER :: igrd ! index for grid type (1,2,3 = T,U,V) 387 LOGICAL :: lluld ! is the variable using the unlimited dimension 388 LOGICAL :: llneed ! 389 LOGICAL :: llread ! 390 TYPE(FLD_N), DIMENSION(1), TARGET :: bn_tem, bn_sal, bn_u3d, bn_v3d ! must be an array to be used with fld_fill 391 TYPE(FLD_N), DIMENSION(1), TARGET :: bn_ssh, bn_u2d, bn_v2d ! informations about the fields to be read 392 TYPE(FLD_N), DIMENSION(1), TARGET :: bn_a_i, bn_h_i, bn_h_s, bn_t_i, bn_t_s, bn_tsu, bn_s_i, bn_aip, bn_hip 393 TYPE(FLD_N), DIMENSION(:), POINTER :: bn_alias ! must be an array to be used with fld_fill 394 TYPE(FLD ), DIMENSION(:), POINTER :: bf_alias 395 ! 439 396 NAMELIST/nambdy_dta/ cn_dir, bn_tem, bn_sal, bn_u3d, bn_v3d, bn_ssh, bn_u2d, bn_v2d 440 #if defined key_si3 441 NAMELIST/nambdy_dta/ bn_a_i, bn_h_i, bn_h_s 442 #endif 443 NAMELIST/nambdy_dta/ ln_full_vel 397 NAMELIST/nambdy_dta/ bn_a_i, bn_h_i, bn_h_s, bn_t_i, bn_t_s, bn_tsu, bn_s_i, bn_aip, bn_hip 398 NAMELIST/nambdy_dta/ rn_ice_tem, rn_ice_sal, rn_ice_age, rn_ice_apnd, rn_ice_hpnd 399 NAMELIST/nambdy_dta/ ln_full_vel, ln_zinterp 444 400 !!--------------------------------------------------------------------------- 445 401 ! … … 449 405 IF(lwp) WRITE(numout,*) '' 450 406 451 ! Set nn_dta 452 DO jbdy = 1, nb_bdy 453 nn_dta(jbdy) = MAX( nn_dyn2d_dta (jbdy) & 454 & , nn_dyn3d_dta (jbdy) & 455 & , nn_tra_dta (jbdy) & 456 #if defined key_si3 457 & , nn_ice_dta (jbdy) & 458 #endif 459 ) 460 IF(nn_dta(jbdy) > 1) nn_dta(jbdy) = 1 461 END DO 462 463 ! Work out upper bound of how many fields there are to read in and allocate arrays 464 ! --------------------------------------------------------------------------- 465 ALLOCATE( nb_bdy_fld(nb_bdy) ) 466 nb_bdy_fld(:) = 0 467 DO jbdy = 1, nb_bdy 468 IF( cn_dyn2d(jbdy) /= 'none' .AND. ( nn_dyn2d_dta(jbdy) == 1 .or. nn_dyn2d_dta(jbdy) == 3 ) ) THEN 469 nb_bdy_fld(jbdy) = nb_bdy_fld(jbdy) + 3 470 ENDIF 471 IF( cn_dyn3d(jbdy) /= 'none' .AND. nn_dyn3d_dta(jbdy) == 1 ) THEN 472 nb_bdy_fld(jbdy) = nb_bdy_fld(jbdy) + 2 473 ENDIF 474 IF( cn_tra(jbdy) /= 'none' .AND. nn_tra_dta(jbdy) == 1 ) THEN 475 nb_bdy_fld(jbdy) = nb_bdy_fld(jbdy) + 2 476 ENDIF 477 #if defined key_si3 478 IF( cn_ice(jbdy) /= 'none' .AND. nn_ice_dta(jbdy) == 1 ) THEN 479 nb_bdy_fld(jbdy) = nb_bdy_fld(jbdy) + 3 480 ENDIF 481 #endif 482 IF(lwp) WRITE(numout,*) 'Maximum number of files to open =', nb_bdy_fld(jbdy) 483 END DO 484 485 nb_bdy_fld_sum = SUM( nb_bdy_fld ) 486 487 ALLOCATE( bf(nb_bdy_fld_sum), STAT=ierror ) 407 ALLOCATE( bf(jpbdyfld,nb_bdy), STAT=ierror ) 488 408 IF( ierror > 0 ) THEN 489 409 CALL ctl_stop( 'bdy_dta: unable to allocate bf structure' ) ; RETURN 490 410 ENDIF 491 ALLOCATE( blf_i(nb_bdy_fld_sum), STAT=ierror ) 492 IF( ierror > 0 ) THEN 493 CALL ctl_stop( 'bdy_dta: unable to allocate blf_i structure' ) ; RETURN 494 ENDIF 495 ALLOCATE( nbmap_ptr(nb_bdy_fld_sum), STAT=ierror ) 496 IF( ierror > 0 ) THEN 497 CALL ctl_stop( 'bdy_dta: unable to allocate nbmap_ptr structure' ) ; RETURN 498 ENDIF 499 ALLOCATE( ilen1(nb_bdy_fld_sum), ilen3(nb_bdy_fld_sum) ) 500 ALLOCATE( ibdy(nb_bdy_fld_sum) ) 501 ALLOCATE( igrid(nb_bdy_fld_sum) ) 502 411 bf(:,:)%clrootname = 'NOT USED' ! default definition used as a flag in fld_read to do nothing. 412 bf(:,:)%lzint = .FALSE. ! default definition 413 bf(:,:)%ltotvel = .FALSE. ! default definition 414 503 415 ! Read namelists 504 416 ! -------------- 505 417 REWIND(numnam_cfg) 506 jfld = 0 507 DO jbdy = 1, nb_bdy 508 IF( nn_dta(jbdy) == 1 ) THEN 509 REWIND(numnam_ref) 510 READ ( numnam_ref, nambdy_dta, IOSTAT = ios, ERR = 901) 511 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy_dta in reference namelist', lwp ) 418 DO jbdy = 1, nb_bdy 419 420 WRITE(ctmp1, '(a,i2)') 'BDY number ', jbdy 421 WRITE(ctmp2, '(a,i2)') 'block nambdy_dta number ', jbdy 422 423 ! There is only one nambdy_dta block in namelist_ref -> use it for each bdy so we do a rewind 424 REWIND(numnam_ref) 425 READ ( numnam_ref, nambdy_dta, IOSTAT = ios, ERR = 901) 426 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy_dta in reference namelist' ) 427 428 ! by-pass nambdy_dta reading if no input data used in this bdy 429 IF( ( dta_bdy(jbdy)%lneed_dyn2d .AND. MOD(nn_dyn2d_dta(jbdy),2) == 1 ) & 430 & .OR. ( dta_bdy(jbdy)%lneed_dyn3d .AND. nn_dyn3d_dta(jbdy) == 1 ) & 431 & .OR. ( dta_bdy(jbdy)%lneed_tra .AND. nn_tra_dta(jbdy) == 1 ) & 432 & .OR. ( dta_bdy(jbdy)%lneed_ice .AND. nn_ice_dta(jbdy) == 1 ) ) THEN 433 ! WARNING: we don't do a rewind here, each bdy reads its own nambdy_dta block one after another 512 434 READ ( numnam_cfg, nambdy_dta, IOSTAT = ios, ERR = 902 ) 513 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nambdy_dta in configuration namelist', lwp ) 514 IF(lwm) WRITE( numond, nambdy_dta ) 515 516 cn_dir_array(jbdy) = cn_dir 517 ln_full_vel_array(jbdy) = ln_full_vel 518 519 nblen => idx_bdy(jbdy)%nblen 520 nblenrim => idx_bdy(jbdy)%nblenrim 521 dta => dta_bdy(jbdy) 522 dta%nread(2) = 0 523 524 ! Only read in necessary fields for this set. 525 ! Important that barotropic variables come first. 526 IF( nn_dyn2d_dta(jbdy) == 1 .or. nn_dyn2d_dta(jbdy) == 3 ) THEN 527 528 IF( dta%ll_ssh ) THEN 529 if(lwp) write(numout,*) '++++++ reading in ssh field' 530 jfld = jfld + 1 531 blf_i(jfld) = bn_ssh 532 ibdy(jfld) = jbdy 533 igrid(jfld) = 1 534 ilen1(jfld) = nblen(igrid(jfld)) 535 ilen3(jfld) = 1 536 dta%nread(2) = dta%nread(2) + 1 537 ENDIF 538 539 IF( dta%ll_u2d .and. .not. ln_full_vel_array(jbdy) ) THEN 540 if(lwp) write(numout,*) '++++++ reading in u2d field' 541 jfld = jfld + 1 542 blf_i(jfld) = bn_u2d 543 ibdy(jfld) = jbdy 544 igrid(jfld) = 2 545 ilen1(jfld) = nblen(igrid(jfld)) 546 ilen3(jfld) = 1 547 dta%nread(2) = dta%nread(2) + 1 548 ENDIF 549 550 IF( dta%ll_v2d .and. .not. ln_full_vel_array(jbdy) ) THEN 551 if(lwp) write(numout,*) '++++++ reading in v2d field' 552 jfld = jfld + 1 553 blf_i(jfld) = bn_v2d 554 ibdy(jfld) = jbdy 555 igrid(jfld) = 3 556 ilen1(jfld) = nblen(igrid(jfld)) 557 ilen3(jfld) = 1 558 dta%nread(2) = dta%nread(2) + 1 559 ENDIF 560 561 ENDIF 562 563 ! read 3D velocities if baroclinic velocities require OR if 564 ! barotropic velocities required and ln_full_vel set to .true. 565 IF( nn_dyn3d_dta(jbdy) == 1 .OR. & 566 & ( ln_full_vel_array(jbdy) .AND. ( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 ) ) ) THEN 567 568 IF( dta%ll_u3d .OR. ( ln_full_vel_array(jbdy) .and. dta%ll_u2d ) ) THEN 569 if(lwp) write(numout,*) '++++++ reading in u3d field' 570 jfld = jfld + 1 571 blf_i(jfld) = bn_u3d 572 ibdy(jfld) = jbdy 573 igrid(jfld) = 2 574 ilen1(jfld) = nblen(igrid(jfld)) 575 ilen3(jfld) = jpk 576 IF( ln_full_vel_array(jbdy) .and. dta%ll_u2d ) dta%nread(2) = dta%nread(2) + 1 577 ENDIF 578 579 IF( dta%ll_v3d .OR. ( ln_full_vel_array(jbdy) .and. dta%ll_v2d ) ) THEN 580 if(lwp) write(numout,*) '++++++ reading in v3d field' 581 jfld = jfld + 1 582 blf_i(jfld) = bn_v3d 583 ibdy(jfld) = jbdy 584 igrid(jfld) = 3 585 ilen1(jfld) = nblen(igrid(jfld)) 586 ilen3(jfld) = jpk 587 IF( ln_full_vel_array(jbdy) .and. dta%ll_v2d ) dta%nread(2) = dta%nread(2) + 1 588 ENDIF 589 590 ENDIF 591 592 ! temperature and salinity 593 IF( nn_tra_dta(jbdy) == 1 ) THEN 594 595 IF( dta%ll_tem ) THEN 596 if(lwp) write(numout,*) '++++++ reading in tem field' 597 jfld = jfld + 1 598 blf_i(jfld) = bn_tem 599 ibdy(jfld) = jbdy 600 igrid(jfld) = 1 601 ilen1(jfld) = nblen(igrid(jfld)) 602 ilen3(jfld) = jpk 603 ENDIF 604 605 IF( dta%ll_sal ) THEN 606 if(lwp) write(numout,*) '++++++ reading in sal field' 607 jfld = jfld + 1 608 blf_i(jfld) = bn_sal 609 ibdy(jfld) = jbdy 610 igrid(jfld) = 1 611 ilen1(jfld) = nblen(igrid(jfld)) 612 ilen3(jfld) = jpk 613 ENDIF 614 615 ENDIF 435 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nambdy_dta in configuration namelist' ) 436 IF(lwm) WRITE( numond, nambdy_dta ) 437 ENDIF 438 439 ! get the number of ice categories in bdy data file (use a_i information to do this) 440 ipl = jpl ! default definition 441 IF( dta_bdy(jbdy)%lneed_ice ) THEN ! if we need ice bdy data 442 IF( nn_ice_dta(jbdy) == 1 ) THEN ! if we get ice bdy data from netcdf file 443 CALL fld_fill( bf(jp_bdya_i,jbdy:jbdy), bn_a_i, cn_dir, 'bdy_dta', 'a_i'//' '//ctmp1, ctmp2 ) ! use namelist info 444 CALL fld_clopn( bf(jp_bdya_i,jbdy), nyear, nmonth, nday ) ! not a problem when we call it again after 445 idvar = iom_varid( bf(jp_bdya_i,jbdy)%num, bf(jp_bdya_i,jbdy)%clvar, kndims=indims, kdimsz=i4dimsz, lduld=lluld ) 446 IF( indims == 4 .OR. ( indims == 3 .AND. .NOT. lluld ) ) THEN ; ipl = i4dimsz(3) ! xylt or xyl 447 ELSE ; ipl = 1 ! xy or xyt 448 ENDIF 449 bf(jp_bdya_i,jbdy)%clrootname = 'NOT USED' ! reset to default value as this subdomain may not need to read this bdy 450 ENDIF 451 ENDIF 616 452 617 453 #if defined key_si3 618 ! sea ice 619 IF( nn_ice_dta(jbdy) == 1 ) THEN 620 ! Test for types of ice input (1cat or Xcat) 621 ! Build file name to find dimensions 622 clname=TRIM( cn_dir )//TRIM(bn_a_i%clname) 623 IF( .NOT. bn_a_i%ln_clim ) THEN 624 WRITE(clname, '(a,"_y",i4.4)' ) TRIM( clname ), nyear ! add year 625 IF( bn_a_i%cltype /= 'yearly' ) WRITE(clname, '(a,"m" ,i2.2)' ) TRIM( clname ), nmonth ! add month 626 ELSE 627 IF( bn_a_i%cltype /= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( clname ), nmonth ! add month 628 ENDIF 629 IF( bn_a_i%cltype == 'daily' .OR. bn_a_i%cltype(1:4) == 'week' ) & 630 & WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname ), nday ! add day 454 IF( .NOT.ln_pnd ) THEN 455 rn_ice_apnd = 0. ; rn_ice_hpnd = 0. 456 CALL ctl_warn( 'rn_ice_apnd & rn_ice_hpnd = 0 when no ponds' ) 457 ENDIF 458 #endif 459 460 ! temp, salt, age and ponds of incoming ice 461 rice_tem (jbdy) = rn_ice_tem 462 rice_sal (jbdy) = rn_ice_sal 463 rice_age (jbdy) = rn_ice_age 464 rice_apnd(jbdy) = rn_ice_apnd 465 rice_hpnd(jbdy) = rn_ice_hpnd 466 467 468 DO jfld = 1, jpbdyfld 469 470 ! ===================== 471 ! ssh 472 ! ===================== 473 IF( jfld == jp_bdyssh ) THEN 474 cl3 = 'ssh' 475 igrd = 1 ! T point 476 ipk = 1 ! surface data 477 llneed = dta_bdy(jbdy)%lneed_ssh ! dta_bdy(jbdy)%ssh will be needed 478 llread = MOD(nn_dyn2d_dta(jbdy),2) == 1 ! get data from NetCDF file 479 bf_alias => bf(jp_bdyssh,jbdy:jbdy) ! alias for ssh structure of bdy number jbdy 480 bn_alias => bn_ssh ! alias for ssh structure of nambdy_dta 481 iszdim = idx_bdy(jbdy)%nblenrim(igrd) ! length of this bdy on this MPI processus : used only on the rim 482 ENDIF 483 ! ===================== 484 ! dyn2d 485 ! ===================== 486 IF( jfld == jp_bdyu2d ) THEN 487 cl3 = 'u2d' 488 igrd = 2 ! U point 489 ipk = 1 ! surface data 490 llneed = dta_bdy(jbdy)%lneed_dyn2d ! dta_bdy(jbdy)%ssh will be needed 491 llread = .NOT. ln_full_vel .AND. MOD(nn_dyn2d_dta(jbdy),2) == 1 ! don't get u2d from u3d and read NetCDF file 492 bf_alias => bf(jp_bdyu2d,jbdy:jbdy) ! alias for u2d structure of bdy number jbdy 493 bn_alias => bn_u2d ! alias for u2d structure of nambdy_dta 494 IF( ln_full_vel ) THEN ; iszdim = idx_bdy(jbdy)%nblen(igrd) ! will be computed from u3d -> need on the full bdy 495 ELSE ; iszdim = idx_bdy(jbdy)%nblenrim(igrd) ! used only on the rim 496 ENDIF 497 ENDIF 498 IF( jfld == jp_bdyv2d ) THEN 499 cl3 = 'v2d' 500 igrd = 3 ! V point 501 ipk = 1 ! surface data 502 llneed = dta_bdy(jbdy)%lneed_dyn2d ! dta_bdy(jbdy)%ssh will be needed 503 llread = .NOT. ln_full_vel .AND. MOD(nn_dyn2d_dta(jbdy),2) == 1 ! don't get v2d from v3d and read NetCDF file 504 bf_alias => bf(jp_bdyv2d,jbdy:jbdy) ! alias for v2d structure of bdy number jbdy 505 bn_alias => bn_v2d ! alias for v2d structure of nambdy_dta 506 IF( ln_full_vel ) THEN ; iszdim = idx_bdy(jbdy)%nblen(igrd) ! will be computed from v3d -> need on the full bdy 507 ELSE ; iszdim = idx_bdy(jbdy)%nblenrim(igrd) ! used only on the rim 508 ENDIF 509 ENDIF 510 ! ===================== 511 ! dyn3d 512 ! ===================== 513 IF( jfld == jp_bdyu3d ) THEN 514 cl3 = 'u3d' 515 igrd = 2 ! U point 516 ipk = jpk ! 3d data 517 llneed = dta_bdy(jbdy)%lneed_dyn3d .OR. & ! dta_bdy(jbdy)%u3d will be needed 518 & ( dta_bdy(jbdy)%lneed_dyn2d .AND. ln_full_vel ) ! u3d needed to compute u2d 519 llread = nn_dyn3d_dta(jbdy) == 1 ! get data from NetCDF file 520 bf_alias => bf(jp_bdyu3d,jbdy:jbdy) ! alias for u3d structure of bdy number jbdy 521 bn_alias => bn_u3d ! alias for u3d structure of nambdy_dta 522 iszdim = idx_bdy(jbdy)%nblen(igrd) ! length of this bdy on this MPI processus 523 ENDIF 524 IF( jfld == jp_bdyv3d ) THEN 525 cl3 = 'v3d' 526 igrd = 3 ! V point 527 ipk = jpk ! 3d data 528 llneed = dta_bdy(jbdy)%lneed_dyn3d .OR. & ! dta_bdy(jbdy)%v3d will be needed 529 & ( dta_bdy(jbdy)%lneed_dyn2d .AND. ln_full_vel ) ! v3d needed to compute v2d 530 llread = nn_dyn3d_dta(jbdy) == 1 ! get data from NetCDF file 531 bf_alias => bf(jp_bdyv3d,jbdy:jbdy) ! alias for v3d structure of bdy number jbdy 532 bn_alias => bn_v3d ! alias for v3d structure of nambdy_dta 533 iszdim = idx_bdy(jbdy)%nblen(igrd) ! length of this bdy on this MPI processus 534 ENDIF 535 536 ! ===================== 537 ! tra 538 ! ===================== 539 IF( jfld == jp_bdytem ) THEN 540 cl3 = 'tem' 541 igrd = 1 ! T point 542 ipk = jpk ! 3d data 543 llneed = dta_bdy(jbdy)%lneed_tra ! dta_bdy(jbdy)%tem will be needed 544 llread = nn_tra_dta(jbdy) == 1 ! get data from NetCDF file 545 bf_alias => bf(jp_bdytem,jbdy:jbdy) ! alias for ssh structure of bdy number jbdy 546 bn_alias => bn_tem ! alias for ssh structure of nambdy_dta 547 iszdim = idx_bdy(jbdy)%nblen(igrd) ! length of this bdy on this MPI processus 548 ENDIF 549 IF( jfld == jp_bdysal ) THEN 550 cl3 = 'sal' 551 igrd = 1 ! T point 552 ipk = jpk ! 3d data 553 llneed = dta_bdy(jbdy)%lneed_tra ! dta_bdy(jbdy)%sal will be needed 554 llread = nn_tra_dta(jbdy) == 1 ! get data from NetCDF file 555 bf_alias => bf(jp_bdysal,jbdy:jbdy) ! alias for ssh structure of bdy number jbdy 556 bn_alias => bn_sal ! alias for ssh structure of nambdy_dta 557 iszdim = idx_bdy(jbdy)%nblen(igrd) ! length of this bdy on this MPI processus 558 ENDIF 559 560 ! ===================== 561 ! ice 562 ! ===================== 563 IF( jfld == jp_bdya_i .OR. jfld == jp_bdyh_i .OR. jfld == jp_bdyh_s .OR. & 564 & jfld == jp_bdyt_i .OR. jfld == jp_bdyt_s .OR. jfld == jp_bdytsu .OR. & 565 & jfld == jp_bdys_i .OR. jfld == jp_bdyaip .OR. jfld == jp_bdyhip ) THEN 566 igrd = 1 ! T point 567 ipk = ipl ! jpl-cat data 568 llneed = dta_bdy(jbdy)%lneed_ice ! ice will be needed 569 llread = nn_ice_dta(jbdy) == 1 ! get data from NetCDF file 570 iszdim = idx_bdy(jbdy)%nblen(igrd) ! length of this bdy on this MPI processus 571 ENDIF 572 IF( jfld == jp_bdya_i ) THEN 573 cl3 = 'a_i' 574 bf_alias => bf(jp_bdya_i,jbdy:jbdy) ! alias for a_i structure of bdy number jbdy 575 bn_alias => bn_a_i ! alias for a_i structure of nambdy_dta 576 ENDIF 577 IF( jfld == jp_bdyh_i ) THEN 578 cl3 = 'h_i' 579 bf_alias => bf(jp_bdyh_i,jbdy:jbdy) ! alias for h_i structure of bdy number jbdy 580 bn_alias => bn_h_i ! alias for h_i structure of nambdy_dta 581 ENDIF 582 IF( jfld == jp_bdyh_s ) THEN 583 cl3 = 'h_s' 584 bf_alias => bf(jp_bdyh_s,jbdy:jbdy) ! alias for h_s structure of bdy number jbdy 585 bn_alias => bn_h_s ! alias for h_s structure of nambdy_dta 586 ENDIF 587 IF( jfld == jp_bdyt_i ) THEN 588 cl3 = 't_i' 589 bf_alias => bf(jp_bdyt_i,jbdy:jbdy) ! alias for t_i structure of bdy number jbdy 590 bn_alias => bn_t_i ! alias for t_i structure of nambdy_dta 591 ENDIF 592 IF( jfld == jp_bdyt_s ) THEN 593 cl3 = 't_s' 594 bf_alias => bf(jp_bdyt_s,jbdy:jbdy) ! alias for t_s structure of bdy number jbdy 595 bn_alias => bn_t_s ! alias for t_s structure of nambdy_dta 596 ENDIF 597 IF( jfld == jp_bdytsu ) THEN 598 cl3 = 'tsu' 599 bf_alias => bf(jp_bdytsu,jbdy:jbdy) ! alias for tsu structure of bdy number jbdy 600 bn_alias => bn_tsu ! alias for tsu structure of nambdy_dta 601 ENDIF 602 IF( jfld == jp_bdys_i ) THEN 603 cl3 = 's_i' 604 bf_alias => bf(jp_bdys_i,jbdy:jbdy) ! alias for s_i structure of bdy number jbdy 605 bn_alias => bn_s_i ! alias for s_i structure of nambdy_dta 606 ENDIF 607 IF( jfld == jp_bdyaip ) THEN 608 cl3 = 'aip' 609 bf_alias => bf(jp_bdyaip,jbdy:jbdy) ! alias for aip structure of bdy number jbdy 610 bn_alias => bn_aip ! alias for aip structure of nambdy_dta 611 ENDIF 612 IF( jfld == jp_bdyhip ) THEN 613 cl3 = 'hip' 614 bf_alias => bf(jp_bdyhip,jbdy:jbdy) ! alias for hip structure of bdy number jbdy 615 bn_alias => bn_hip ! alias for hip structure of nambdy_dta 616 ENDIF 617 618 IF( llneed .AND. iszdim > 0 ) THEN ! dta_bdy(jbdy)%xxx will be needed 619 ! ! -> must be associated with an allocated target 620 ALLOCATE( bf_alias(1)%fnow( iszdim, 1, ipk ) ) ! allocate the target 631 621 ! 632 CALL iom_open ( clname, inum ) 633 id1 = iom_varid( inum, bn_a_i%clvar, kdimsz=kdimsz, kndims=kndims, ldstop = .FALSE. ) 634 CALL iom_close ( inum ) 635 636 IF ( kndims == 4 ) THEN 637 nice_cat = kdimsz(4) ! Xcat input 638 ELSE 639 nice_cat = 1 ! 1cat input 640 ENDIF 641 ! End test 642 643 IF( dta%ll_a_i ) THEN 644 jfld = jfld + 1 645 blf_i(jfld) = bn_a_i 646 ibdy(jfld) = jbdy 647 igrid(jfld) = 1 648 ilen1(jfld) = nblen(igrid(jfld)) 649 ilen3(jfld) = nice_cat 650 ENDIF 651 652 IF( dta%ll_h_i ) THEN 653 jfld = jfld + 1 654 blf_i(jfld) = bn_h_i 655 ibdy(jfld) = jbdy 656 igrid(jfld) = 1 657 ilen1(jfld) = nblen(igrid(jfld)) 658 ilen3(jfld) = nice_cat 659 ENDIF 660 661 IF( dta%ll_h_s ) THEN 662 jfld = jfld + 1 663 blf_i(jfld) = bn_h_s 664 ibdy(jfld) = jbdy 665 igrid(jfld) = 1 666 ilen1(jfld) = nblen(igrid(jfld)) 667 ilen3(jfld) = nice_cat 668 ENDIF 669 670 ENDIF 671 #endif 672 ! Recalculate field counts 673 !------------------------- 674 IF( jbdy == 1 ) THEN 675 nb_bdy_fld_sum = 0 676 nb_bdy_fld(jbdy) = jfld 677 nb_bdy_fld_sum = jfld 678 ELSE 679 nb_bdy_fld(jbdy) = jfld - nb_bdy_fld_sum 680 nb_bdy_fld_sum = nb_bdy_fld_sum + nb_bdy_fld(jbdy) 681 ENDIF 682 683 dta%nread(1) = nb_bdy_fld(jbdy) 684 685 ENDIF ! nn_dta == 1 686 ENDDO ! jbdy 687 688 DO jfld = 1, nb_bdy_fld_sum 689 ALLOCATE( bf(jfld)%fnow(ilen1(jfld),1,ilen3(jfld)) ) 690 IF( blf_i(jfld)%ln_tint ) ALLOCATE( bf(jfld)%fdta(ilen1(jfld),1,ilen3(jfld),2) ) 691 nbmap_ptr(jfld)%ptr => idx_bdy(ibdy(jfld))%nbmap(:,igrid(jfld)) 692 nbmap_ptr(jfld)%ll_unstruc = ln_coords_file(ibdy(jfld)) 693 ENDDO 694 695 ! fill bf with blf_i and control print 696 !------------------------------------- 697 jstart = 1 698 DO jbdy = 1, nb_bdy 699 jend = jstart - 1 + nb_bdy_fld(jbdy) 700 CALL fld_fill( bf(jstart:jend), blf_i(jstart:jend), cn_dir_array(jbdy), 'bdy_dta', & 701 & 'open boundary conditions', 'nambdy_dta' ) 702 jstart = jend + 1 703 ENDDO 704 705 DO jfld = 1, nb_bdy_fld_sum 706 bf(jfld)%igrd = igrid(jfld) 707 bf(jfld)%ibdy = ibdy(jfld) 708 ENDDO 709 710 ! Initialise local boundary data arrays 711 ! nn_xxx_dta=0 : allocate space - will be filled from initial conditions later 712 ! nn_xxx_dta=1 : point to "fnow" arrays 713 !------------------------------------- 714 715 jfld = 0 716 DO jbdy=1, nb_bdy 717 718 nblen => idx_bdy(jbdy)%nblen 719 dta => dta_bdy(jbdy) 720 721 if(lwp) then 722 write(numout,*) '++++++ dta%ll_ssh = ',dta%ll_ssh 723 write(numout,*) '++++++ dta%ll_u2d = ',dta%ll_u2d 724 write(numout,*) '++++++ dta%ll_v2d = ',dta%ll_v2d 725 write(numout,*) '++++++ dta%ll_u3d = ',dta%ll_u3d 726 write(numout,*) '++++++ dta%ll_v3d = ',dta%ll_v3d 727 write(numout,*) '++++++ dta%ll_tem = ',dta%ll_tem 728 write(numout,*) '++++++ dta%ll_sal = ',dta%ll_sal 729 endif 730 731 IF ( nn_dyn2d_dta(jbdy) == 0 .or. nn_dyn2d_dta(jbdy) == 2 ) THEN 732 if(lwp) write(numout,*) '++++++ dta%ssh/u2d/u3d allocated space' 733 IF( dta%ll_ssh ) ALLOCATE( dta%ssh(nblen(1)) ) 734 IF( dta%ll_u2d ) ALLOCATE( dta%u2d(nblen(2)) ) 735 IF( dta%ll_v2d ) ALLOCATE( dta%v2d(nblen(3)) ) 736 ENDIF 737 IF ( nn_dyn2d_dta(jbdy) == 1 .or. nn_dyn2d_dta(jbdy) == 3 ) THEN 738 IF( dta%ll_ssh ) THEN 739 if(lwp) write(numout,*) '++++++ dta%ssh pointing to fnow' 740 jfld = jfld + 1 741 dta%ssh => bf(jfld)%fnow(:,1,1) 742 ENDIF 743 IF ( dta%ll_u2d ) THEN 744 IF ( ln_full_vel_array(jbdy) ) THEN 745 if(lwp) write(numout,*) '++++++ dta%u2d allocated space' 746 ALLOCATE( dta%u2d(nblen(2)) ) 747 ELSE 748 if(lwp) write(numout,*) '++++++ dta%u2d pointing to fnow' 749 jfld = jfld + 1 750 dta%u2d => bf(jfld)%fnow(:,1,1) 751 ENDIF 752 ENDIF 753 IF ( dta%ll_v2d ) THEN 754 IF ( ln_full_vel_array(jbdy) ) THEN 755 if(lwp) write(numout,*) '++++++ dta%v2d allocated space' 756 ALLOCATE( dta%v2d(nblen(3)) ) 757 ELSE 758 if(lwp) write(numout,*) '++++++ dta%v2d pointing to fnow' 759 jfld = jfld + 1 760 dta%v2d => bf(jfld)%fnow(:,1,1) 761 ENDIF 762 ENDIF 763 ENDIF 764 765 IF ( nn_dyn3d_dta(jbdy) == 0 ) THEN 766 if(lwp) write(numout,*) '++++++ dta%u3d/v3d allocated space' 767 IF( dta%ll_u3d ) ALLOCATE( dta_bdy(jbdy)%u3d(nblen(2),jpk) ) 768 IF( dta%ll_v3d ) ALLOCATE( dta_bdy(jbdy)%v3d(nblen(3),jpk) ) 769 ENDIF 770 IF ( nn_dyn3d_dta(jbdy) == 1 .or. & 771 & ( ln_full_vel_array(jbdy) .and. ( nn_dyn2d_dta(jbdy) == 1 .or. nn_dyn2d_dta(jbdy) == 3 ) ) ) THEN 772 IF ( dta%ll_u3d .or. ( ln_full_vel_array(jbdy) .and. dta%ll_u2d ) ) THEN 773 if(lwp) write(numout,*) '++++++ dta%u3d pointing to fnow' 774 jfld = jfld + 1 775 dta_bdy(jbdy)%u3d => bf(jfld)%fnow(:,1,:) 776 ENDIF 777 IF ( dta%ll_v3d .or. ( ln_full_vel_array(jbdy) .and. dta%ll_v2d ) ) THEN 778 if(lwp) write(numout,*) '++++++ dta%v3d pointing to fnow' 779 jfld = jfld + 1 780 dta_bdy(jbdy)%v3d => bf(jfld)%fnow(:,1,:) 781 ENDIF 782 ENDIF 783 784 IF( nn_tra_dta(jbdy) == 0 ) THEN 785 if(lwp) write(numout,*) '++++++ dta%tem/sal allocated space' 786 IF( dta%ll_tem ) ALLOCATE( dta_bdy(jbdy)%tem(nblen(1),jpk) ) 787 IF( dta%ll_sal ) ALLOCATE( dta_bdy(jbdy)%sal(nblen(1),jpk) ) 788 ELSE 789 IF( dta%ll_tem ) THEN 790 if(lwp) write(numout,*) '++++++ dta%tem pointing to fnow' 791 jfld = jfld + 1 792 dta_bdy(jbdy)%tem => bf(jfld)%fnow(:,1,:) 793 ENDIF 794 IF( dta%ll_sal ) THEN 795 if(lwp) write(numout,*) '++++++ dta%sal pointing to fnow' 796 jfld = jfld + 1 797 dta_bdy(jbdy)%sal => bf(jfld)%fnow(:,1,:) 798 ENDIF 799 ENDIF 800 801 #if defined key_si3 802 IF (cn_ice(jbdy) /= 'none') THEN 803 IF( nn_ice_dta(jbdy) == 0 ) THEN 804 ALLOCATE( dta_bdy(jbdy)%a_i(nblen(1),jpl) ) 805 ALLOCATE( dta_bdy(jbdy)%h_i(nblen(1),jpl) ) 806 ALLOCATE( dta_bdy(jbdy)%h_s(nblen(1),jpl) ) 807 ELSE 808 IF ( nice_cat == jpl ) THEN ! case input cat = jpl 809 jfld = jfld + 1 810 dta_bdy(jbdy)%a_i => bf(jfld)%fnow(:,1,:) 811 jfld = jfld + 1 812 dta_bdy(jbdy)%h_i => bf(jfld)%fnow(:,1,:) 813 jfld = jfld + 1 814 dta_bdy(jbdy)%h_s => bf(jfld)%fnow(:,1,:) 815 ELSE ! case input cat = 1 OR (/=1 and /=jpl) 816 jfld_ait(jbdy) = jfld + 1 817 jfld_htit(jbdy) = jfld + 2 818 jfld_htst(jbdy) = jfld + 3 819 jfld = jfld + 3 820 ALLOCATE( dta_bdy(jbdy)%a_i(nblen(1),jpl) ) 821 ALLOCATE( dta_bdy(jbdy)%h_i(nblen(1),jpl) ) 822 ALLOCATE( dta_bdy(jbdy)%h_s(nblen(1),jpl) ) 823 dta_bdy(jbdy)%a_i(:,:) = 0._wp 824 dta_bdy(jbdy)%h_i(:,:) = 0._wp 825 dta_bdy(jbdy)%h_s(:,:) = 0._wp 826 ENDIF 827 828 ENDIF 829 ENDIF 830 #endif 622 IF( llread ) THEN ! get data from NetCDF file 623 CALL fld_fill( bf_alias, bn_alias, cn_dir, 'bdy_dta', cl3//' '//ctmp1, ctmp2 ) ! use namelist info 624 IF( bf_alias(1)%ln_tint ) ALLOCATE( bf_alias(1)%fdta( iszdim, 1, ipk, 2 ) ) 625 bf_alias(1)%imap => idx_bdy(jbdy)%nbmap(1:iszdim,igrd) ! associate the mapping used for this bdy 626 bf_alias(1)%igrd = igrd ! used only for vertical integration of 3D arrays 627 bf_alias(1)%ibdy = jbdy ! " " " " " " " " 628 bf_alias(1)%ltotvel = ln_full_vel ! T if u3d is full velocity 629 bf_alias(1)%lzint = ln_zinterp ! T if it requires a vertical interpolation 630 ENDIF 631 632 ! associate the pointer and get rid of the dimensions with a size equal to 1 633 IF( jfld == jp_bdyssh ) dta_bdy(jbdy)%ssh => bf_alias(1)%fnow(:,1,1) 634 IF( jfld == jp_bdyu2d ) dta_bdy(jbdy)%u2d => bf_alias(1)%fnow(:,1,1) 635 IF( jfld == jp_bdyv2d ) dta_bdy(jbdy)%v2d => bf_alias(1)%fnow(:,1,1) 636 IF( jfld == jp_bdyu3d ) dta_bdy(jbdy)%u3d => bf_alias(1)%fnow(:,1,:) 637 IF( jfld == jp_bdyv3d ) dta_bdy(jbdy)%v3d => bf_alias(1)%fnow(:,1,:) 638 IF( jfld == jp_bdytem ) dta_bdy(jbdy)%tem => bf_alias(1)%fnow(:,1,:) 639 IF( jfld == jp_bdysal ) dta_bdy(jbdy)%sal => bf_alias(1)%fnow(:,1,:) 640 IF( jfld == jp_bdya_i ) THEN 641 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%a_i => bf_alias(1)%fnow(:,1,:) 642 ELSE ; ALLOCATE( dta_bdy(jbdy)%a_i(iszdim,jpl) ) 643 ENDIF 644 ENDIF 645 IF( jfld == jp_bdyh_i ) THEN 646 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%h_i => bf_alias(1)%fnow(:,1,:) 647 ELSE ; ALLOCATE( dta_bdy(jbdy)%h_i(iszdim,jpl) ) 648 ENDIF 649 ENDIF 650 IF( jfld == jp_bdyh_s ) THEN 651 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%h_s => bf_alias(1)%fnow(:,1,:) 652 ELSE ; ALLOCATE( dta_bdy(jbdy)%h_s(iszdim,jpl) ) 653 ENDIF 654 ENDIF 655 IF( jfld == jp_bdyt_i ) THEN 656 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%t_i => bf_alias(1)%fnow(:,1,:) 657 ELSE ; ALLOCATE( dta_bdy(jbdy)%t_i(iszdim,jpl) ) 658 ENDIF 659 ENDIF 660 IF( jfld == jp_bdyt_s ) THEN 661 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%t_s => bf_alias(1)%fnow(:,1,:) 662 ELSE ; ALLOCATE( dta_bdy(jbdy)%t_s(iszdim,jpl) ) 663 ENDIF 664 ENDIF 665 IF( jfld == jp_bdytsu ) THEN 666 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%tsu => bf_alias(1)%fnow(:,1,:) 667 ELSE ; ALLOCATE( dta_bdy(jbdy)%tsu(iszdim,jpl) ) 668 ENDIF 669 ENDIF 670 IF( jfld == jp_bdys_i ) THEN 671 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%s_i => bf_alias(1)%fnow(:,1,:) 672 ELSE ; ALLOCATE( dta_bdy(jbdy)%s_i(iszdim,jpl) ) 673 ENDIF 674 ENDIF 675 IF( jfld == jp_bdyaip ) THEN 676 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%aip => bf_alias(1)%fnow(:,1,:) 677 ELSE ; ALLOCATE( dta_bdy(jbdy)%aip(iszdim,jpl) ) 678 ENDIF 679 ENDIF 680 IF( jfld == jp_bdyhip ) THEN 681 IF( ipk == jpl ) THEN ; dta_bdy(jbdy)%hip => bf_alias(1)%fnow(:,1,:) 682 ELSE ; ALLOCATE( dta_bdy(jbdy)%hip(iszdim,jpl) ) 683 ENDIF 684 ENDIF 685 ENDIF 686 687 END DO ! jpbdyfld 831 688 ! 832 689 END DO ! jbdy 833 690 ! 834 691 END SUBROUTINE bdy_dta_init 835 692 836 693 !!============================================================================== 837 694 END MODULE bdydta
Note: See TracChangeset
for help on using the changeset viewer.