SUBROUTINE OMEGA_DIAG 1,6OMEGDI1A.30
1 ( GSM3F403.67
*CALL ARGFLDPT 
GSM3F403.68
& U,V,OMEGA,SEC_P_LATITUDE,COS_U_LATITUDE, GSM3F403.69
2 PSTAR,PSTAR_OLD,DELTA_AK,DELTA_BK, OMEGDI1A.32
3 AK,BK,AKH,BKH,U_FIELD,P_FIELD,P_LEVELS, OMEGDI1A.33
4 ROW_LENGTH,LATITUDE_STEP_INVERSE, OMEGDI1A.34
5 LONGITUDE_STEP_INVERSE,ADVECTION_TIMESTEP) OMEGDI1A.35
OMEGDI1A.36
IMPLICIT NONE OMEGDI1A.37
OMEGDI1A.38
*CALL TYPFLDPT GSM3F403.70
GSM3F403.71
INTEGER OMEGDI1A.39
& ROW_LENGTH !IN NUMBER OF POINTS PER ROW OMEGDI1A.40
&, P_LEVELS !IN NUMBER OF PRESSURE LEVELS OF DATA OMEGDI1A.41
&, P_FIELD !IN NUMBER OF POINTS IN PRESSURE FIELD OMEGDI1A.42
&, U_FIELD !IN NUMBER OF POINTS IN VELOCITY FIELD OMEGDI1A.43
OMEGDI1A.44
REAL OMEGDI1A.45
& U(U_FIELD,P_LEVELS) !IN U VELOCITY OMEGDI1A.46
&,V(U_FIELD,P_LEVELS) !IN V VELOCITY OMEGDI1A.47
&,PSTAR(P_FIELD) !IN SURFACE PRESSURE AT P POINTS. OMEGDI1A.48
&,PSTAR_OLD(P_FIELD) !IN SURFACE PRESSURE AT PREVIOUS TIMESTEP OMEGDI1A.49
&,SEC_P_LATITUDE(P_FIELD)!IN 1/COS(LAT) AT P POINTS OMEGDI1A.50
&,COS_U_LATITUDE(U_FIELD)!IN COS(LAT) AT U POINTS OMEGDI1A.51
&,LONGITUDE_STEP_INVERSE !IN 1/LONGITUDE INCREMENT OMEGDI1A.52
&,LATITUDE_STEP_INVERSE !IN 1/LATITUDE INCREMENT OMEGDI1A.53
&,BK(P_LEVELS) !IN HOLDS COEFFICIENT WHICH OMEGDI1A.54
& ! MULTIPLIES PSTAR IN HYBRID CO-ORDS OMEGDI1A.55
& ! AT LEVELS K OMEGDI1A.56
&,AK(P_LEVELS) !IN HOLDS FIRST COEFFICIENT OMEGDI1A.57
& ! IN HYBRID CO-ORDS AT LEVELS K OMEGDI1A.58
&,BKH(P_LEVELS+1) !IN HOLDS COEFFICIENT WHICH OMEGDI1A.59
& ! MULTIPLIES PSTAR IN HYBRID CO-ORDS OMEGDI1A.60
& ! AT half LEVELS K OMEGDI1A.61
&,AKH(P_LEVELS+1) !IN HOLDS FIRST COEFFICIENT OMEGDI1A.62
& ! IN HYBRID CO-ORDS AT half LEVELS K OMEGDI1A.63
&,DELTA_AK(P_LEVELS) !IN AK(K+1/2)-AK(K-1/2) OMEGDI1A.64
&,DELTA_BK(P_LEVELS) !IN BK(K+1/2)-BK(K-1/2) OMEGDI1A.65
&,ADVECTION_TIMESTEP !IN OMEGDI1A.66
OMEGDI1A.67
OMEGDI1A.68
REAL OMEGDI1A.69
& OMEGA(U_FIELD,P_LEVELS) !OUT. OMEGA AT U POINTS. OMEGDI1A.70
C*--------------------------------------------------------------------- OMEGDI1A.71
OMEGDI1A.72
C*L 9 LOCAL ARRAYS NEEDED. ----------------------------------------- OMEGDI1A.73
OMEGDI1A.74
REAL OMEGDI1A.75
& ETADOT(P_FIELD,P_LEVELS+1) ! HOLDS ETADOT OR PARTS THEREOF. OMEGDI1A.76
&, WORK1(P_FIELD) ! \ OMEGDI1A.77
&, WORK2(P_FIELD) ! > GENERAL WORK ARRAYS. OMEGDI1A.78
&, WORK3(P_FIELD) ! / OMEGDI1A.79
&, PSTAR_UV(U_FIELD) ! PSTAR AT U POINTS. OMEGDI1A.80
&, DELTA_P(U_FIELD) ! DEPTH OF LAYER AT U POINTS. OMEGDI1A.81
&, WP(U_FIELD) ! TERM IN OMEGA EQUATION. OMEGDI1A.82
&, DELTA_AKH(P_LEVELS+1) ! AK(K) -AK(K-1) OMEGDI1A.83
&, DELTA_BKH(P_LEVELS+1) ! BK(K) - BK(K-1) OMEGDI1A.84
C*--------------------------------------------------------------------- OMEGDI1A.85
OMEGDI1A.86
C DEFINE COUNT VARIABLES FOR DO LOOPS ETC. OMEGDI1A.87
INTEGER OMEGDI1A.88
& I,J,K OMEGDI1A.89
&, START,END OMEGDI1A.90
C DEFINE LOCAL SCALARS OMEGDI1A.91
REAL OMEGDI1A.92
& SCALAR OMEGDI1A.93
OMEGDI1A.94
*IF DEF,GLOBAL OMEGDI1A.95
REAL SUM_N(P_LEVELS),SUM_S(P_LEVELS) GSM3F403.72
*ENDIF OMEGDI1A.97
OMEGDI1A.98
*IF DEF,MPP GSM3F403.73
REAL GSM3F403.74
& SUM_N_ARRAY(ROW_LENGTH,P_LEVELS) GSM3F403.75
& ,SUM_S_ARRAY(ROW_LENGTH,P_LEVELS) GSM3F403.76
GSM3F403.77
INTEGER GSM3F403.78
& info GSM3F403.79
GSM3F403.80
*ENDIF GSM3F403.81
C*L EXTERNAL SUBROUTINE CALLS:--------------------------------------- OMEGDI1A.99
EXTERNAL P_TO_UV OMEGDI1A.100
C*--------------------------------------------------------------------- OMEGDI1A.101
*CALL C_A OMEGDI1A.102
OMEGDI1A.103
CL MAXIMUM VECTOR LENGTH ASSUMED IS P_FIELD OMEGDI1A.104
CL--------------------------------------------------------------------- OMEGDI1A.105
CL INTERNAL STRUCTURE. OMEGDI1A.106
CL--------------------------------------------------------------------- OMEGDI1A.107
CL OMEGDI1A.108
CL--------------------------------------------------------------------- OMEGDI1A.109
CL SECTION 0. SET UP PSTAR AT U POINTS AND DELTA_AKH,_BKH OMEGDI1A.110
CL--------------------------------------------------------------------- OMEGDI1A.111
CL OMEGDI1A.112
OMEGDI1A.113
CL CALL P_TO_UV TO INTERPOLATE PSTAR TO U POINTS. OMEGDI1A.114
OMEGDI1A.115
CALL P_TO_UV(PSTAR,PSTAR_UV,P_FIELD,U_FIELD,ROW_LENGTH, OMEGDI1A.116
& P_FIELD/ROW_LENGTH) OMEGDI1A.117
*IF DEF,MPP GSM3F403.82
CALL SWAPBOUNDS(PSTAR_UV,ROW_LENGTH,tot_U_ROWS, GSM3F403.83
& EW_Halo,NS_Halo,1) GSM3F403.84
*ENDIF GSM3F403.85
OMEGDI1A.118
C SET UP DELTA_AKH,_BKH OMEGDI1A.119
DO K=2,P_LEVELS OMEGDI1A.121
DELTA_AKH(K) = AK(K) - AK(K-1) OMEGDI1A.122
DELTA_BKH(K) = BK(K) - BK(K-1) OMEGDI1A.123
END DO OMEGDI1A.124
DELTA_AKH(1) = 0. OMEGDI1A.125
DELTA_AKH(P_LEVELS+1) = 0. OMEGDI1A.126
DELTA_BKH(1) = 0. OMEGDI1A.127
DELTA_BKH(P_LEVELS+1) = 0. OMEGDI1A.128
OMEGDI1A.129
C LOOP OVER LEVELS OMEGDI1A.130
DO 100 K=1,P_LEVELS OMEGDI1A.131
OMEGDI1A.132
CL--------------------------------------------------------------------- OMEGDI1A.133
CL SECTION 1. CALCULATE DIVERGENCE AS IN EQUATION (30). OMEGDI1A.134
CL--------------------------------------------------------------------- OMEGDI1A.135
OMEGDI1A.136
DO I=1,U_FIELD OMEGDI1A.137
DELTA_P(I) = DELTA_AK(K)+DELTA_BK(K)*PSTAR_UV(I) OMEGDI1A.138
END DO OMEGDI1A.139
OMEGDI1A.140
C CALCULATE DU/D(LAMDA), STORE IN ETADOT. OMEGDI1A.141
DO I=2,U_FIELD OMEGDI1A.142
WORK1(I) = LONGITUDE_STEP_INVERSE*(U(I,K)*DELTA_P(I) OMEGDI1A.143
& -U(I-1,K)*DELTA_P(I-1)) OMEGDI1A.144
END DO OMEGDI1A.145
OMEGDI1A.146
C CALCULATE DV/D(PHI), store in work2. OMEGDI1A.147
DO I=ROW_LENGTH+1,U_FIELD OMEGDI1A.148
WORK2(I) = LATITUDE_STEP_INVERSE*((V(I-ROW_LENGTH,K) GSM3F403.86
& *COS_U_LATITUDE(I-ROW_LENGTH))*DELTA_P(I-ROW_LENGTH) GSM3F403.87
& -(V(I,K)*COS_U_LATITUDE(I))*DELTA_P(I)) GSM3F403.88
END DO OMEGDI1A.152
OMEGDI1A.153
*IF DEF,GLOBAL OMEGDI1A.154
C CALCULATE AVERAGE OF DV/D(PHI), STORE IN WORK3. OMEGDI1A.155
DO I=ROW_LENGTH+2,U_FIELD OMEGDI1A.156
WORK3(I) = WORK2(I) + WORK2(I-1) OMEGDI1A.157
END DO OMEGDI1A.158
OMEGDI1A.159
*IF -DEF,MPP GSM3F403.89
C NOW DO FIRST POINT ON EACH SLICE FOR DU/D(LAMDA) AND DV/D(PHI) OMEGDI1A.160
WORK1(1) = LONGITUDE_STEP_INVERSE * (U(1,K)*DELTA_P(1) OMEGDI1A.161
& - U(ROW_LENGTH,K)*DELTA_P(ROW_LENGTH)) OMEGDI1A.162
DO I=ROW_LENGTH+1,U_FIELD,ROW_LENGTH OMEGDI1A.163
WORK1(I) = LONGITUDE_STEP_INVERSE * (U(I,K)*DELTA_P(I) OMEGDI1A.164
& - U(I+ROW_LENGTH-1,K)*DELTA_P(I+ROW_LENGTH-1)) OMEGDI1A.165
WORK3(I)=WORK2(I)+WORK2(I-1+ROW_LENGTH) OMEGDI1A.166
END DO OMEGDI1A.167
*ENDIF GSM3F403.90
C CALCULATE DIVERGENCES. STORE IN ETADOT. OMEGDI1A.169
OMEGDI1A.170
*IF DEF,MPP GSM3F403.91
DO I=START_POINT_NO_HALO+EW_Halo,END_P_POINT_NO_HALO-EW_Halo GSM3F403.92
*ELSE GSM3F403.93
DO I=ROW_LENGTH+1,U_FIELD OMEGDI1A.171
*ENDIF GSM3F403.94
ETADOT(I,K)= (SEC_P_LATITUDE(I)*.5)*(WORK1(I) GSM3F403.95
& + WORK1(I-ROW_LENGTH) OMEGDI1A.173
& + WORK3(I)) OMEGDI1A.174
END DO OMEGDI1A.175
*ELSE OMEGDI1A.176
WORK2(1) = 0. OMEGDI1A.177
OMEGDI1A.178
C CALCULATE DIVERGENCES. STORE IN ETADOT OMEGDI1A.179
OMEGDI1A.180
*IF DEF,MPP GSM3F403.96
DO I=START_POINT_NO_HALO+EW_Halo,END_P_POINT_NO_HALO-EW_Halo GSM3F403.97
*ELSE GSM3F403.98
DO I=ROW_LENGTH+2,U_FIELD-1 OMEGDI1A.181
*ENDIF GSM3F403.99
ETADOT(I,K)= (SEC_P_LATITUDE(I)*.5)*(WORK1(I) GSM3F403.100
& + WORK1(I-ROW_LENGTH) OMEGDI1A.183
& + WORK2(I) + WORK2(I-1)) OMEGDI1A.184
END DO OMEGDI1A.185
OMEGDI1A.186
C ZERO DIVERGENCES ON BOUNDARIES. OMEGDI1A.187
*IF DEF,MPP GSM3F403.101
IF (at_left_of_LPG) THEN GSM3F403.102
DO I=ROW_LENGTH+1+EW_Halo,P_FIELD,ROW_LENGTH GSM3F403.103
ETADOT(I,K) = 0. GSM3F403.104
ENDDO GSM3F403.105
ENDIF GSM3F403.106
IF (at_right_of_LPG) THEN GSM3F403.107
DO I=2*ROW_LENGTH-EW_Halo,P_FIELD,ROW_LENGTH GSM3F403.108
ETADOT(I,K) = 0. GSM3F403.109
ENDDO GSM3F403.110
ENDIF GSM3F403.111
*ELSE GSM3F403.112
DO I=ROW_LENGTH+1,U_FIELD,ROW_LENGTH OMEGDI1A.188
ETADOT(I,K) = 0. OMEGDI1A.189
ETADOT(I+ROW_LENGTH-1,K) = 0. OMEGDI1A.190
END DO OMEGDI1A.191
*ENDIF GSM3F403.113
*ENDIF OMEGDI1A.192
OMEGDI1A.193
OMEGDI1A.194
*IF DEF,GLOBAL OMEGDI1A.195
C CALCULATE DIVERGENCE AT POLES BY SUMMING DV/D(LAT) AROUND POLAR OMEGDI1A.196
C CIRCLE AND AVERAGING. OMEGDI1A.197
SCALAR = (8.*LATITUDE_STEP_INVERSE) * (LATITUDE_STEP_INVERSE GSM3F403.114
*IF DEF,MPP GSM3F403.115
& /GLOBAL_ROW_LENGTH) GSM3F403.116
*ELSE GSM3F403.117
& /ROW_LENGTH) GSM3F403.118
*ENDIF GSM3F403.119
SUM_N(K) = 0. GSM3F403.120
SUM_S(K) = 0. GSM3F403.121
*IF DEF,MPP GSM3F403.122
IF (at_top_of_LPG) THEN GSM3F403.123
DO I=TOP_ROW_START+EW_Halo,START_POINT_NO_HALO-1-EW_Halo GSM3F403.124
SUM_N_ARRAY(I-TOP_ROW_START,K)= GSM3F403.125
& -(V(I,K)*DELTA_P(I))*(COS_U_LATITUDE(I)*SCALAR) GSM3F403.126
ENDDO GSM3F403.127
ENDIF GSM3F403.128
IF (at_base_of_LPG) THEN GSM3F403.129
DO I=U_BOT_ROW_START+EW_Halo,LAST_U_FLD_PT-EW_Halo GSM3F403.130
SUM_S_ARRAY(I-U_BOT_ROW_START,K)= GSM3F403.131
& (V(I,K)*SCALAR)*(DELTA_P(I)*COS_U_LATITUDE(I)) GSM3F403.132
ENDDO GSM3F403.133
ENDIF GSM3F403.134
GSM3F403.135
*ELSE GSM3F403.136
DO I=1,ROW_LENGTH OMEGDI1A.202
SUM_N(K) = SUM_N(K) - (V(I,K)*DELTA_P(I)) GSM3F403.137
& *(COS_U_LATITUDE(I)*SCALAR) GSM3F403.138
SUM_S(K) = SUM_S(K) + (V(U_FIELD+I-ROW_LENGTH,K)*SCALAR) GSM3F403.139
& *(delta_p(u_field+i-row_length) GSM3F403.140
& *cos_u_latitude(u_field+i-row_length)) GSM3F403.141
END DO OMEGDI1A.207
*ENDIF GSM3F403.142
*ENDIF GSM3F403.143
100 CONTINUE GSM3F403.144
GSM3F403.145
*IF DEF,GLOBAL GSM3F403.146
*IF DEF,MPP GSM3F403.147
IF (at_top_of_LPG) THEN GSM3F403.148
CALL GCG_RVECSUMR(ROW_LENGTH,ROW_LENGTH-2*EW_Halo, GSM3F403.149
& EW_Halo,P_LEVELS,SUM_N_ARRAY,GC_ROW_GROUP,info,SUM_N) GSM3F403.150
ENDIF GSM3F403.151
IF (at_base_of_LPG) THEN GSM3F403.152
CALL GCG_RVECSUMR(ROW_LENGTH,ROW_LENGTH-2*EW_Halo, GSM3F403.153
& EW_Halo,P_LEVELS,SUM_S_ARRAY,GC_ROW_GROUP,info,SUM_S) GSM3F403.154
ENDIF GSM3F403.155
*ENDIF GSM3F403.156
GSM3F403.157
DO K=1,P_LEVELS GSM3F403.158
C DIVERGENCE AT POLES GSM3F403.159
*IF DEF,MPP GSM3F403.160
IF (at_top_of_LPG) THEN GSM3F403.161
ETADOT(START_POINT_NO_HALO+EW_Halo-1,K) = SUM_N(K) GSM3F403.162
ENDIF GSM3F403.163
IF (at_base_of_LPG) THEN GSM3F403.164
ETADOT(END_P_POINT_NO_HALO-EW_Halo+1,K) = SUM_S(K) GSM3F403.165
ENDIF GSM3F403.166
*ELSE GSM3F403.167
ETADOT(ROW_LENGTH,K) = SUM_N(K) GSM3F403.168
ETADOT(P_FIELD-ROW_LENGTH+1,K) = SUM_S(K) GSM3F403.169
*ENDIF OMEGDI1A.210
C END LOOP OVER LEVELS OMEGDI1A.211
ENDDO GSM3F403.170
*ENDIF GSM3F403.171
OMEGDI1A.213
CL OMEGDI1A.214
CL--------------------------------------------------------------------- OMEGDI1A.215
CL SECTION 2. CALCULATE VERTICAL VELOCITY. EQUATION (29). OMEGDI1A.216
CL--------------------------------------------------------------------- OMEGDI1A.217
OMEGDI1A.218
*IF DEF,GLOBAL OMEGDI1A.219
C IF GLOBAL PERFORM SECTION 2 OVER POLES AS WELL AS OTHER POINTS. OMEGDI1A.220
C ONLY ONE POLAR POINT NEEDS TO BE DONE. OMEGDI1A.221
*IF DEF,MPP GSM3F403.172
START = START_POINT_NO_HALO+EW_Halo GSM3F403.173
END = END_P_POINT_NO_HALO-EW_Halo GSM3F403.174
C At poles, extra point holds polar values GSM3F403.175
IF (at_top_of_LPG) START = START-1 GSM3F403.176
IF (at_base_of_LPG) END = END+1 GSM3F403.177
*ELSE GSM3F403.178
START = ROW_LENGTH OMEGDI1A.222
END = P_FIELD-ROW_LENGTH + 1 GSM3F403.179
*ENDIF GSM3F403.180
*ELSE GSM3F403.181
*IF DEF,MPP GSM3F403.182
START = START_POINT_NO_HALO+EW_Halo GSM3F403.183
END = END_P_POINT_NO_HALO-EW_Halo GSM3F403.184
*ELSE OMEGDI1A.224
START = ROW_LENGTH + 1 AMN1F400.3
END = P_FIELD-ROW_LENGTH GSM3F403.185
*ENDIF GSM3F403.186
*ENDIF OMEGDI1A.227
OMEGDI1A.228
C --------------------------------------------------------------------- OMEGDI1A.229
CL SECTION 2.1 SUM DIVERGENCES THROUGHOUT ATMOSPHERE. OMEGDI1A.230
C --------------------------------------------------------------------- OMEGDI1A.231
OMEGDI1A.232
C BY CODING THE SUMMATION AS FOLLOWS THE VALUES PUT INTO EACH LEVEL OMEGDI1A.233
C OF ETADOT ARE THE ONES NEEDED FOR THE SECOND SUMMATION TERM OMEGDI1A.234
C IN EQUATION 29, WHILE THE TOTAL SUM IS HELD IN ETADOT( ,1) OMEGDI1A.235
OMEGDI1A.236
DO 210 K=P_LEVELS-1,1,-1 OMEGDI1A.237
DO I=START,END OMEGDI1A.238
ETADOT(I,K)= ETADOT(I,K)+ETADOT(I,K+1) OMEGDI1A.239
END DO OMEGDI1A.240
210 CONTINUE OMEGDI1A.241
OMEGDI1A.242
C --------------------------------------------------------------------- OMEGDI1A.243
CL SECTION 2.2 CALCULATE MASS-WEIGHTED VERTICAL VELOCITY. OMEGDI1A.244
C --------------------------------------------------------------------- OMEGDI1A.245
OMEGDI1A.246
C DP/D(PSTAR) IS NOTHING MORE THAN THE BK COEFFICENT. OMEGDI1A.247
OMEGDI1A.248
DO 220 K= P_LEVELS,2,-1 OMEGDI1A.249
DO I= START, END OMEGDI1A.250
ETADOT(I,K)= ETADOT(I,K) - BKH(K) * ETADOT(I,1) OMEGDI1A.251
END DO OMEGDI1A.252
220 CONTINUE OMEGDI1A.253
OMEGDI1A.254
DO 230 K=1,P_LEVELS OMEGDI1A.255
*IF DEF,GLOBAL OMEGDI1A.256
C IF GLOBAL MODEL SET ALL POINTS AT POLES TO THE UNIQUE VALUE. OMEGDI1A.257
*IF DEF,MPP GSM3F403.187
C MPP: Above, polar values stored in START and END points GSM3F403.188
IF (at_top_of_LPG) THEN GSM3F403.189
DO I=TOP_ROW_START+EW_Halo,START-1 GSM3F403.190
ETADOT(I,K)=ETADOT(START,K) GSM3F403.191
ENDDO GSM3F403.192
ENDIF GSM3F403.193
IF (at_base_of_LPG) THEN GSM3F403.194
DO I=END+1,LAST_P_FLD_PT-EW_Halo GSM3F403.195
ETADOT(I,K)=ETADOT(END,K) GSM3F403.196
ENDDO GSM3F403.197
ENDIF GSM3F403.198
*ELSE GSM3F403.199
CDIR$ IVDEP OMEGDI1A.258
! Fujitsu vectorization directive GRB0F405.429
!OCL NOVREC GRB0F405.430
DO I=1,ROW_LENGTH-1 OMEGDI1A.259
ETADOT(I,K)= ETADOT(ROW_LENGTH,K) OMEGDI1A.260
ETADOT(END+I,K)= ETADOT(END,K) OMEGDI1A.261
END DO OMEGDI1A.262
*ENDIF GSM3F403.200
*ELSE OMEGDI1A.263
C SET NORTHERN AND SOUTHERN BOUNDARIES TO ZERO. OMEGDI1A.264
*IF DEF,MPP GSM3F403.201
IF (at_top_of_LPG) THEN GSM3F403.202
DO I=TOP_ROW_START+EW_Halo,START-1 GSM3F403.203
ETADOT(I,K) = 0. GSM3F403.204
ENDDO GSM3F403.205
ENDIF GSM3F403.206
IF (at_base_of_LPG) THEN GSM3F403.207
DO I=END+1,LAST_P_FLD_PT-EW_Halo GSM3F403.208
ETADOT(I,K) = 0. GSM3F403.209
ENDDO GSM3F403.210
ENDIF GSM3F403.211
*ELSE GSM3F403.212
DO I=1,ROW_LENGTH OMEGDI1A.265
ETADOT(I,K) = 0. OMEGDI1A.266
ETADOT(P_FIELD+1-I,K) = 0. OMEGDI1A.267
END DO OMEGDI1A.268
*ENDIF OMEGDI1A.269
*ENDIF GSM3F403.213
230 CONTINUE OMEGDI1A.270
OMEGDI1A.271
CL INTERPOLATE ETADOT TO U POINTS. OMEGDI1A.272
CL LOOP OVER LEVELS OMEGDI1A.273
*IF DEF,MPP GSM3F403.214
CALL SWAPBOUNDS(ETADOT(1,2),ROW_LENGTH,tot_P_ROWS,EW_Halo, GSM3F403.215
& NS_Halo,P_LEVELS-1) GSM3F403.216
*ENDIF GSM3F403.217
DO K=2,P_LEVELS OMEGDI1A.274
CALL P_TO_UV(ETADOT(1,K),WORK1,P_FIELD,U_FIELD,ROW_LENGTH, OMEGDI1A.275
& P_FIELD/ROW_LENGTH) OMEGDI1A.276
DO I=1,U_FIELD OMEGDI1A.277
ETADOT(I,K) = WORK1(I) OMEGDI1A.278
END DO OMEGDI1A.279
END DO OMEGDI1A.280
OMEGDI1A.281
CL SET ETADOT AT TOP AND BOTTOM TO ZERO. OMEGDI1A.282
DO I=1,P_FIELD OMEGDI1A.283
ETADOT(I,1) = 0. OMEGDI1A.284
ETADOT(I,P_LEVELS+1) = 0. OMEGDI1A.285
END DO OMEGDI1A.286
OMEGDI1A.287
CL LOOP OVER LEVELS OMEGDI1A.288
DO 300 K=1,P_LEVELS OMEGDI1A.289
CL OMEGDI1A.290
CL--------------------------------------------------------------------- OMEGDI1A.291
CL SECTION 3. CALCULATE DP/DT OMEGDI1A.292
CL--------------------------------------------------------------------- OMEGDI1A.293
OMEGDI1A.294
IF(BK(K).EQ.0.) THEN OMEGDI1A.295
C A CONSTANT PRESSURE LEVEL SO DP/DT IS ZERO. OMEGDI1A.296
DO I=1,U_FIELD OMEGDI1A.297
WORK3(I) = 0. OMEGDI1A.298
END DO OMEGDI1A.299
ELSE OMEGDI1A.300
C CALCULATE DP/DT. OMEGDI1A.301
SCALAR = BK(k)/ADVECTION_TIMESTEP OMEGDI1A.302
DO I=1,P_FIELD OMEGDI1A.303
WORK2(I) = (PSTAR(I)-PSTAR_OLD(I))*SCALAR OMEGDI1A.304
END DO OMEGDI1A.305
C INTERPOLATE DP/DT TO U POINTS. OMEGDI1A.306
CALL P_TO_UV(WORK2,WORK3,P_FIELD,U_FIELD,ROW_LENGTH, OMEGDI1A.307
& P_FIELD/ROW_LENGTH) OMEGDI1A.308
END IF OMEGDI1A.309
OMEGDI1A.310
CL--------------------------------------------------------------------- OMEGDI1A.311
CL SECTION 4. CALCULATE U.GRAD P OMEGDI1A.312
CL--------------------------------------------------------------------- OMEGDI1A.313
OMEGDI1A.314
C---------------------------------------------------------------------- OMEGDI1A.315
CL SECTION 4.1 CALCULATE U DP/D(LAMDA) OMEGDI1A.316
C---------------------------------------------------------------------- OMEGDI1A.317
OMEGDI1A.318
C CALCULATE U DP/D(LAMDA) BETWEEN P POINTS OMEGDI1A.319
SCALAR=(.5*LONGITUDE_STEP_INVERSE)*(BK(K)/A) GSM3F403.218
*IF DEF,MPP GSM3F403.219
DO I=START_POINT_NO_HALO+EW_Halo,END_P_POINT_NO_Halo-EW_Halo GSM3F403.220
*ELSE GSM3F403.221
DO I=ROW_LENGTH+1,P_FIELD-ROW_LENGTH-1 GSM3F403.222
*ENDIF GSM3F403.223
WORK1(I) = ((U(I,K)+U(I-ROW_LENGTH,K))* GSM3F403.224
& (PSTAR(I+1)-PSTAR(I)))*(SEC_P_LATITUDE(I)* GSM3F403.225
& SCALAR) GSM3F403.226
END DO OMEGDI1A.324
C SET POLAR/BOUNDARY VALUES TO ZERO. OMEGDI1A.325
*IF DEF,MPP GSM3F403.227
IF (at_top_of_LPG) THEN GSM3F403.228
DO I=TOP_ROW_START+EW_Halo,START_POINT_NO_HALO-1-EW_Halo GSM3F403.229
WORK1(I) = 0. GSM3F403.230
ENDDO GSM3F403.231
ENDIF GSM3F403.232
IF (at_base_of_LPG) THEN GSM3F403.233
DO I=P_BOT_ROW_START+EW_Halo,LAST_P_FLD_PT-EW_Halo GSM3F403.234
WORK1(I) = 0. GSM3F403.235
ENDDO GSM3F403.236
ELSE GSM3F403.237
C Calculate bottom row halo values to allow WP to be calculated GSM3F403.238
DO I=P_BOT_ROW_START+EW_Halo,LAST_P_VALID_PT-EW_Halo GSM3F403.239
WORK1(I) = ((U(I,K)+U(I-ROW_LENGTH,K))* GSM3F403.240
& (PSTAR(I+1)-PSTAR(I)))*(SEC_P_LATITUDE(I) GSM3F403.241
& *SCALAR) GSM3F403.242
ENDDO GSM3F403.243
ENDIF GSM3F403.244
*ELSE GSM3F403.245
DO I=1,ROW_LENGTH OMEGDI1A.326
WORK1(I) = 0. OMEGDI1A.327
WORK1(P_FIELD+1-I) = 0. OMEGDI1A.328
END DO OMEGDI1A.329
*ENDIF GSM3F403.246
*IF -DEF,MPP GSM3F403.247
*IF DEF,GLOBAL OMEGDI1A.331
C GLOBAL MODEL SO RECALCULATE VALUE AT END-POINT OMEGDI1A.332
DO I=ROW_LENGTH*2,U_FIELD,ROW_LENGTH OMEGDI1A.333
WORK1(I) = ((U(I,K)+U(I-ROW_LENGTH,K))* GSM3F403.248
& (PSTAR(I+1-ROW_LENGTH)-PSTAR(I)))*(SEC_P_LATITUDE(I) GSM3F403.249
& *SCALAR) GSM3F403.250
END DO OMEGDI1A.337
*ELSE OMEGDI1A.338
WORK1(U_FIELD) = 0. OMEGDI1A.339
*ENDIF OMEGDI1A.340
*ENDIF GSM3F403.251
OMEGDI1A.341
C CALCULATE U DP/D(LONGITUDE) AT U POINTS OMEGDI1A.342
*IF DEF,MPP GSM3F403.252
IF (at_top_of_LPG) THEN GSM3F403.253
DO I=TOP_ROW_START+EW_Halo,START_POINT_NO_HALO-1+EW_Halo GSM3F403.254
WP(I) = .5*WORK1(I+ROW_LENGTH) GSM3F403.255
ENDDO GSM3F403.256
ENDIF GSM3F403.257
DO I=START_POINT_NO_HALO+EW_Halo,LAST_U_FLD_PT-EW_Halo GSM3F403.258
*ELSE GSM3F403.259
DO I=1,U_FIELD OMEGDI1A.344
*ENDIF GSM3F403.260
WP(I) = .5*(WORK1(I)+WORK1(I+ROW_LENGTH)) OMEGDI1A.345
END DO OMEGDI1A.346
OMEGDI1A.347
C---------------------------------------------------------------------- OMEGDI1A.348
CL SECTION 4.2 CALCULATE V DP/D(PHI) AND HENCE U.GRAD P OMEGDI1A.349
C---------------------------------------------------------------------- OMEGDI1A.350
OMEGDI1A.351
C CALCULATE V DP/D(PHI) BETWEEN P POINTS. OMEGDI1A.352
SCALAR=(.5*LATITUDE_STEP_INVERSE)*(BK(K)/A) GSM3F403.261
*IF DEF,MPP GSM3F403.262
DO I=TOP_ROW_START+1,LAST_U_FLD_PT GSM3F403.263
*ELSE GSM3F403.264
DO I=2,U_FIELD OMEGDI1A.354
*ENDIF GSM3F403.265
WORK2(I) = ((V(I,K)+V(I-1,K))* GSM3F403.266
& (PSTAR(I)-PSTAR(I+ROW_LENGTH))) GSM3F403.267
& *SCALAR GSM3F403.268
END DO OMEGDI1A.358
*IF -DEF,MPP GSM3F403.269
*IF DEF,GLOBAL OMEGDI1A.360
OMEGDI1A.361
C GLOBAL MODEL SO RECALCULATE FIRST POINT VALUES FOR V DP/D(LAT) OMEGDI1A.362
OMEGDI1A.363
DO I=1,U_FIELD,ROW_LENGTH OMEGDI1A.364
WORK2(I) = ((V(I,k)+V(I-1+ROW_LENGTH,k))*(PSTAR(I)- GSM3F403.270
& PSTAR(I+ROW_LENGTH)))*SCALAR GSM3F403.271
END DO OMEGDI1A.367
*ELSE OMEGDI1A.368
WORK2(1) = 0. OMEGDI1A.369
OMEGDI1A.370
*ENDIF OMEGDI1A.371
*ENDIF GSM3F403.272
OMEGDI1A.372
CL CALCULATE V DP/D(LAT) AT U POINTS AND ADD TO WP. OMEGDI1A.373
*IF DEF,MPP GSM3F403.273
DO I=TOP_ROW_START+1,LAST_U_FLD_PT-EW_Halo GSM3F403.274
*ELSE GSM3F403.275
DO I=1,U_FIELD-1 OMEGDI1A.374
*ENDIF GSM3F403.276
WORK1(I) = .5*(WORK2(I)+WORK2(I+1)) OMEGDI1A.375
END DO OMEGDI1A.376
OMEGDI1A.377
*IF -DEF,MPP GSM3F403.277
*IF DEF,GLOBAL OMEGDI1A.378
C RE-CALCULATE END POINT VALUE. OMEGDI1A.379
DO I=ROW_LENGTH,U_FIELD,ROW_LENGTH OMEGDI1A.380
WORK1(I) = .5*(WORK2(I)+WORK2(I+1-ROW_LENGTH)) OMEGDI1A.381
END DO OMEGDI1A.382
*ENDIF OMEGDI1A.383
*ENDIF GSM3F403.278
GSM3F403.279
*IF DEF,MPP GSM3F403.280
DO I=FIRST_FLD_PT+EW_Halo,LAST_U_FLD_PT-EW_Halo GSM3F403.281
*ELSE GSM3F403.282
DO I=1,U_FIELD OMEGDI1A.384
*ENDIF GSM3F403.283
WP(I) = WP(I) +WORK1(I) OMEGDI1A.385
END DO OMEGDI1A.386
OMEGDI1A.387
CL--------------------------------------------------------------------- OMEGDI1A.388
CL SECTION 5. CALCULATE OMEGA AS IN EQUATION (44). OMEGDI1A.389
CL--------------------------------------------------------------------- OMEGDI1A.390
OMEGDI1A.391
*IF DEF,MPP GSM3F403.284
DO I=FIRST_FLD_PT+EW_Halo,LAST_U_FLD_PT-EW_Halo GSM3F403.285
*ELSE GSM3F403.286
*IF DEF,GLOBAL OMEGDI1A.392
DO I=1,U_FIELD OMEGDI1A.393
*ELSE OMEGDI1A.394
DO I=2,U_FIELD-1 OMEGDI1A.395
*ENDIF OMEGDI1A.396
*ENDIF GSM3F403.287
OMEGA(I,K)= (WP(I)+WORK3(I))+((.5*(ETADOT(I,K+1)* GSM3F403.288
& (DELTA_AKH(K+1)+DELTA_BKH(K+1)*PSTAR_UV(I)) OMEGDI1A.398
& +ETADOT(I,K) OMEGDI1A.399
& *(DELTA_AKH(K)+DELTA_BKH(K)*PSTAR_UV(I)))) GSM3F403.289
& /(A*(DELTA_AK(K)+DELTA_BK(K)*PSTAR_UV(I)))) GSM3F403.290
END DO OMEGDI1A.402
*IF DEF,MPP GSM3F403.291
C Initialise unused rows of OMEGA GSM3F403.292
DO I=FIRST_FLD_PT+EW_Halo-1,1,-1 GSM3F403.293
OMEGA(I,K)=OMEGA(I+ROW_LENGTH,K) GSM3F403.294
ENDDO GSM3F403.295
DO I=LAST_U_FLD_PT-EW_Halo+1,U_FIELD GSM3F403.296
OMEGA(I,K)=OMEGA(I-ROW_LENGTH,K) GSM3F403.297
ENDDO GSM3F403.298
*ENDIF GSM3F403.299
*IF -DEF,GLOBAL OMEGDI1A.404
OMEGDI1A.405
CL LIMITED AREA MODEL SET VERTICAL VELOCITY ON BOUNDARY TO ZERO. OMEGDI1A.406
*IF -DEF,MPP GSM3F403.300
DO I=1,U_FIELD,ROW_LENGTH OMEGDI1A.408
OMEGA(I,K) = 0. OMEGDI1A.409
OMEGA(I+ROW_LENGTH-1,K) = 0. OMEGDI1A.410
OMEGA(I+ROW_LENGTH-2,K) = 0. OMEGDI1A.411
END DO OMEGDI1A.412
*ELSE GSM3F403.301
IF (at_left_of_LPG) THEN GSM3F403.302
DO I=ROW_LENGTH+1+EW_Halo,LAST_U_FLD_PT,ROW_LENGTH GSM3F403.303
OMEGA(I,K) = 0. GSM3F403.304
ENDDO GSM3F403.305
ENDIF GSM3F403.306
IF (at_right_of_LPG) THEN GSM3F403.307
DO I=ROW_LENGTH-EW_Halo,LAST_U_FLD_PT,ROW_LENGTH GSM3F403.308
OMEGA(I,K) = 0. GSM3F403.309
OMEGA(I-1,K) = 0. GSM3F403.310
ENDDO GSM3F403.311
ENDIF GSM3F403.312
*ENDIF GSM3F403.313
*ENDIF OMEGDI1A.414
OMEGDI1A.415
CL END LOOP OVER LEVELS. OMEGDI1A.416
300 CONTINUE OMEGDI1A.417
*IF DEF,MPP GSM3F403.314
CALL SWAPBOUNDS(OMEGA(1,1),ROW_LENGTH,tot_U_ROWS,EW_Halo, GSM3F403.315
& NS_Halo,P_LEVELS) GSM3F403.316
*ENDIF GSM3F403.317
OMEGDI1A.418
CL END OF ROUTINE OMEGA_DIAG OMEGDI1A.419
OMEGDI1A.420
RETURN OMEGDI1A.421
END OMEGDI1A.422
*ENDIF OMEGDI1A.423