- Timestamp:
- 2013-11-20T17:28:04+01:00 (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90
r4148 r4292 21 21 !! bdy_init : Initialization of unstructured open boundaries 22 22 !!---------------------------------------------------------------------- 23 USE wrk_nemo ! Memory Allocation 23 24 USE timing ! Timing 24 25 USE oce ! ocean dynamics and tracers variables … … 79 80 INTEGER :: jpbdtau, jpbdtas ! - - 80 81 INTEGER :: ib_bdy1, ib_bdy2, ib1, ib2 ! - - 82 INTEGER :: i_offset, j_offset ! - - 81 83 INTEGER, POINTER :: nbi, nbj, nbr ! short cuts 82 REAL , POINTER :: flagu, flagv ! - - 84 REAL(wp), POINTER :: flagu, flagv ! - - 85 REAL(wp), POINTER, DIMENSION(:,:) :: pmask ! pointer to 2D mask fields 83 86 REAL(wp) :: zefl, zwfl, znfl, zsfl ! local scalars 84 87 INTEGER, DIMENSION (2) :: kdimsz … … 90 93 INTEGER :: com_east_b, com_west_b, com_south_b, com_north_b ! Flags for boundaries receiving 91 94 INTEGER :: iw_b(4), ie_b(4), is_b(4), in_b(4) ! Arrays for neighbours coordinates 95 REAL(wp), POINTER, DIMENSION(:,:) :: zfmask ! temporary fmask array excluding coastal boundary condition (shlat) 92 96 93 97 !! 94 NAMELIST/nambdy/ nb_bdy, ln_coords_file, cn_coords_file, &95 & ln_mask_file, cn_mask_file, nn_dyn2d, nn_dyn2d_dta,&96 & nn_dyn3d, nn_dyn3d_dta, nn_tra, nn_tra_dta,&97 & ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, 98 #if defined key_lim299 & nn_ice_lim2, nn_ice_lim2_dta,&98 NAMELIST/nambdy/ nb_bdy, ln_coords_file, cn_coords_file, & 99 & ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta, & 100 & cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta, & 101 & ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, & 102 #if ( defined key_lim2 || defined key_lim3 ) 103 & cn_ice_lim, nn_ice_lim_dta, & 100 104 #endif 101 105 & ln_vol, nn_volctl, nn_rimwidth … … 156 160 157 161 IF(lwp) WRITE(numout,*) 'Boundary conditions for barotropic solution: ' 158 SELECT CASE( nn_dyn2d(ib_bdy) ) 159 CASE(jp_none) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 160 CASE(jp_frs) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 161 CASE(jp_flather) ; IF(lwp) WRITE(numout,*) ' Flather radiation condition' 162 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_dyn2d' ) 162 SELECT CASE( cn_dyn2d(ib_bdy) ) 163 CASE('none') 164 IF(lwp) WRITE(numout,*) ' no open boundary condition' 165 dta_bdy(ib_bdy)%ll_ssh = .false. 166 dta_bdy(ib_bdy)%ll_u2d = .false. 167 dta_bdy(ib_bdy)%ll_v2d = .false. 168 CASE('frs') 169 IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 170 dta_bdy(ib_bdy)%ll_ssh = .false. 171 dta_bdy(ib_bdy)%ll_u2d = .true. 172 dta_bdy(ib_bdy)%ll_v2d = .true. 173 CASE('flather') 174 IF(lwp) WRITE(numout,*) ' Flather radiation condition' 175 dta_bdy(ib_bdy)%ll_ssh = .true. 176 dta_bdy(ib_bdy)%ll_u2d = .true. 177 dta_bdy(ib_bdy)%ll_v2d = .true. 178 CASE('orlanski') 179 IF(lwp) WRITE(numout,*) ' Orlanski (fully oblique) radiation condition with adaptive nudging' 180 dta_bdy(ib_bdy)%ll_ssh = .false. 181 dta_bdy(ib_bdy)%ll_u2d = .true. 182 dta_bdy(ib_bdy)%ll_v2d = .true. 183 CASE('orlanski_npo') 184 IF(lwp) WRITE(numout,*) ' Orlanski (NPO) radiation condition with adaptive nudging' 185 dta_bdy(ib_bdy)%ll_ssh = .false. 186 dta_bdy(ib_bdy)%ll_u2d = .true. 187 dta_bdy(ib_bdy)%ll_v2d = .true. 188 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for cn_dyn2d' ) 163 189 END SELECT 164 IF( nn_dyn2d(ib_bdy) .gt. 0) THEN190 IF( cn_dyn2d(ib_bdy) /= 'none' ) THEN 165 191 SELECT CASE( nn_dyn2d_dta(ib_bdy) ) ! 166 192 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' initial state used for bdy data' … … 177 203 178 204 IF(lwp) WRITE(numout,*) 'Boundary conditions for baroclinic velocities: ' 179 SELECT CASE( nn_dyn3d(ib_bdy) ) 180 CASE(jp_none) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 181 CASE(jp_frs) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 182 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Specified value' 183 CASE( 3 ) ; IF(lwp) WRITE(numout,*) ' Zero baroclinic velocities (runoff case)' 184 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_dyn3d' ) 205 SELECT CASE( cn_dyn3d(ib_bdy) ) 206 CASE('none') 207 IF(lwp) WRITE(numout,*) ' no open boundary condition' 208 dta_bdy(ib_bdy)%ll_u3d = .false. 209 dta_bdy(ib_bdy)%ll_v3d = .false. 210 CASE('frs') 211 IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 212 dta_bdy(ib_bdy)%ll_u3d = .true. 213 dta_bdy(ib_bdy)%ll_v3d = .true. 214 CASE('specified') 215 IF(lwp) WRITE(numout,*) ' Specified value' 216 dta_bdy(ib_bdy)%ll_u3d = .true. 217 dta_bdy(ib_bdy)%ll_v3d = .true. 218 CASE('zero') 219 IF(lwp) WRITE(numout,*) ' Zero baroclinic velocities (runoff case)' 220 dta_bdy(ib_bdy)%ll_u3d = .false. 221 dta_bdy(ib_bdy)%ll_v3d = .false. 222 CASE('orlanski') 223 IF(lwp) WRITE(numout,*) ' Orlanski (fully oblique) radiation condition with adaptive nudging' 224 dta_bdy(ib_bdy)%ll_u3d = .true. 225 dta_bdy(ib_bdy)%ll_v3d = .true. 226 CASE('orlanski_npo') 227 IF(lwp) WRITE(numout,*) ' Orlanski (NPO) radiation condition with adaptive nudging' 228 dta_bdy(ib_bdy)%ll_u3d = .true. 229 dta_bdy(ib_bdy)%ll_v3d = .true. 230 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for cn_dyn3d' ) 185 231 END SELECT 186 IF( nn_dyn3d(ib_bdy) .gt. 0) THEN232 IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 187 233 SELECT CASE( nn_dyn3d_dta(ib_bdy) ) ! 188 234 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' initial state used for bdy data' … … 193 239 194 240 IF ( ln_dyn3d_dmp(ib_bdy) ) THEN 195 IF ( nn_dyn3d(ib_bdy).EQ.0) THEN241 IF ( cn_dyn3d(ib_bdy) == 'none' ) THEN 196 242 IF(lwp) WRITE(numout,*) 'No open boundary condition for baroclinic velocities: ln_dyn3d_dmp is set to .false.' 197 243 ln_dyn3d_dmp(ib_bdy)=.false. 198 ELSEIF ( nn_dyn3d(ib_bdy).EQ.1) THEN244 ELSEIF ( cn_dyn3d(ib_bdy) == 'frs' ) THEN 199 245 CALL ctl_stop( 'Use FRS OR relaxation' ) 200 246 ELSE … … 202 248 IF(lwp) WRITE(numout,*) ' Damping time scale: ',rn_time_dmp(ib_bdy),' days' 203 249 IF((lwp).AND.rn_time_dmp(ib_bdy)<0) CALL ctl_stop( 'Time scale must be positive' ) 250 dta_bdy(ib_bdy)%ll_u3d = .true. 251 dta_bdy(ib_bdy)%ll_v3d = .true. 204 252 ENDIF 205 253 ELSE … … 209 257 210 258 IF(lwp) WRITE(numout,*) 'Boundary conditions for temperature and salinity: ' 211 SELECT CASE( nn_tra(ib_bdy) ) 212 CASE(jp_none) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 213 CASE(jp_frs) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 214 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Specified value' 215 CASE( 3 ) ; IF(lwp) WRITE(numout,*) ' Neumann conditions' 216 CASE( 4 ) ; IF(lwp) WRITE(numout,*) ' Runoff conditions : Neumann for T and specified to 0.1 for salinity' 217 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_tra' ) 259 SELECT CASE( cn_tra(ib_bdy) ) 260 CASE('none') 261 IF(lwp) WRITE(numout,*) ' no open boundary condition' 262 dta_bdy(ib_bdy)%ll_tem = .false. 263 dta_bdy(ib_bdy)%ll_sal = .false. 264 CASE('frs') 265 IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 266 dta_bdy(ib_bdy)%ll_tem = .true. 267 dta_bdy(ib_bdy)%ll_sal = .true. 268 CASE('specified') 269 IF(lwp) WRITE(numout,*) ' Specified value' 270 dta_bdy(ib_bdy)%ll_tem = .true. 271 dta_bdy(ib_bdy)%ll_sal = .true. 272 CASE('neumann') 273 IF(lwp) WRITE(numout,*) ' Neumann conditions' 274 dta_bdy(ib_bdy)%ll_tem = .false. 275 dta_bdy(ib_bdy)%ll_sal = .false. 276 CASE('runoff') 277 IF(lwp) WRITE(numout,*) ' Runoff conditions : Neumann for T and specified to 0.1 for salinity' 278 dta_bdy(ib_bdy)%ll_tem = .false. 279 dta_bdy(ib_bdy)%ll_sal = .false. 280 CASE('orlanski') 281 IF(lwp) WRITE(numout,*) ' Orlanski (fully oblique) radiation condition with adaptive nudging' 282 dta_bdy(ib_bdy)%ll_tem = .true. 283 dta_bdy(ib_bdy)%ll_sal = .true. 284 CASE('orlanski_npo') 285 IF(lwp) WRITE(numout,*) ' Orlanski (NPO) radiation condition with adaptive nudging' 286 dta_bdy(ib_bdy)%ll_tem = .true. 287 dta_bdy(ib_bdy)%ll_sal = .true. 288 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for cn_tra' ) 218 289 END SELECT 219 IF( nn_tra(ib_bdy) .gt. 0) THEN290 IF( cn_tra(ib_bdy) /= 'none' ) THEN 220 291 SELECT CASE( nn_tra_dta(ib_bdy) ) ! 221 292 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' initial state used for bdy data' … … 226 297 227 298 IF ( ln_tra_dmp(ib_bdy) ) THEN 228 IF ( nn_tra(ib_bdy).EQ.0) THEN299 IF ( cn_tra(ib_bdy) == 'none' ) THEN 229 300 IF(lwp) WRITE(numout,*) 'No open boundary condition for tracers: ln_tra_dmp is set to .false.' 230 301 ln_tra_dmp(ib_bdy)=.false. 231 ELSEIF ( nn_tra(ib_bdy).EQ.1) THEN302 ELSEIF ( cn_tra(ib_bdy) == 'frs' ) THEN 232 303 CALL ctl_stop( 'Use FRS OR relaxation' ) 233 304 ELSE 234 305 IF(lwp) WRITE(numout,*) ' + T/S relaxation zone' 235 306 IF(lwp) WRITE(numout,*) ' Damping time scale: ',rn_time_dmp(ib_bdy),' days' 307 IF(lwp) WRITE(numout,*) ' Outflow damping time scale: ',rn_time_dmp_out(ib_bdy),' days' 236 308 IF((lwp).AND.rn_time_dmp(ib_bdy)<0) CALL ctl_stop( 'Time scale must be positive' ) 309 dta_bdy(ib_bdy)%ll_tem = .true. 310 dta_bdy(ib_bdy)%ll_sal = .true. 237 311 ENDIF 238 312 ELSE … … 243 317 #if defined key_lim2 244 318 IF(lwp) WRITE(numout,*) 'Boundary conditions for sea ice: ' 245 SELECT CASE( nn_ice_lim2(ib_bdy) ) 246 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' no open boundary condition' 247 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 248 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_tra' ) 319 SELECT CASE( cn_ice_lim(ib_bdy) ) 320 CASE('none') 321 IF(lwp) WRITE(numout,*) ' no open boundary condition' 322 dta_bdy(ib_bdy)%ll_frld = .false. 323 dta_bdy(ib_bdy)%ll_hicif = .false. 324 dta_bdy(ib_bdy)%ll_hsnif = .false. 325 CASE('frs') 326 IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 327 dta_bdy(ib_bdy)%ll_frld = .true. 328 dta_bdy(ib_bdy)%ll_hicif = .true. 329 dta_bdy(ib_bdy)%ll_hsnif = .true. 330 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for cn_ice_lim' ) 249 331 END SELECT 250 IF( nn_ice_lim2(ib_bdy) .gt. 0) THEN251 SELECT CASE( nn_ice_lim 2_dta(ib_bdy) ) !332 IF( cn_ice_lim(ib_bdy) /= 'none' ) THEN 333 SELECT CASE( nn_ice_lim_dta(ib_bdy) ) ! 252 334 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' initial state used for bdy data' 253 335 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' boundary data taken from file' 254 CASE DEFAULT ; CALL ctl_stop( 'nn_ice_lim2_dta must be 0 or 1' ) 336 CASE DEFAULT ; CALL ctl_stop( 'nn_ice_lim_dta must be 0 or 1' ) 337 END SELECT 338 ENDIF 339 IF(lwp) WRITE(numout,*) 340 #elif defined key_lim3 341 IF(lwp) WRITE(numout,*) 'Boundary conditions for sea ice: ' 342 SELECT CASE( cn_ice_lim(ib_bdy) ) 343 CASE('none') 344 IF(lwp) WRITE(numout,*) ' no open boundary condition' 345 dta_bdy(ib_bdy)%ll_a_i = .false. 346 dta_bdy(ib_bdy)%ll_ht_i = .false. 347 dta_bdy(ib_bdy)%ll_ht_s = .false. 348 CASE('frs') 349 IF(lwp) WRITE(numout,*) ' Flow Relaxation Scheme' 350 dta_bdy(ib_bdy)%ll_a_i = .true. 351 dta_bdy(ib_bdy)%ll_ht_i = .true. 352 dta_bdy(ib_bdy)%ll_ht_s = .true. 353 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for cn_ice_lim' ) 354 END SELECT 355 IF( cn_ice_lim(ib_bdy) /= 'none' ) THEN 356 SELECT CASE( nn_ice_lim_dta(ib_bdy) ) ! 357 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' initial state used for bdy data' 358 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' boundary data taken from file' 359 CASE DEFAULT ; CALL ctl_stop( 'nn_ice_lim_dta must be 0 or 1' ) 255 360 END SELECT 256 361 ENDIF … … 740 845 IF(lwp) THEN ! Since all procs read global data only need to do this check on one proc... 741 846 IF( nbrdta(ib,igrd,ib_bdy) < nbrdta(ibm1,igrd,ib_bdy) ) THEN 742 CALL ctl_stop('bdy_init : ERROR : boundary data in file & 743 must be defined in order of distance from edge nbr.', & 744 'A utility for re-ordering boundary coordinates and data & 745 files exists in the TOOLS/OBC directory') 847 CALL ctl_stop('bdy_init : ERROR : boundary data in file must be defined in order of distance from edge nbr.', & 848 'A utility for re-ordering boundary coordinates and data files exists in the TOOLS/OBC directory') 746 849 ENDIF 747 850 ENDIF … … 766 869 ALLOCATE( idx_bdy(ib_bdy)%nbr(ilen1,jpbgrd) ) 767 870 ALLOCATE( idx_bdy(ib_bdy)%nbd(ilen1,jpbgrd) ) 871 ALLOCATE( idx_bdy(ib_bdy)%nbdout(ilen1,jpbgrd) ) 768 872 ALLOCATE( idx_bdy(ib_bdy)%nbmap(ilen1,jpbgrd) ) 769 873 ALLOCATE( idx_bdy(ib_bdy)%nbw(ilen1,jpbgrd) ) 770 ALLOCATE( idx_bdy(ib_bdy)%flagu(ilen1 ) )771 ALLOCATE( idx_bdy(ib_bdy)%flagv(ilen1 ) )874 ALLOCATE( idx_bdy(ib_bdy)%flagu(ilen1,jpbgrd) ) 875 ALLOCATE( idx_bdy(ib_bdy)%flagv(ilen1,jpbgrd) ) 772 876 773 877 ! Dispatch mapping indices and discrete distances on each processor … … 937 1041 ENDDO 938 1042 ENDDO 1043 939 1044 ! definition of the i- and j- direction local boundaries arrays 940 1045 ! used for sending the boudaries … … 990 1095 nbr => idx_bdy(ib_bdy)%nbr(ib,igrd) 991 1096 idx_bdy(ib_bdy)%nbd(ib,igrd) = 1. / ( rn_time_dmp(ib_bdy) * rday ) & 1097 & *(FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2. ! quadratic 1098 idx_bdy(ib_bdy)%nbdout(ib,igrd) = 1. / ( rn_time_dmp_out(ib_bdy) * rday ) & 992 1099 & *(FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2. ! quadratic 993 1100 END DO … … 1092 1199 ENDDO 1093 1200 1201 ! For the flagu/flagv calculation below we require a version of fmask without 1202 ! the land boundary condition (shlat) included: 1203 CALL wrk_alloc(jpi,jpj,zfmask) 1204 DO ij = 2, jpjm1 1205 DO ii = 2, jpim1 1206 zfmask(ii,ij) = tmask(ii,ij ,1) * tmask(ii+1,ij ,1) & 1207 & * tmask(ii,ij+1,1) * tmask(ii+1,ij+1,1) 1208 END DO 1209 END DO 1210 1094 1211 ! Lateral boundary conditions 1212 CALL lbc_lnk( zfmask , 'F', 1. ) 1095 1213 CALL lbc_lnk( fmask , 'F', 1. ) ; CALL lbc_lnk( bdytmask(:,:), 'T', 1. ) 1096 1214 CALL lbc_lnk( bdyumask(:,:), 'U', 1. ) ; CALL lbc_lnk( bdyvmask(:,:), 'V', 1. ) … … 1098 1216 DO ib_bdy = 1, nb_bdy ! Indices and directions of rim velocity components 1099 1217 1100 idx_bdy(ib_bdy)%flagu(: ) = 0.e01101 idx_bdy(ib_bdy)%flagv(: ) = 0.e01218 idx_bdy(ib_bdy)%flagu(:,:) = 0.e0 1219 idx_bdy(ib_bdy)%flagv(:,:) = 0.e0 1102 1220 icount = 0 1103 1221 1104 !flagu = -1 : u component is normal to the dynamical boundary but its direction is outward 1105 !flagu = 0 : u is tangential 1106 !flagu = 1 : u is normal to the boundary and is direction is inward 1222 ! Calculate relationship of U direction to the local orientation of the boundary 1223 ! flagu = -1 : u component is normal to the dynamical boundary and its direction is outward 1224 ! flagu = 0 : u is tangential 1225 ! flagu = 1 : u is normal to the boundary and is direction is inward 1107 1226 1108 igrd = 2 ! u-component 1109 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1110 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 1111 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 1112 zefl = bdytmask(nbi ,nbj) 1113 zwfl = bdytmask(nbi+1,nbj) 1114 IF( zefl + zwfl == 2 ) THEN 1115 icount = icount + 1 1116 ELSE 1117 idx_bdy(ib_bdy)%flagu(ib)=-zefl+zwfl 1118 ENDIF 1227 DO igrd = 1,jpbgrd 1228 SELECT CASE( igrd ) 1229 CASE( 1 ) 1230 pmask => umask(:,:,1) 1231 i_offset = 0 1232 CASE( 2 ) 1233 pmask => bdytmask 1234 i_offset = 1 1235 CASE( 3 ) 1236 pmask => zfmask(:,:) 1237 i_offset = 0 1238 END SELECT 1239 icount = 0 1240 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1241 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 1242 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 1243 zefl = pmask(nbi+i_offset-1,nbj) 1244 zwfl = pmask(nbi+i_offset,nbj) 1245 ! This error check only works if you are using the bdyXmask arrays 1246 IF( i_offset == 1 .and. zefl + zwfl == 2 ) THEN 1247 icount = icount + 1 1248 IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(nbi),mjg(nbj) 1249 ELSE 1250 idx_bdy(ib_bdy)%flagu(ib,igrd) = -zefl + zwfl 1251 ENDIF 1252 END DO 1253 IF( icount /= 0 ) THEN 1254 IF(lwp) WRITE(numout,*) 1255 IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid(igrd),' grid points,', & 1256 ' are not boundary points (flagu calculation). Check nbi, nbj, indices for boundary set ',ib_bdy 1257 IF(lwp) WRITE(numout,*) ' ========== ' 1258 IF(lwp) WRITE(numout,*) 1259 nstop = nstop + 1 1260 ENDIF 1119 1261 END DO 1120 1262 1121 !flagv = -1 : u component is normal to the dynamical boundary but its direction is outward 1122 !flagv = 0 : u is tangential 1123 !flagv = 1 : u is normal to the boundary and is direction is inward 1124 1125 igrd = 3 ! v-component 1126 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1127 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 1128 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 1129 znfl = bdytmask(nbi,nbj ) 1130 zsfl = bdytmask(nbi,nbj+1) 1131 IF( znfl + zsfl == 2 ) THEN 1132 icount = icount + 1 1133 ELSE 1134 idx_bdy(ib_bdy)%flagv(ib) = -znfl + zsfl 1135 END IF 1263 ! Calculate relationship of V direction to the local orientation of the boundary 1264 ! flagv = -1 : v component is normal to the dynamical boundary but its direction is outward 1265 ! flagv = 0 : v is tangential 1266 ! flagv = 1 : v is normal to the boundary and is direction is inward 1267 1268 DO igrd = 1,jpbgrd 1269 SELECT CASE( igrd ) 1270 CASE( 1 ) 1271 pmask => vmask(:,:,1) 1272 j_offset = 0 1273 CASE( 2 ) 1274 pmask => zfmask(:,:) 1275 j_offset = 0 1276 CASE( 3 ) 1277 pmask => bdytmask 1278 j_offset = 1 1279 END SELECT 1280 icount = 0 1281 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1282 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 1283 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 1284 znfl = pmask(nbi,nbj+j_offset-1 ) 1285 zsfl = pmask(nbi,nbj+j_offset) 1286 ! This error check only works if you are using the bdyXmask arrays 1287 IF( j_offset == 1 .and. znfl + zsfl == 2 ) THEN 1288 IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(nbi),mjg(nbj) 1289 icount = icount + 1 1290 ELSE 1291 idx_bdy(ib_bdy)%flagv(ib,igrd) = -znfl + zsfl 1292 END IF 1293 END DO 1294 IF( icount /= 0 ) THEN 1295 IF(lwp) WRITE(numout,*) 1296 IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid(igrd),' grid points,', & 1297 ' are not boundary points (flagv calculation). Check nbi, nbj, indices for boundary set ',ib_bdy 1298 IF(lwp) WRITE(numout,*) ' ========== ' 1299 IF(lwp) WRITE(numout,*) 1300 nstop = nstop + 1 1301 ENDIF 1136 1302 END DO 1137 1303 1138 IF( icount /= 0 ) THEN 1139 IF(lwp) WRITE(numout,*) 1140 IF(lwp) WRITE(numout,*) ' E R R O R : Some data velocity points,', & 1141 ' are not boundary points. Check nbi, nbj, indices for boundary set ',ib_bdy 1142 IF(lwp) WRITE(numout,*) ' ========== ' 1143 IF(lwp) WRITE(numout,*) 1144 nstop = nstop + 1 1145 ENDIF 1146 1147 ENDDO 1304 END DO 1148 1305 1149 1306 ! Compute total lateral surface for volume correction: … … 1157 1314 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 1158 1315 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 1159 flagu => idx_bdy(ib_bdy)%flagu(ib )1316 flagu => idx_bdy(ib_bdy)%flagu(ib,igrd) 1160 1317 bdysurftot = bdysurftot + hu (nbi , nbj) & 1161 1318 & * e2u (nbi , nbj) * ABS( flagu ) & … … 1170 1327 nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 1171 1328 nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 1172 flagv => idx_bdy(ib_bdy)%flagv(ib )1329 flagv => idx_bdy(ib_bdy)%flagv(ib,igrd) 1173 1330 bdysurftot = bdysurftot + hv (nbi, nbj ) & 1174 1331 & * e1v (nbi, nbj ) * ABS( flagv ) & … … 1186 1343 DEALLOCATE(nbidta, nbjdta, nbrdta) 1187 1344 ENDIF 1345 1346 CALL wrk_dealloc(jpi,jpj,zfmask) 1188 1347 1189 1348 IF( nn_timing == 1 ) CALL timing_stop('bdy_init') … … 1580 1739 itest = 0 1581 1740 1582 IF ( nn_dyn2d(ib1)/=nn_dyn2d(ib2)) itest = itest + 11583 IF ( nn_dyn3d(ib1)/=nn_dyn3d(ib2)) itest = itest + 11584 IF ( nn_tra(ib1)/=nn_tra(ib2)) itest = itest + 11741 IF (cn_dyn2d(ib1)/=cn_dyn2d(ib2)) itest = itest + 1 1742 IF (cn_dyn3d(ib1)/=cn_dyn3d(ib2)) itest = itest + 1 1743 IF (cn_tra(ib1)/=cn_tra(ib2)) itest = itest + 1 1585 1744 ! 1586 1745 IF (nn_dyn2d_dta(ib1)/=nn_dyn2d_dta(ib2)) itest = itest + 1
Note: See TracChangeset
for help on using the changeset viewer.
