*IF DEF,A11_1A                                                             TRBDRY1A.2      
C ******************************COPYRIGHT******************************    GTS2F400.10585  
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved.    GTS2F400.10586  
C                                                                          GTS2F400.10587  
C Use, duplication or disclosure of this code is subject to the            GTS2F400.10588  
C restrictions as set forth in the contract.                               GTS2F400.10589  
C                                                                          GTS2F400.10590  
C                Meteorological Office                                     GTS2F400.10591  
C                London Road                                               GTS2F400.10592  
C                BRACKNELL                                                 GTS2F400.10593  
C                Berkshire UK                                              GTS2F400.10594  
C                RG12 2SZ                                                  GTS2F400.10595  
C                                                                          GTS2F400.10596  
C If no contract has been raised with this copy of the code, the use,      GTS2F400.10597  
C duplication or disclosure of it is strictly prohibited.  Permission      GTS2F400.10598  
C to do so must first be obtained in writing from the Head of Numerical    GTS2F400.10599  
C Modelling at the above address.                                          GTS2F400.10600  
C ******************************COPYRIGHT******************************    GTS2F400.10601  
C                                                                          GTS2F400.10602  
C*LL  SUBROUTINE TRBDRY-------------------------------------------------   TRBDRY1A.3      
CLL                                                                        TRBDRY1A.4      
CLL  Purpose: Special routine to add psuedo source terms to boundary       TRBDRY1A.5      
CLL           data in limited area model.                                  TRBDRY1A.6      
CLL  Method:  runs along each boundary setting boundary concentration      TRBDRY1A.7      
CLL           if there is inflow. The boundary concentration is            TRBDRY1A.8      
CLL           set using a call to function BDRYV. The function BDRYV       TRBDRY1A.9      
CLL           is specific to a model configuration: the current            TRBDRY1A.10     
CLL           version (3.4) is specific to UK MES.                         TRBDRY1A.11     
CLL           On outflow boundaries the concentration is set equal         TRBDRY1A.12     
CLL           to the nearest gridpoint inside the boundary.                TRBDRY1A.13     
CLL                                                                        TRBDRY1A.14     
CLL                                                                        TRBDRY1A.15     
CLL Pete Clark  <- programmer of some or all of previous code or changes   TRBDRY1A.16     
CLL                                                                        TRBDRY1A.17     
CLL  Model            Modification history from model version 3.4:         TRBDRY1A.18     
CLL version  Date                                                          TRBDRY1A.19     
CLL  4.2  15/08/96  Add MPP code. Remove unused variables.  RTHBarnes.     ARB1F402.727    
CLL                                                                        TRBDRY1A.20     
CLL  Programming standard: Unified Model Documentation Paper No 3,         TRBDRY1A.21     
CLL                        Version 7, dated 11/3/93.                       TRBDRY1A.22     
CLL                                                                        TRBDRY1A.23     
CLL                                                                        TRBDRY1A.24     
C*L  Arguments:---------------------------------------------------------   TRBDRY1A.25     

      SUBROUTINE TRBDRY(                                                    1,4TRBDRY1A.26     
     & AK,BK,                                                              TRBDRY1A.27     
     & POINTS,PFIELD,UFIELD,ROW_LENGTH,                                    TRBDRY1A.30     
*CALL ARGFLDPT                                                             ARB1F402.728    
     & PSTAR,                                                              TRBDRY1A.31     
     & U,V,                                                                TRBDRY1A.32     
     & TR,TIMESTEP,ERROR                                                   TRBDRY1A.33     
     &)                                                                    TRBDRY1A.34     
      IMPLICIT NONE                                                        TRBDRY1A.35     
      INTEGER                                                              TRBDRY1A.36     
     & POINTS              ! IN No. of gridpoints being processed.         TRBDRY1A.37     
     &,PFIELD              ! IN No. of points in global field (at one      TRBDRY1A.38     
C                          !    vertical level).                           TRBDRY1A.39     
     &,UFIELD              ! IN No. of u points in global field (at one    TRBDRY1A.40     
C                          !    vertical level).                           TRBDRY1A.41     
     &,ROW_LENGTH          ! IN Length of a row.                           TRBDRY1A.42     
