New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 3991 for branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90 – NEMO

Ignore:
Timestamp:
2013-07-29T11:04:44+02:00 (11 years ago)
Author:
davestorkey
Message:

New branch from later branch point on trunk so you can do a clean
diff of all the changes. Copy in changes from dev_r3891_METO1_MERCATOR6_OBC_BDY_merge.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_r3987_METO1_MERCATOR6_OBC_BDY_merge/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90

    r3909 r3991  
    7979      INTEGER  ::   jpbdtau, jpbdtas                       !   -       - 
    8080      INTEGER  ::   ib_bdy1, ib_bdy2, ib1, ib2             !   -       - 
     81      INTEGER  ::   i_offset, j_offset                     !   -       - 
    8182      INTEGER, POINTER  ::  nbi, nbj, nbr                  ! short cuts 
    82       REAL   , POINTER  ::  flagu, flagv                   !    -   - 
     83      REAL(wp), POINTER  ::  flagu, flagv                  !    -   - 
     84      REAL(wp), POINTER, DIMENSION(:,:)       ::   pmask    ! pointer to 2D mask fields 
    8385      REAL(wp) ::   zefl, zwfl, znfl, zsfl                 ! local scalars 
    8486      INTEGER, DIMENSION (2)                  ::   kdimsz 
     
    9294 
    9395      !! 
    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,             & 
     96      NAMELIST/nambdy/ nb_bdy, ln_coords_file, cn_coords_file,                 & 
     97         &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,    & 
     98         &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             &   
     99         &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, & 
    98100#if defined key_lim2 
    99          &             nn_ice_lim2, nn_ice_lim2_dta,                       & 
     101         &             cn_ice_lim2, nn_ice_lim2_dta,                           & 
    100102#endif 
    101103         &             ln_vol, nn_volctl, nn_rimwidth 
     
    128130      ln_mask_file      = .false. 
    129131      cn_mask_file(:)   = '' 
    130       nn_dyn2d(:)       = 0 
     132      cn_dyn2d(:)       = '' 
    131133      nn_dyn2d_dta(:)   = -1  ! uninitialised flag 
    132       nn_dyn3d(:)       = 0 
     134      cn_dyn3d(:)       = '' 
    133135      nn_dyn3d_dta(:)   = -1  ! uninitialised flag 
    134       nn_tra(:)         = 0 
     136      cn_tra(:)         = '' 
    135137      nn_tra_dta(:)     = -1  ! uninitialised flag 
    136138      ln_tra_dmp(:)     = .false. 
     
    138140      rn_time_dmp(:)    = 1. 
    139141#if defined key_lim2 
    140       nn_ice_lim2(:)    = 0 
     142      cn_ice_lim2(:)    = '' 
    141143      nn_ice_lim2_dta(:)= -1  ! uninitialised flag 
    142144#endif 
     
    172174 
    173175        IF(lwp) WRITE(numout,*) 'Boundary conditions for barotropic solution:  ' 
    174         SELECT CASE( nn_dyn2d(ib_bdy) )                   
    175           CASE(jp_none)         ;   IF(lwp) WRITE(numout,*) '      no open boundary condition'         
    176           CASE(jp_frs)          ;   IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
    177           CASE(jp_flather)      ;   IF(lwp) WRITE(numout,*) '      Flather radiation condition' 
    178           CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for nn_dyn2d' ) 
     176        SELECT CASE( cn_dyn2d(ib_bdy) )                   
     177          CASE('none')          
     178             IF(lwp) WRITE(numout,*) '      no open boundary condition'         
     179             dta_bdy(ib_bdy)%ll_ssh = .false. 
     180             dta_bdy(ib_bdy)%ll_u2d = .false. 
     181             dta_bdy(ib_bdy)%ll_v2d = .false. 
     182          CASE('frs')           
     183             IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
     184             dta_bdy(ib_bdy)%ll_ssh = .false. 
     185             dta_bdy(ib_bdy)%ll_u2d = .true. 
     186             dta_bdy(ib_bdy)%ll_v2d = .true. 
     187          CASE('flather')       
     188             IF(lwp) WRITE(numout,*) '      Flather radiation condition' 
     189             dta_bdy(ib_bdy)%ll_ssh = .true. 
     190             dta_bdy(ib_bdy)%ll_u2d = .true. 
     191             dta_bdy(ib_bdy)%ll_v2d = .true. 
     192          CASE('orlanski')      
     193             IF(lwp) WRITE(numout,*) '      Orlanski (fully oblique) radiation condition with adaptive nudging' 
     194             dta_bdy(ib_bdy)%ll_ssh = .false. 
     195             dta_bdy(ib_bdy)%ll_u2d = .true. 
     196             dta_bdy(ib_bdy)%ll_v2d = .true. 
     197          CASE('orlanski_npo')  
     198             IF(lwp) WRITE(numout,*) '      Orlanski (NPO) radiation condition with adaptive nudging' 
     199             dta_bdy(ib_bdy)%ll_ssh = .false. 
     200             dta_bdy(ib_bdy)%ll_u2d = .true. 
     201             dta_bdy(ib_bdy)%ll_v2d = .true. 
     202          CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for cn_dyn2d' ) 
    179203        END SELECT 
    180         IF( nn_dyn2d(ib_bdy) .gt. 0 ) THEN 
     204        IF( cn_dyn2d(ib_bdy) /= 'none' ) THEN 
    181205           SELECT CASE( nn_dyn2d_dta(ib_bdy) )                   !  
    182206              CASE( 0 )      ;   IF(lwp) WRITE(numout,*) '      initial state used for bdy data'         
     
    193217 
    194218        IF(lwp) WRITE(numout,*) 'Boundary conditions for baroclinic velocities:  ' 
    195         SELECT CASE( nn_dyn3d(ib_bdy) )                   
    196           CASE(jp_none)  ;   IF(lwp) WRITE(numout,*) '      no open boundary condition'         
    197           CASE(jp_frs)   ;   IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
    198           CASE( 2 )      ;   IF(lwp) WRITE(numout,*) '      Specified value' 
    199           CASE( 3 )      ;   IF(lwp) WRITE(numout,*) '      Zero baroclinic velocities (runoff case)' 
    200           CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for nn_dyn3d' ) 
     219        SELECT CASE( cn_dyn3d(ib_bdy) )                   
     220          CASE('none') 
     221             IF(lwp) WRITE(numout,*) '      no open boundary condition'         
     222             dta_bdy(ib_bdy)%ll_u3d = .false. 
     223             dta_bdy(ib_bdy)%ll_v3d = .false. 
     224          CASE('frs')        
     225             IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
     226             dta_bdy(ib_bdy)%ll_u3d = .true. 
     227             dta_bdy(ib_bdy)%ll_v3d = .true. 
     228          CASE('specified') 
     229             IF(lwp) WRITE(numout,*) '      Specified value' 
     230             dta_bdy(ib_bdy)%ll_u3d = .true. 
     231             dta_bdy(ib_bdy)%ll_v3d = .true. 
     232          CASE('zero') 
     233             IF(lwp) WRITE(numout,*) '      Zero baroclinic velocities (runoff case)' 
     234             dta_bdy(ib_bdy)%ll_u3d = .false. 
     235             dta_bdy(ib_bdy)%ll_v3d = .false. 
     236          CASE('orlanski') 
     237             IF(lwp) WRITE(numout,*) '      Orlanski (fully oblique) radiation condition with adaptive nudging' 
     238             dta_bdy(ib_bdy)%ll_u3d = .true. 
     239             dta_bdy(ib_bdy)%ll_v3d = .true. 
     240          CASE('orlanski_npo') 
     241             IF(lwp) WRITE(numout,*) '      Orlanski (NPO) radiation condition with adaptive nudging' 
     242             dta_bdy(ib_bdy)%ll_u3d = .true. 
     243             dta_bdy(ib_bdy)%ll_v3d = .true. 
     244          CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for cn_dyn3d' ) 
    201245        END SELECT 
    202         IF( nn_dyn3d(ib_bdy) .gt. 0 ) THEN 
     246        IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 
    203247           SELECT CASE( nn_dyn3d_dta(ib_bdy) )                   !  
    204248              CASE( 0 )      ;   IF(lwp) WRITE(numout,*) '      initial state used for bdy data'         
     
    209253 
    210254        IF ( ln_dyn3d_dmp(ib_bdy) ) THEN 
    211            IF ( nn_dyn3d(ib_bdy).EQ.0 ) THEN 
     255           IF ( cn_dyn3d(ib_bdy) == 'none' ) THEN 
    212256              IF(lwp) WRITE(numout,*) 'No open boundary condition for baroclinic velocities: ln_dyn3d_dmp is set to .false.' 
    213257              ln_dyn3d_dmp(ib_bdy)=.false. 
    214            ELSEIF ( nn_dyn3d(ib_bdy).EQ.1 ) THEN 
     258           ELSEIF ( cn_dyn3d(ib_bdy) == 'frs' ) THEN 
    215259              CALL ctl_stop( 'Use FRS OR relaxation' ) 
    216260           ELSE 
     
    218262              IF(lwp) WRITE(numout,*) '      Damping time scale: ',rn_time_dmp(ib_bdy),' days' 
    219263              IF((lwp).AND.rn_time_dmp(ib_bdy)<0) CALL ctl_stop( 'Time scale must be positive' ) 
     264              dta_bdy(ib_bdy)%ll_u3d = .true. 
     265              dta_bdy(ib_bdy)%ll_v3d = .true. 
    220266           ENDIF 
    221267        ELSE 
     
    225271 
    226272        IF(lwp) WRITE(numout,*) 'Boundary conditions for temperature and salinity:  ' 
    227         SELECT CASE( nn_tra(ib_bdy) )                   
    228           CASE(jp_none)  ;   IF(lwp) WRITE(numout,*) '      no open boundary condition'         
    229           CASE(jp_frs)   ;   IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
    230           CASE( 2 )      ;   IF(lwp) WRITE(numout,*) '      Specified value' 
    231           CASE( 3 )      ;   IF(lwp) WRITE(numout,*) '      Neumann conditions' 
    232           CASE( 4 )      ;   IF(lwp) WRITE(numout,*) '      Runoff conditions : Neumann for T and specified to 0.1 for salinity' 
    233           CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for nn_tra' ) 
     273        SELECT CASE( cn_tra(ib_bdy) )                   
     274          CASE('none') 
     275             IF(lwp) WRITE(numout,*) '      no open boundary condition'         
     276             dta_bdy(ib_bdy)%ll_tem = .false. 
     277             dta_bdy(ib_bdy)%ll_sal = .false. 
     278          CASE('frs') 
     279             IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
     280             dta_bdy(ib_bdy)%ll_tem = .true. 
     281             dta_bdy(ib_bdy)%ll_sal = .true. 
     282          CASE('specified') 
     283             IF(lwp) WRITE(numout,*) '      Specified value' 
     284             dta_bdy(ib_bdy)%ll_tem = .true. 
     285             dta_bdy(ib_bdy)%ll_sal = .true. 
     286          CASE('neumann') 
     287             IF(lwp) WRITE(numout,*) '      Neumann conditions' 
     288             dta_bdy(ib_bdy)%ll_tem = .false. 
     289             dta_bdy(ib_bdy)%ll_sal = .false. 
     290          CASE('runoff') 
     291             IF(lwp) WRITE(numout,*) '      Runoff conditions : Neumann for T and specified to 0.1 for salinity' 
     292             dta_bdy(ib_bdy)%ll_tem = .false. 
     293             dta_bdy(ib_bdy)%ll_sal = .false. 
     294          CASE('orlanski') 
     295             IF(lwp) WRITE(numout,*) '      Orlanski (fully oblique) radiation condition with adaptive nudging' 
     296             dta_bdy(ib_bdy)%ll_tem = .true. 
     297             dta_bdy(ib_bdy)%ll_sal = .true. 
     298          CASE('orlanski_npo') 
     299             IF(lwp) WRITE(numout,*) '      Orlanski (NPO) radiation condition with adaptive nudging' 
     300             dta_bdy(ib_bdy)%ll_tem = .true. 
     301             dta_bdy(ib_bdy)%ll_sal = .true. 
     302          CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for cn_tra' ) 
    234303        END SELECT 
    235         IF( nn_tra(ib_bdy) .gt. 0 ) THEN 
     304        IF( cn_tra(ib_bdy) /= 'none' ) THEN 
    236305           SELECT CASE( nn_tra_dta(ib_bdy) )                   !  
    237306              CASE( 0 )      ;   IF(lwp) WRITE(numout,*) '      initial state used for bdy data'         
     
    242311 
    243312        IF ( ln_tra_dmp(ib_bdy) ) THEN 
    244            IF ( nn_tra(ib_bdy).EQ.0 ) THEN 
     313           IF ( cn_tra(ib_bdy) == 'none' ) THEN 
    245314              IF(lwp) WRITE(numout,*) 'No open boundary condition for tracers: ln_tra_dmp is set to .false.' 
    246315              ln_tra_dmp(ib_bdy)=.false. 
    247            ELSEIF ( nn_tra(ib_bdy).EQ.1 ) THEN 
     316           ELSEIF ( cn_tra(ib_bdy) == 'frs' ) THEN 
    248317              CALL ctl_stop( 'Use FRS OR relaxation' ) 
    249318           ELSE 
    250319              IF(lwp) WRITE(numout,*) '      + T/S relaxation zone' 
    251320              IF(lwp) WRITE(numout,*) '      Damping time scale: ',rn_time_dmp(ib_bdy),' days' 
     321              IF(lwp) WRITE(numout,*) '      Outflow damping time scale: ',rn_time_dmp_out(ib_bdy),' days' 
    252322              IF((lwp).AND.rn_time_dmp(ib_bdy)<0) CALL ctl_stop( 'Time scale must be positive' ) 
     323              dta_bdy(ib_bdy)%ll_tem = .true. 
     324              dta_bdy(ib_bdy)%ll_sal = .true. 
    253325           ENDIF 
    254326        ELSE 
     
    259331#if defined key_lim2 
    260332        IF(lwp) WRITE(numout,*) 'Boundary conditions for sea ice:  ' 
    261         SELECT CASE( nn_ice_lim2(ib_bdy) )                   
    262           CASE( 0 )      ;   IF(lwp) WRITE(numout,*) '      no open boundary condition'         
    263           CASE( 1 )      ;   IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
    264           CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for nn_tra' ) 
     333        SELECT CASE( cn_ice_lim2(ib_bdy) )                   
     334          CASE('none') 
     335             IF(lwp) WRITE(numout,*) '      no open boundary condition'         
     336             dta_bdy(ib_bdy)%ll_frld  = .false. 
     337             dta_bdy(ib_bdy)%ll_hicif = .false. 
     338             dta_bdy(ib_bdy)%ll_hsnif = .false. 
     339          CASE('frs') 
     340             IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme' 
     341             dta_bdy(ib_bdy)%ll_frld  = .true. 
     342             dta_bdy(ib_bdy)%ll_hicif = .true. 
     343             dta_bdy(ib_bdy)%ll_hsnif = .true. 
     344          CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for cn_ice_lim2' ) 
    265345        END SELECT 
    266         IF( nn_ice_lim2(ib_bdy) .gt. 0 ) THEN  
     346        IF( cn_ice_lim2(ib_bdy) /= 'none' ) THEN  
    267347           SELECT CASE( nn_ice_lim2_dta(ib_bdy) )                   !  
    268348              CASE( 0 )      ;   IF(lwp) WRITE(numout,*) '      initial state used for bdy data'         
     
    744824               IF(lwp) THEN         ! Since all procs read global data only need to do this check on one proc... 
    745825                  IF( nbrdta(ib,igrd,ib_bdy) < nbrdta(ibm1,igrd,ib_bdy) ) THEN 
    746                      CALL ctl_stop('bdy_init : ERROR : boundary data in file  & 
    747                                     must be defined in order of distance from edge nbr.', & 
    748                                    'A utility for re-ordering boundary coordinates and data & 
    749                                     files exists in the TOOLS/OBC directory') 
     826                     CALL ctl_stop('bdy_init : ERROR : boundary data in file must be defined in order of distance from edge nbr.', & 
     827                                   'A utility for re-ordering boundary coordinates and data files exists in the TOOLS/OBC directory') 
    750828                  ENDIF     
    751829               ENDIF 
     
    770848         ALLOCATE( idx_bdy(ib_bdy)%nbr(ilen1,jpbgrd) ) 
    771849         ALLOCATE( idx_bdy(ib_bdy)%nbd(ilen1,jpbgrd) ) 
     850         ALLOCATE( idx_bdy(ib_bdy)%nbdout(ilen1,jpbgrd) ) 
    772851         ALLOCATE( idx_bdy(ib_bdy)%nbmap(ilen1,jpbgrd) ) 
    773852         ALLOCATE( idx_bdy(ib_bdy)%nbw(ilen1,jpbgrd) ) 
    774          ALLOCATE( idx_bdy(ib_bdy)%flagu(ilen1) ) 
    775          ALLOCATE( idx_bdy(ib_bdy)%flagv(ilen1) ) 
     853         ALLOCATE( idx_bdy(ib_bdy)%flagu(ilen1,jpbgrd) ) 
     854         ALLOCATE( idx_bdy(ib_bdy)%flagv(ilen1,jpbgrd) ) 
    776855 
    777856         ! Dispatch mapping indices and discrete distances on each processor 
     
    9411020            ENDDO 
    9421021         ENDDO  
     1022 
    9431023         ! definition of the i- and j- direction local boundaries arrays 
    9441024         ! used for sending the boudaries 
     
    9941074               nbr => idx_bdy(ib_bdy)%nbr(ib,igrd) 
    9951075               idx_bdy(ib_bdy)%nbd(ib,igrd) = 1. / ( rn_time_dmp(ib_bdy) * rday ) &  
     1076               & *(FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2.   ! quadratic 
     1077               idx_bdy(ib_bdy)%nbdout(ib,igrd) = 1. / ( rn_time_dmp_out(ib_bdy) * rday ) &  
    9961078               & *(FLOAT(nn_rimwidth(ib_bdy)+1-nbr)/FLOAT(nn_rimwidth(ib_bdy)))**2.   ! quadratic 
    9971079            END DO 
     
    10961178      ENDDO 
    10971179 
     1180      ! bdyfmask required for flagu, flagv calculations below even though F-points  
     1181      ! not defined for BDY grid.  
     1182      DO ij = 2, jpjm1 
     1183         DO ii = 2, jpim1 
     1184            bdyfmask(ii,ij) = bdytmask(ii,ij  ) * bdytmask(ii+1,ij  )   & 
     1185           &                * bdytmask(ii,ij+1) * bdytmask(ii+1,ij+1) 
     1186         END DO       
     1187      END DO 
     1188 
    10981189      ! Lateral boundary conditions 
    10991190      CALL lbc_lnk( fmask        , 'F', 1. )   ;   CALL lbc_lnk( bdytmask(:,:), 'T', 1. ) 
     
    11021193      DO ib_bdy = 1, nb_bdy       ! Indices and directions of rim velocity components 
    11031194 
    1104          idx_bdy(ib_bdy)%flagu(:) = 0.e0 
    1105          idx_bdy(ib_bdy)%flagv(:) = 0.e0 
     1195         idx_bdy(ib_bdy)%flagu(:,:) = 0.e0 
     1196         idx_bdy(ib_bdy)%flagv(:,:) = 0.e0 
    11061197         icount = 0  
    11071198 
    1108          !flagu = -1 : u component is normal to the dynamical boundary but its direction is outward 
    1109          !flagu =  0 : u is tangential 
    1110          !flagu =  1 : u is normal to the boundary and is direction is inward 
     1199         ! Calculate relationship of U direction to the local orientation of the boundary 
     1200         ! flagu = -1 : u component is normal to the dynamical boundary and its direction is outward 
     1201         ! flagu =  0 : u is tangential 
     1202         ! flagu =  1 : u is normal to the boundary and is direction is inward 
    11111203   
    1112          igrd = 2      ! u-component  
    1113          DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
    1114             nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 
    1115             nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1116             zefl = bdytmask(nbi  ,nbj) 
    1117             zwfl = bdytmask(nbi+1,nbj) 
    1118             IF( zefl + zwfl == 2 ) THEN 
    1119                icount = icount + 1 
    1120             ELSE 
    1121                idx_bdy(ib_bdy)%flagu(ib)=-zefl+zwfl 
    1122             ENDIF 
     1204         DO igrd = 1,jpbgrd  
     1205            SELECT CASE( igrd ) 
     1206               CASE( 1 ) 
     1207                  cgrid = 'T'  
     1208                  pmask => umask(:,:,1) 
     1209                  i_offset = 0 
     1210               CASE( 2 )  
     1211                  cgrid = 'U' 
     1212                  pmask => bdytmask 
     1213                  i_offset = 1 
     1214               CASE( 3 )  
     1215                  cgrid = 'V' 
     1216                  pmask => fmask(:,:,1) 
     1217                  i_offset = 0 
     1218            END SELECT  
     1219            icount = 0 
     1220            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
     1221               nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 
     1222               nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 
     1223               zefl = pmask(nbi+i_offset-1,nbj) 
     1224               zwfl = pmask(nbi+i_offset,nbj) 
     1225               IF( zefl + zwfl == 2 * i_offset ) THEN 
     1226                  icount = icount + 1 
     1227                  IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(nbi),mjg(nbj) 
     1228               ELSE 
     1229                  idx_bdy(ib_bdy)%flagu(ib,igrd)=-zefl+zwfl 
     1230               ENDIF 
     1231            END DO 
     1232            IF( icount /= 0 ) THEN 
     1233               IF(lwp) WRITE(numout,*) 
     1234               IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid,' grid points,',   & 
     1235                  ' are not boundary points (flagu calculation). Check nbi, nbj, indices for boundary set ',ib_bdy 
     1236               IF(lwp) WRITE(numout,*) ' ========== ' 
     1237               IF(lwp) WRITE(numout,*) 
     1238               nstop = nstop + 1 
     1239            ENDIF  
    11231240         END DO 
    11241241 
    1125          !flagv = -1 : u component is normal to the dynamical boundary but its direction is outward 
    1126          !flagv =  0 : u is tangential 
    1127          !flagv =  1 : u is normal to the boundary and is direction is inward 
    1128  
    1129          igrd = 3      ! v-component 
    1130          DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
    1131             nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 
    1132             nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1133             znfl = bdytmask(nbi,nbj  ) 
    1134             zsfl = bdytmask(nbi,nbj+1) 
    1135             IF( znfl + zsfl == 2 ) THEN 
    1136                icount = icount + 1 
    1137             ELSE 
    1138                idx_bdy(ib_bdy)%flagv(ib) = -znfl + zsfl 
    1139             END IF 
     1242         ! Calculate relationship of V direction to the local orientation of the boundary 
     1243         ! flagv = -1 : u component is normal to the dynamical boundary but its direction is outward 
     1244         ! flagv =  0 : u is tangential 
     1245         ! flagv =  1 : u is normal to the boundary and is direction is inward 
     1246 
     1247         DO igrd = 1,jpbgrd  
     1248            SELECT CASE( igrd ) 
     1249               CASE( 1 ) 
     1250                  cgrid = 'T' 
     1251                  pmask => vmask(:,:,1) 
     1252                  j_offset = 0 
     1253               CASE( 2 ) 
     1254                  cgrid = 'U' 
     1255                  pmask => fmask(:,:,1) 
     1256                  j_offset = 0 
     1257               CASE( 3 ) 
     1258                  cgrid = 'V' 
     1259                  pmask => bdytmask 
     1260                  j_offset = 1 
     1261            END SELECT  
     1262            icount = 0 
     1263            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
     1264               nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 
     1265               nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 
     1266               znfl = pmask(nbi,nbj+j_offset-1  ) 
     1267               zsfl = pmask(nbi,nbj+j_offset) 
     1268               IF( znfl + zsfl == 2 * j_offset ) THEN 
     1269                  IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(nbi),mjg(nbj) 
     1270                  icount = icount + 1 
     1271               ELSE 
     1272                  idx_bdy(ib_bdy)%flagv(ib,igrd) = -znfl + zsfl 
     1273               END IF 
     1274            END DO 
     1275            IF( icount /= 0 ) THEN 
     1276               IF(lwp) WRITE(numout,*) 
     1277               IF(lwp) WRITE(numout,*) ' E R R O R : Some ',cgrid,' grid points,',   & 
     1278                  ' are not boundary points (flagv calculation). Check nbi, nbj, indices for boundary set ',ib_bdy 
     1279               IF(lwp) WRITE(numout,*) ' ========== ' 
     1280               IF(lwp) WRITE(numout,*) 
     1281               nstop = nstop + 1 
     1282            ENDIF  
    11401283         END DO 
    11411284 
    1142          IF( icount /= 0 ) THEN 
    1143             IF(lwp) WRITE(numout,*) 
    1144             IF(lwp) WRITE(numout,*) ' E R R O R : Some data velocity points,',   & 
    1145                ' are not boundary points. Check nbi, nbj, indices for boundary set ',ib_bdy 
    1146             IF(lwp) WRITE(numout,*) ' ========== ' 
    1147             IF(lwp) WRITE(numout,*) 
    1148             nstop = nstop + 1 
    1149          ENDIF  
    1150      
    1151       ENDDO 
     1285      END DO 
    11521286 
    11531287      ! Compute total lateral surface for volume correction: 
     
    11611295               nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 
    11621296               nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1163                flagu => idx_bdy(ib_bdy)%flagu(ib) 
     1297               flagu => idx_bdy(ib_bdy)%flagu(ib,igrd) 
    11641298               bdysurftot = bdysurftot + hu     (nbi  , nbj)                           & 
    11651299                  &                    * e2u    (nbi  , nbj) * ABS( flagu ) & 
     
    11741308               nbi => idx_bdy(ib_bdy)%nbi(ib,igrd) 
    11751309               nbj => idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1176                flagv => idx_bdy(ib_bdy)%flagv(ib) 
     1310               flagv => idx_bdy(ib_bdy)%flagv(ib,igrd) 
    11771311               bdysurftot = bdysurftot + hv     (nbi, nbj  )                           & 
    11781312                  &                    * e1v    (nbi, nbj  ) * ABS( flagv ) & 
     
    15841718      itest = 0 
    15851719 
    1586       IF (nn_dyn2d(ib1)/=nn_dyn2d(ib2)) itest = itest + 1 
    1587       IF (nn_dyn3d(ib1)/=nn_dyn3d(ib2)) itest = itest + 1 
    1588       IF (nn_tra(ib1)/=nn_tra(ib2)) itest = itest + 1 
     1720      IF (cn_dyn2d(ib1)/=cn_dyn2d(ib2)) itest = itest + 1 
     1721      IF (cn_dyn3d(ib1)/=cn_dyn3d(ib2)) itest = itest + 1 
     1722      IF (cn_tra(ib1)/=cn_tra(ib2)) itest = itest + 1 
    15891723      ! 
    15901724      IF (nn_dyn2d_dta(ib1)/=nn_dyn2d_dta(ib2)) itest = itest + 1 
Note: See TracChangeset for help on using the changeset viewer.