*IF DEF,A16_1A DIATHA1A.2
C ******************************COPYRIGHT****************************** GTS2F400.2161
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.2162
C GTS2F400.2163
C Use, duplication or disclosure of this code is subject to the GTS2F400.2164
C restrictions as set forth in the contract. GTS2F400.2165
C GTS2F400.2166
C Meteorological Office GTS2F400.2167
C London Road GTS2F400.2168
C BRACKNELL GTS2F400.2169
C Berkshire UK GTS2F400.2170
C RG12 2SZ GTS2F400.2171
C GTS2F400.2172
C If no contract has been raised with this copy of the code, the use, GTS2F400.2173
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.2174
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.2175
C Modelling at the above address. GTS2F400.2176
C ******************************COPYRIGHT****************************** GTS2F400.2177
C GTS2F400.2178
CLL SUBROUTINE THADV--------------------------------------------------- DIATHA1A.3
CLL DIATHA1A.4
CLL PURPOSE: Calculates thermal advection at pressure levels DIATHA1A.5
CLL Tested under compiler CFT77 DIATHA1A.6
CLL Tested under OS version 5.1 DIATHA1A.7
CLL DIATHA1A.8
CLL D.Robinson <- programmer of some or all of previous code or changes DIATHA1A.9
CLL DIATHA1A.10
CLL Model Modification history from model version 3.0: DIATHA1A.11
CLL version Date DIATHA1A.12
CLL 4.4 23/06/97 Initialise TH_TP. D. Robinson. GDR2F404.1
!LL 4.5 15/04/98 Added start-end arguments to V_INT and GSM1F405.610
!LL V_INT_T routines. S.D.Mullerworth GSM1F405.611
CLL DIATHA1A.13
CLL Logical components covered: D424 DIATHA1A.14
CLL DIATHA1A.15
CLL Project TASK: DIATHA1A.16
CLL DIATHA1A.17
CLL Programming standard: U M DOC Paper NO. 4, DIATHA1A.18
CLL DIATHA1A.19
CLL External documentation: DIATHA1A.20
CLL DIATHA1A.21
CLLEND------------------------------------------------------------------ DIATHA1A.22
C DIATHA1A.23
C*L ARGUMENTS:--------------------------------------------------------- DIATHA1A.24
SUBROUTINE DIA_THADV( 2,4DIATHA1A.25
C data in DIATHA1A.26
& U,V,THETA,PPT,PUV,PSTAR,PRESS_REQD,P_EXNER_HALF, DIATHA1A.27
C data out DIATHA1A.28
& TH_TP, DIATHA1A.29
C constants in DIATHA1A.30
& P_FIELD,U_FIELD,P_LEVELS,ROW_LENGTH,P_ROWS,SEC_U_LATITUDE, DIATHA1A.31
& EW_SPACE,NS_SPACE,AKH,BKH) DIATHA1A.32
C* DIATHA1A.33
C*L DIATHA1A.34
C----------------------------------------------------------------------- DIATHA1A.35
IMPLICIT NONE DIATHA1A.36
C----------------------------------------------------------------------- DIATHA1A.37
EXTERNAL V_INT_T,V_INT,UV_TO_P_FULL DIATHA1A.38
C----------------------------------------------------------------------- DIATHA1A.39
INTEGER DIATHA1A.40
* P_FIELD ! IN NO OF POINTS ON P/T GRID DIATHA1A.41
*,U_FIELD ! IN NO OF POINTS ON U/V GRID DIATHA1A.42
*,P_LEVELS ! IN NO OF MODEL LEVELS DIATHA1A.43
*,P_ROWS ! IN NO OF ROWS ON P/T GRID DIATHA1A.44
*,ROW_LENGTH ! IN NO OF COLUMNS DIATHA1A.45
C----------------------------------------------------------------------- DIATHA1A.46
REAL DIATHA1A.47
* U(U_FIELD,P_LEVELS) ! IN U FIELD AT FULL LEVELS DIATHA1A.48
*,V(U_FIELD,P_LEVELS) ! IN V FIELD AT FULL LEVELS DIATHA1A.49
*,THETA(P_FIELD,P_LEVELS) ! IN POTENTIAL TEMPERATURE " " DIATHA1A.50
*,PRESS_REQD ! IN PRESSURE AT WHICH TO DIATHA1A.51
* ! CALCULATE THERMAL ADVECTION DIATHA1A.52
*,TH_TP(P_FIELD) ! OUT THERMAL ADVECTION AT DIATHA1A.53
* ! REQUIRED LEVEL DIATHA1A.54
*,PPT(P_FIELD,P_LEVELS) ! IN PRESS FIELD AT P/T POINTS DIATHA1A.55
*,PUV(U_FIELD,P_LEVELS) ! IN PRESS FIELD AT U/V POINTS DIATHA1A.56
*,PSTAR(P_FIELD) ! IN SURFACE PRESSURE FIELD DIATHA1A.57
*,P_EXNER_HALF(P_FIELD,P_LEVELS+1) ! IN EXNER PRESSURE AT MODEL DIATHA1A.58
* ! HALF LEVELS DIATHA1A.59
*,SEC_U_LATITUDE(U_FIELD) ! IN 1/COS(LAT) AT U/V POINTS DIATHA1A.60
*,EW_SPACE ! IN LONGITUDE GRID SPACING DIATHA1A.61
*,NS_SPACE ! IN LATITUDE GRID SPACING DIATHA1A.62
* ! BOTH THE ABOVE IN DEGREES DIATHA1A.63
*,AKH,BKH ! IN Hybrid Coords. A and B DIATHA1A.64
* ! values on half levels. DIATHA1A.65
C* DIATHA1A.66
C*L DIATHA1A.67
C----------------------------------------------------------------------- DIATHA1A.68
C Local Variables DIATHA1A.69
C----------------------------------------------------------------------- DIATHA1A.70
INTEGER DIATHA1A.71
* I,J,K ! LOOP COUNTERS DIATHA1A.72
*,T_REF ! REFERENCE LEVEL FOR VERTICAL INTERPOLATION DIATHA1A.73
* ! OF TEMPERATURES DIATHA1A.74
C----------------------------------------------------------------------- DIATHA1A.75
REAL DIATHA1A.76
* LONG_SI_OVER_TWO_A ! LONGITUDE STEP INVERSE (IN RADIANS) DIATHA1A.77
* ! DIVIDED BY (2*EARTH'S RADIUS(A)) DIATHA1A.78
*,LAT_SI_OVER_TWO_A ! LATITUDE STEP INVERSE (IN RADIANS) DIATHA1A.79
* ! DIVIDED BY (2*EARTH'S RADIUS(A)) DIATHA1A.80
*,T_P(P_FIELD) ! TEMPERATURE FIELD ON REQUIRED LEVEL DIATHA1A.81
*,U_P(U_FIELD) ! U FIELD ON REQUIRED LEVEL DIATHA1A.82
*,V_P(U_FIELD) ! V FIELD ON REQUIRED LEVEL DIATHA1A.83
*,TH_UV(U_FIELD) ! THERMAL ADVECTION AT U/V POINTS DIATHA1A.84
*,PPT_REQ(P_FIELD) ! PRESSURE SURFACE REQUIRED ARRAY DIATHA1A.85
*,DUMMY1 ! DUMMY ARGUMENT1 DIATHA1A.86
*,DUMMY2 ! DUMMY ARGUMENT2 DIATHA1A.87
*,TERM1 ! 1st TERM USED IN CALCULATION DIATHA1A.88
*,TERM2 ! 2nd TERM USED IN CALCULATION DIATHA1A.89
C----------------------------------------------------------------------- DIATHA1A.90
C Constants required : A=radius of Earth, DIATHA1A.91
C RECIP_PI_OVER_180=180/Pi DIATHA1A.92
C----------------------------------------------------------------------- DIATHA1A.93
*CALL C_A 
DIATHA1A.94
*CALL C_PI DIATHA1A.95
T_REF=2 ! SET TO SECOND MODEL LEVEL DIATHA1A.96
C*---------------------------------------------------------------------- DIATHA1A.97
C----------------------------------------------------------------------- GDR2F404.2
CL Initialise TH_TP GDR2F404.3
C----------------------------------------------------------------------- GDR2F404.4
DO I=1,P_FIELD GDR2F404.5
TH_TP(I)=0.0 GDR2F404.6
ENDDO GDR2F404.7
C----------------------------------------------------------------------- DIATHA1A.98
CL Set up Pressure required array DIATHA1A.99
C----------------------------------------------------------------------- DIATHA1A.100
DO I=1,P_FIELD DIATHA1A.101
PPT_REQ(I)=PRESS_REQD DIATHA1A.102
ENDDO DIATHA1A.103
C----------------------------------------------------------------------- DIATHA1A.104
CL Interpolate the temperature field onto the required level DIATHA1A.105
C----------------------------------------------------------------------- DIATHA1A.106
CALL V_INT_T(T_P,PPT_REQ,PPT(1,T_REF),PSTAR,P_EXNER_HALF, DIATHA1A.107
&THETA,P_FIELD,P_LEVELS,T_REF,AKH,BKH,1,P_FIELD) GSM1F405.612
C----------------------------------------------------------------------- DIATHA1A.109
CL Interpolate the U field onto the required level DIATHA1A.110
C----------------------------------------------------------------------- DIATHA1A.111
CALL V_INT(PUV,PPT_REQ,U,U_P, DIATHA1A.112
&U_FIELD,P_LEVELS,DUMMY1,DUMMY2,.FALSE.,1,U_FIELD) GSM1F405.613
C----------------------------------------------------------------------- DIATHA1A.114
CL Interpolate the V field onto the required level DIATHA1A.115
C----------------------------------------------------------------------- DIATHA1A.116
CALL V_INT(PUV,PPT_REQ,V,V_P, DIATHA1A.117
& U_FIELD,P_LEVELS,DUMMY1,DUMMY2,.FALSE.,1,U_FIELD) GSM1F405.614
C----------------------------------------------------------------------- DIATHA1A.119
CL Calculation of thermal advection DIATHA1A.120
C----------------------------------------------------------------------- DIATHA1A.121
C DIATHA1A.122
C *T4 _ *T1 Y = Latitude DIATHA1A.123
C ^ X = Longitude DIATHA1A.124
C dY dY= Latitude interval in radians DIATHA1A.125
C ^ dX= Longitude interval in radians DIATHA1A.126
C <---dX--*U,V---> A = Radius of Earth in metres DIATHA1A.127
C ^ T1-4= Temperatures in Kelvin (T_P) DIATHA1A.128
C ^ U,V= Wind components in m/s (U_P,V_P) DIATHA1A.129
C ^ DIATHA1A.130
C *T3 _ *T2 DIATHA1A.131
C DIATHA1A.132
C Term1 Term2 DIATHA1A.133
C d(Temp) ( U*(T1+T2-T3-T4) V*(T4+T1-T2-T3)) DIATHA1A.134
C Thermal advection = ------- = -( --------------- + ---------------) DIATHA1A.135
C Degrees/second d(Time) ( Cos(Y)*dX*2*A 2*dY*A ) DIATHA1A.136
C DIATHA1A.137
C----------------------------------------------------------------------- DIATHA1A.138
CL Calculate 1/(dX*2*A) and 1/(dY*2*A) DIATHA1A.139
C----------------------------------------------------------------------- DIATHA1A.140
LONG_SI_OVER_TWO_A = 0.5*RECIP_PI_OVER_180/EW_SPACE/A DIATHA1A.141
LAT_SI_OVER_TWO_A = 0.5*RECIP_PI_OVER_180/NS_SPACE/A DIATHA1A.142
C----------------------------------------------------------------------- DIATHA1A.143
DO I=1,U_FIELD-1 DIATHA1A.144
C----------------------------------------------------------------------- DIATHA1A.145
C Calculate terms for all cases DIATHA1A.146
C----------------------------------------------------------------------- DIATHA1A.147
TERM1=U_P(I)*(T_P(I+1)+T_P(I+1+ROW_LENGTH)-T_P(I+ROW_LENGTH) DIATHA1A.148
* -T_P(I))*SEC_U_LATITUDE(I)*LONG_SI_OVER_TWO_A DIATHA1A.149
C DIATHA1A.150
TERM2=V_P(I)*(T_P(I)+T_P(I+1)-T_P(I+1+ROW_LENGTH) DIATHA1A.151
* -T_P(I+ROW_LENGTH))*LAT_SI_OVER_TWO_A DIATHA1A.152
TH_UV(I)=-(TERM1+TERM2) DIATHA1A.153
ENDDO DIATHA1A.154
*IF DEF,GLOBAL DIATHA1A.155
C----------------------------------------------------------------------- DIATHA1A.156
C Recalculate terms for end of each row in the Global model DIATHA1A.157
C----------------------------------------------------------------------- DIATHA1A.158
DO I=ROW_LENGTH,U_FIELD,ROW_LENGTH DIATHA1A.159
TERM1=U_P(I)*(T_P(I+1-ROW_LENGTH)+T_P(I+1)-T_P(I+ROW_LENGTH) DIATHA1A.160
* -T_P(I))*SEC_U_LATITUDE(I)*LONG_SI_OVER_TWO_A DIATHA1A.161
C DIATHA1A.162
TERM2=V_P(I)*(T_P(I)+T_P(I+1-ROW_LENGTH)-T_P(I+1) DIATHA1A.163
* -T_P(I+ROW_LENGTH))*LAT_SI_OVER_TWO_A DIATHA1A.164
TH_UV(I)=-(TERM1+TERM2) DIATHA1A.165
ENDDO DIATHA1A.166
*ELSE DIATHA1A.167
C----------------------------------------------------------------------- DIATHA1A.168
CL For regional model first and last values in each row recalculated DIATHA1A.169
CL by extrapolation of adjacent two points DIATHA1A.170
C----------------------------------------------------------------------- DIATHA1A.171
DO I=ROW_LENGTH,U_FIELD,ROW_LENGTH DIATHA1A.172
J=I-ROW_LENGTH+1 DIATHA1A.173
TH_UV(I)=2.*TH_UV(I-1)-TH_UV(I-2) DIATHA1A.174
TH_UV(J)=2.*TH_UV(J+1)-TH_UV(J+2) DIATHA1A.175
ENDDO DIATHA1A.176
*ENDIF DIATHA1A.177
DIATHA1A.178
C----------------------------------------------------------------------- DIATHA1A.179
CL Interpolate thermal advection from U/V to P/T points DIATHA1A.180
C----------------------------------------------------------------------- DIATHA1A.181
CALL UV_TO_P_FULL DIATHA1A.182
1(TH_UV,TH_TP,U_FIELD,P_FIELD,ROW_LENGTH,P_ROWS) DIATHA1A.183
C----------------------------------------------------------------------- DIATHA1A.184
DIATHA1A.185
RETURN DIATHA1A.186
END DIATHA1A.187
*ENDIF DIATHA1A.188