! All TYPFLDPT arguments are intent IN                                     ARB1F402.729    
*CALL TYPFLDPT                                                             ARB1F402.730    
      REAL                                                                 TRBDRY1A.43     
     & AK,BK             ! IN Layer  ak and bk                             TRBDRY1A.44     
     &,PSTAR(PFIELD)     ! IN Surface pressure                             TRBDRY1A.47     
     &,U(UFIELD),V(UFIELD) ! IN U and V component of wind                  TRBDRY1A.48     
     &,TR(PFIELD)        ! INOUT Tracer field (kg per kg air).             ARB1F402.731    
     &,TIMESTEP          ! IN Timestep in seconds                          TRBDRY1A.51     
      INTEGER ERROR      ! OUT Error return code.                          TRBDRY1A.52     
*CALL C_PI                                                                 TRBDRY1A.53     
C                                                                          TRBDRY1A.55     
C*L  Workspace usage----------------------------------------------------   TRBDRY1A.56     
C*L  External subroutine called ----------------------------------------   TRBDRY1A.57     
C     None                                                                 TRBDRY1A.58     
C     EXTERNAL                                                             TRBDRY1A.59     
C* Local, including SAVE'd, storage------------------------------------    TRBDRY1A.60     
C                                                                          TRBDRY1A.61     
C  (a) Scalars effectively expanded to workspace by the Cray (using        TRBDRY1A.62     
C      vector registers).                                                  TRBDRY1A.63     
      REAL                                                                 TRBDRY1A.64     
     & DC,WDIR,WSPEED,PRESS                                                ARB1F402.732    
      REAL                                                                 TRBDRY1A.66     
     & BDRYV ! FUNCTION giving boundary value.                             TRBDRY1A.67     
C                                                                          TRBDRY1A.68     
C  (b) Others.                                                             TRBDRY1A.69     
      INTEGER I,IU   ! Loop counters                                       ARB1F402.733    
C-----------------------------------------------------------------------   TRBDRY1A.72     
C  Check input arguments for potential over-writing problems.              TRBDRY1A.73     
C-----------------------------------------------------------------------   TRBDRY1A.74     
      ERROR=0                                                              TRBDRY1A.75     
      IF(POINTS.GT.PFIELD)THEN                                             TRBDRY1A.76     
        ERROR=1                                                            TRBDRY1A.77     
        WRITE(6,*)'Error in TRBDRY: POINTS greater than PFIELD.'           TRBDRY1A.78     
        GOTO 9999                                                          TRBDRY1A.79     
      ENDIF                                                                TRBDRY1A.80     
C                                                                          TRBDRY1A.84     
*IF DEF,MPP                                                                ARB1F402.734    
      IF (at_top_of_LPG) THEN                                              ARB1F402.735    
*ENDIF                                                                     ARB1F402.736    
C-----------------------------------------------------------------------   TRBDRY1A.85     
CL Loop across top row.                                                    TRBDRY1A.86     
C-----------------------------------------------------------------------   TRBDRY1A.87     
C                                                                          TRBDRY1A.88     
      DO  I = TOP_ROW_START,TOP_ROW_START+ROW_LENGTH-1                     ARB1F402.737    
        IU=I                                                               TRBDRY1A.93     
        IF(V(IU) .LT. 0.0)THEN                                             TRBDRY1A.94     
          PRESS=AK+BK*PSTAR(I)                                             TRBDRY1A.95     
          WDIR=ATAN2(V(IU),U(IU))*RECIP_PI_OVER_180                        TRBDRY1A.96     
          WSPEED=SQRT(U(IU)*U(IU) + V(IU)*V(IU))                           TRBDRY1A.97     
          DC = BDRYV(WDIR,WSPEED,PRESS)                                    TRBDRY1A.98     
          TR(I) = DC                                                       TRBDRY1A.99     
        ELSE                                                               TRBDRY1A.100    
          TR(I) = TR(I+ROW_LENGTH)                                         TRBDRY1A.101    
        ENDIF                                                              TRBDRY1A.102    
      ENDDO ! Loop over points                                             TRBDRY1A.103    
