*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