C                                                                          TRBDRY1A.104    
*IF DEF,MPP                                                                ARB1F402.738    
      ENDIF                                                                ARB1F402.739    
*ENDIF                                                                     ARB1F402.740    
C                                                                          TRBDRY1A.105    
*IF DEF,MPP                                                                ARB1F402.741    
      IF (at_base_of_LPG) THEN                                             ARB1F402.742    
*ENDIF                                                                     ARB1F402.743    
C-----------------------------------------------------------------------   TRBDRY1A.106    
CL Loop across bottom row.                                                 TRBDRY1A.107    
C-----------------------------------------------------------------------   TRBDRY1A.108    
C                                                                          TRBDRY1A.109    
      DO  I = P_BOT_ROW_START,P_BOT_ROW_START+ROW_LENGTH-1                 ARB1F402.744    
        IU = I-ROW_LENGTH                                                  ARB1F402.745    
        IF(V(IU) .GT. 0.0)THEN                                             TRBDRY1A.115    
          PRESS=AK+BK*PSTAR(I)                                             TRBDRY1A.116    
          WDIR=ATAN2(V(IU),U(IU))*RECIP_PI_OVER_180                        TRBDRY1A.117    
          WSPEED=SQRT(U(IU)*U(IU) + V(IU)*V(IU))                           TRBDRY1A.118    
          DC = BDRYV(WDIR,WSPEED,PRESS)                                    TRBDRY1A.119    
          TR(I) = DC                                                       TRBDRY1A.120    
        ELSE                                                               TRBDRY1A.121    
          TR(I) = TR(I-ROW_LENGTH)                                         TRBDRY1A.122    
        ENDIF                                                              TRBDRY1A.123    
      ENDDO ! Loop over points                                             TRBDRY1A.124    
C                                                                          TRBDRY1A.125    
*IF DEF,MPP                                                                ARB1F402.746    
      ENDIF                                                                ARB1F402.747    
*ENDIF                                                                     ARB1F402.748    
C                                                                          ARB1F402.749    
*IF DEF,MPP                                                                ARB1F402.750    
      IF (at_left_of_LPG) THEN                                             ARB1F402.751    
*ENDIF                                                                     ARB1F402.752    
C-----------------------------------------------------------------------   TRBDRY1A.126    
CL Loop across left column                                                 TRBDRY1A.127    
C-----------------------------------------------------------------------   TRBDRY1A.128    
C                                                                          TRBDRY1A.129    
      DO  I = TOP_ROW_START+FIRST_ROW_PT-1,                                ARB1F402.753    
     &        P_BOT_ROW_START+FIRST_ROW_PT-1,ROW_LENGTH                    ARB1F402.754    
        IF (I .eq. P_BOT_ROW_START+FIRST_ROW_PT-1) THEN                    ARB1F402.755    
*IF DEF,MPP                                                                ARB1F402.756    
         IF (at_base_of_LPG) THEN                                          ARB1F402.757    
*ENDIF                                                                     ARB1F402.758    
          IU = I-ROW_LENGTH                                                ARB1F402.759    
*IF DEF,MPP                                                                ARB1F402.760    
         END IF                                                            ARB1F402.761    
*ENDIF                                                                     ARB1F402.762    
        ELSE                                                               TRBDRY1A.136    
          IU = I                                                           ARB1F402.763    
        END IF                                                             ARB1F402.764    
        IF(U(IU) .GT. 0.0)THEN                                             TRBDRY1A.139    
          PRESS=AK+BK*PSTAR(I)                                             TRBDRY1A.140    
          WDIR=ATAN2(V(IU),U(IU))*RECIP_PI_OVER_180                        TRBDRY1A.141    
          WSPEED=SQRT(U(IU)*U(IU) + V(IU)*V(IU))                           TRBDRY1A.142    
          DC = BDRYV(WDIR,WSPEED,PRESS)                                    TRBDRY1A.143    
          TR(I) = DC                                                       TRBDRY1A.144    
        ELSE                                                               TRBDRY1A.145    
          TR(I) = TR(I+1)                                                  TRBDRY1A.146    
        ENDIF                                                              TRBDRY1A.147    
      ENDDO ! Loop over points                                             TRBDRY1A.148    
C                                                                          TRBDRY1A.149    
*IF DEF,MPP                                                                ARB1F402.765    
      ENDIF                                                                ARB1F402.766    
*ENDIF                                                                     ARB1F402.767    
C                                                                          ARB1F402.768    
*IF DEF,MPP                                                                ARB1F402.769    
      IF (at_right_of_LPG) THEN                                            ARB1F402.770    
*ENDIF                                                                     ARB1F402.771    
C-----------------------------------------------------------------------   TRBDRY1A.150    
CL Loop across right column                                                TRBDRY1A.151    
C-----------------------------------------------------------------------   TRBDRY1A.152    
C                                                                          TRBDRY1A.153    
      DO  I = TOP_ROW_START+LAST_ROW_PT-1,                                 ARB1F402.772    
     &        P_BOT_ROW_START+LAST_ROW_PT-1,ROW_LENGTH                     ARB1F402.773    
        IF (I .eq. P_BOT_ROW_START+LAST_ROW_PT-1) THEN                     ARB1F402.774    
*IF DEF,MPP                                                                ARB1F402.775    
         IF (at_base_of_LPG) THEN                                          ARB1F402.776    
*ENDIF                                                                     ARB1F402.777    
          IU = I-ROW_LENGTH                                                ARB1F402.778    
*IF DEF,MPP                                                                ARB1F402.779    
         END IF                                                            ARB1F402.780    
*ENDIF                                                                     ARB1F402.781    
        ELSE                                                               TRBDRY1A.160    
          IU = I                                                           ARB1F402.782    
        END IF                                                             ARB1F402.783    
        IF(U(IU) .LT. 0.0)THEN                                             TRBDRY1A.163    
          PRESS=AK+BK*PSTAR(I)                                             TRBDRY1A.164    
          WDIR=ATAN2(V(IU),U(IU))*RECIP_PI_OVER_180                        TRBDRY1A.165    
          WSPEED=SQRT(U(IU)*U(IU) + V(IU)*V(IU))                           TRBDRY1A.166    
          DC = BDRYV(WDIR,WSPEED,PRESS)                                    TRBDRY1A.167    
          TR(I) = DC                                                       TRBDRY1A.168    
        ELSE                                                               TRBDRY1A.169    
          TR(I) = TR(I-1)                                                  TRBDRY1A.170    
        ENDIF                                                              TRBDRY1A.171    
      ENDDO ! Loop over points                                             TRBDRY1A.172    
C                                                                          ARB1F402.784    
*IF DEF,MPP                                                                ARB1F402.785    
      ENDIF                                                                ARB1F402.786    
*ENDIF                                                                     ARB1F402.787    
C                                                                          TRBDRY1A.173    
 9999 CONTINUE ! Error exit                                                TRBDRY1A.174    
      RETURN                                                               TRBDRY1A.175    
      END                                                                  TRBDRY1A.176    
C*LL  FUNCTION BDRYV----------------------------------------------------   TRBDRY1A.177    
CLL                                                                        TRBDRY1A.178    
CLL  Purpose: Special routine to add psuedo source terms to boundary       TRBDRY1A.179    
CLL           data in limited area.                                        TRBDRY1A.180    
CLL  PARAMETERS ARE SPECIFIC TO UK MESOSCALE MODEL                         TRBDRY1A.181    
CLL  Method:  The boundary concentrations are computed using a             TRBDRY1A.182    
CLL           simple model of transport from sources outside the           TRBDRY1A.183    
CLL           model. Analysis of the source distribution outside           TRBDRY1A.184    
CLL           the UK MES shows that it can be well represented by          TRBDRY1A.185    
CLL           a line source at constant radius from the centre of          TRBDRY1A.186    
CLL           the model, with a source distribution given by the           TRBDRY1A.187    
CLL           sum of two Gaussians. Concentrations from these are          TRBDRY1A.188    
CLL           computed assuming transport using the local windspeed u      TRBDRY1A.189    
CLL           or 1 m/s, whichever is stronger, over a distance             TRBDRY1A.190    
CLL           determined from the centroid of the source distribution, x   TRBDRY1A.191    
CLL           with a linear transformation rate k from emission to         TRBDRY1A.192    
CLL           aerosol, dry deposition at a rate determined from the        TRBDRY1A.193    
CLL           dry deposition velocity vd and mean mixed layer depth h.     TRBDRY1A.194    
CLL           Thus the max concentration is given by                       TRBDRY1A.195    
CLL               Q/(uh)*k/(k+vd/h)*(1-exp(-k*x/u))                        TRBDRY1A.196    
CLL           The source term is assumed to decrease with level            TRBDRY1A.197    
CLL           pressure. See forthcoming documentation for details.         TRBDRY1A.198    
CLL                                                                        TRBDRY1A.199    
CLL Pete Clark  <- programmer of some or all of previous code or changes   TRBDRY1A.200    
CLL                                                                        TRBDRY1A.201    
CLL  Model            Modification history from model version 3.4:         TRBDRY1A.202    
CLL version  Date                                                          TRBDRY1A.203    
CLL                                                                        TRBDRY1A.204    
CLL  Programming standard: Unified Model Documentation Paper No 3,         TRBDRY1A.205    
CLL                        Version 7, dated 11/3/93.                       TRBDRY1A.206    
CLL                                                                        TRBDRY1A.207    
CLL                                                                        TRBDRY1A.208    
C*L  Arguments:---------------------------------------------------------   TRBDRY1A.209    

      REAL FUNCTION BDRYV(                                                  4TRBDRY1A.210    
     & WDIR                                                                TRBDRY1A.211    
     &,WSPEED                                                              TRBDRY1A.212    
     &,PRESS                                                               TRBDRY1A.213    
     &)                                                                    TRBDRY1A.214    
      REAL                                                                 TRBDRY1A.215    
     & WDIR          ! IN Wind direction : Cartesian degrees               TRBDRY1A.216    
     &,WSPEED        ! IN Wind speed m/s                                   TRBDRY1A.217    
     &,PRESS         ! IN Pressure                                         TRBDRY1A.218    
C*                                                                         TRBDRY1A.219    
C* Local, including SAVE'd, storage------------------------------------    TRBDRY1A.220    
C                                                                          TRBDRY1A.221    
C  (a) Scalars effectively expanded to workspace by the Cray (using        TRBDRY1A.222    
C      vector registers).                                                  TRBDRY1A.223    
      REAL                                                                 TRBDRY1A.224    
     & ZANGLE1,WIDTH1 ! Centre and width of first source Gaussian.         TRBDRY1A.225    
     & ZANGLE2,WIDTH2 ! Centre and width of second source Gaussian.        TRBDRY1A.226    
     &,CMAX,CZER ! Max concentration and 'background'.                     TRBDRY1A.227    
     &,WDIRN,RECIPROOT2PI                                                  TRBDRY1A.228    
     &,MIXD,TRAVEL    ! Average mixed layer depth and travel distance.     TRBDRY1A.229    
     &,QMAX1          ! Peak height of first source Gaussian.              TRBDRY1A.230    
     &,QMAX2          ! Peak height of second source Gaussian.             TRBDRY1A.231    
     &,VD             ! Dry deposition velocity.                           TRBDRY1A.232    
     &,K,K1 ! Transformation parameters.                                   TRBDRY1A.233    
     &,PH ! Pressure height scale.                                         TRBDRY1A.234    
     &,KRAT,KT                                                             TRBDRY1A.235    
      PARAMETER(ZANGLE1=178.0)                                             TRBDRY1A.236    
      PARAMETER(WIDTH1=5.0)                                                TRBDRY1A.237    
      PARAMETER(ZANGLE2=173.0)                                             TRBDRY1A.238    
      PARAMETER(WIDTH2=25.0)                                               TRBDRY1A.239    
      PARAMETER(RECIPROOT2PI=0.3989422803)                                 TRBDRY1A.240    
      PARAMETER(QMAX1=4.7E5,QMAX2=9.0E4,MIXD=800.0,TRAVEL=7.7E5)           TRBDRY1A.241    
      PARAMETER(CZER=6.0)                                                  TRBDRY1A.242    
      PARAMETER(VD=5.0E-3)                                                 TRBDRY1A.243    
      PARAMETER(K=3.0E-6, K1=K+VD/MIXD, KRAT=K/K1,KT=-K1*TRAVEL)           TRBDRY1A.244    
      PARAMETER(PH=3.0E4)                                                  TRBDRY1A.245    
C                                                                          TRBDRY1A.246    
*CALL C_PI                                                                 TRBDRY1A.247    
!                                                                          TRBDRY1A.248    
!     Max concentration = Q/(uh)*k/(k+vd/h)*(1-exp(-k*x/u))                TRBDRY1A.249    
!                                                                          TRBDRY1A.250    
      CMAX=1.0/MAX(WSPEED,1.0)/MIXD                                        TRBDRY1A.251    
      CMAX=CMAX*KRAT*(1-EXP(KT/MAX(WSPEED,1.0)))                           TRBDRY1A.252    
      WDIRN = WDIR - ZANGLE1                                               TRBDRY1A.253    
      IF (WDIRN .LT. -180.0) WDIRN=WDIRN+360.0                             TRBDRY1A.254    
      IF (WDIRN .GT.  180.0) WDIRN=WDIRN-360.0                             TRBDRY1A.255    
      WDIRN=WDIRN/WIDTH1                                                   TRBDRY1A.256    
      BDRYV= QMAX1 * EXP(-WDIRN*WDIRN/2.0)                                 TRBDRY1A.257    
      WDIRN = WDIR*RECIP_PI_OVER_180 - ZANGLE2                             TRBDRY1A.258    
      IF (WDIRN .LT. -180.0) WDIRN=WDIRN+360.0                             TRBDRY1A.259    
      IF (WDIRN .GT.  180.0) WDIRN=WDIRN-360.0                             TRBDRY1A.260    
      WDIRN=WDIRN/WIDTH2                                                   TRBDRY1A.261    
      BDRYV= BDRYV + QMAX2 * EXP(-WDIRN*WDIRN/2.0)                         TRBDRY1A.262    
!                                                                          TRBDRY1A.263    
!     Add 'background' value.                                              TRBDRY1A.264    
!                                                                          TRBDRY1A.265    
      BDRYV= BDRYV * CMAX + CZER                                           TRBDRY1A.266    
!                                                                          TRBDRY1A.267    
!     Reduce concentration with pressure altitude.                         TRBDRY1A.268    
!                                                                          TRBDRY1A.269    
      BDRYV= BDRYV * EXP(-(1.E5-PRESS)/PH)                                 TRBDRY1A.270    
      RETURN                                                               TRBDRY1A.271    
      END                                                                  TRBDRY1A.272    
*ENDIF                                                                     TRBDRY1A.273