1 |
molod |
1.2 |
subroutine plume2dyn(qplume,Nxplume,Lmplume,uref,vref,flag, |
2 |
|
|
. idim1,idim2,jdim1,jdim2,Lmout,Nsx,Nsy,bi,bj,qdyn1,qdyn2) |
3 |
molod |
1.1 |
C*********************************************************************** |
4 |
|
|
C Purpose: |
5 |
molod |
1.2 |
C To interpolate an arbitrary quantity from higher resolution plume |
6 |
molod |
1.1 |
C grid to the model's dynamics grid |
7 |
|
|
C Algorithm: |
8 |
molod |
1.2 |
C Plumes -> Dynamics computes the plumes mean value, and in the case |
9 |
|
|
C of a vector field, preserves the direction of a vector |
10 |
|
|
C given in (uref,vref) |
11 |
molod |
1.1 |
C |
12 |
|
|
C Input: |
13 |
molod |
1.2 |
C qplume... [idim2,jdim2,im,Lmplume,bi] Quantity on Input Grid |
14 |
molod |
1.1 |
C idimin... Longitude Dimension of Input |
15 |
|
|
C Lmplume.. Vertical Dimension of Input |
16 |
molod |
1.2 |
C uref .... [im,jm,Lmout,bi,bj] Reference u-component of velocity |
17 |
|
|
C vref .... [im,jm,Lmout,bi,bj] Reference v-component of velocity |
18 |
|
|
C flag .... Flag to indicate vector (1) or scalar (0) interpolation |
19 |
|
|
C idim1,2.. Beginning and ending i-values of output grid |
20 |
|
|
C jdim1,2.. Beginning and ending j-values of output grid |
21 |
|
|
C Lmout.... Vertical Dimension of Output |
22 |
molod |
1.1 |
C Nsx...... Number of processes in x-direction |
23 |
|
|
C Nsy...... Number of processes in y-direction |
24 |
|
|
C bi....... Index of process number in x-direction |
25 |
|
|
C bj....... Index of process number in x-direction |
26 |
|
|
C |
27 |
|
|
C Output: |
28 |
molod |
1.2 |
C qdyn1..... [im,jm,Lmout,bi,bj] Field at output grid (dynamics) |
29 |
|
|
C qdyn2..... [im,jm,Lmout,bi,bj] Field at output grid (dynamics) |
30 |
molod |
1.1 |
C |
31 |
|
|
C Notes: |
32 |
molod |
1.2 |
C 1) Assume (for now) that the number of vertical levels is the |
33 |
|
|
C same on both the input and output grids |
34 |
molod |
1.1 |
C*********************************************************************** |
35 |
|
|
implicit none |
36 |
|
|
#include "CPP_OPTIONS.h" |
37 |
|
|
|
38 |
molod |
1.2 |
integer Nxplume, Lmplume, Lmout, Nsx, Nsy |
39 |
|
|
integer idim1, idim2, jdim1, jdim2, bi, bj, flag |
40 |
|
|
_RL qplume(idim2,jdim2,Nxplume,Lmplume,Nsx) |
41 |
|
|
_RL uref(idim1:idim2,jdim1:jdim2,Lmout,Nsx,Nsy) |
42 |
|
|
_RL vref(idim1:idim2,jdim1:jdim2,Lmout,Nsx,Nsy) |
43 |
|
|
_RL qdyn1(idim1:idim2,jdim1:jdim2,Lmout,Nsx,Nsy) |
44 |
|
|
_RL qdyn2(idim1:idim2,jdim1:jdim2,Lmout,Nsx,Nsy) |
45 |
molod |
1.1 |
|
46 |
molod |
1.2 |
integer i,j,L,iplume |
47 |
|
|
_RL qplumeav(idim1,jdim2,Lmplume) |
48 |
|
|
_RL sqrtarg |
49 |
|
|
|
50 |
|
|
C First step - compute the average of qplume over Nxplume |
51 |
|
|
do j = jdim1,jdim2 |
52 |
|
|
do i = idim1,idim2 |
53 |
|
|
do L = 1,Lmplume |
54 |
|
|
qplumeav(i,j,L) = 0. |
55 |
|
|
do iplume = 1,Nxplume |
56 |
|
|
qplumeav(i,j,L)=qplumeav(i,j,L)+qplume(i,j,iplume,L,bi)/Nxplume |
57 |
|
|
enddo |
58 |
|
|
enddo |
59 |
|
|
enddo |
60 |
|
|
enddo |
61 |
molod |
1.1 |
|
62 |
molod |
1.2 |
C Now check the flag -- if a scalar, we are done - just assign |
63 |
|
|
C the average to all the i and j points of the output grid. |
64 |
|
|
C If a vector, there is some more work to do in order to preserve |
65 |
|
|
C the angle given by uref and vref |
66 |
molod |
1.1 |
|
67 |
molod |
1.2 |
if (flag.eq.0) then |
68 |
molod |
1.1 |
do j = jdim1,jdim2 |
69 |
molod |
1.2 |
do i = idim1,idim2 |
70 |
|
|
do L = 1,Lmplume |
71 |
|
|
qdyn1(i,j,L,bi,bj) = qplumeav(i,j,L) |
72 |
|
|
enddo |
73 |
|
|
enddo |
74 |
|
|
enddo |
75 |
|
|
elseif (flag.eq.1) then |
76 |
|
|
do j = jdim1,jdim2 |
77 |
|
|
do i = idim1,idim2 |
78 |
|
|
do L = 1,Lmplume |
79 |
molod |
1.3 |
if(vref(i,j,L,bi,bj).ne.0.) then |
80 |
|
|
sqrtarg = (qplumeav(i,j,L)*qplumeav(i,j,L)) / |
81 |
molod |
1.2 |
. ( ( (uref(i,j,L,bi,bj)*uref(i,j,L,bi,bj)) / |
82 |
|
|
. (vref(i,j,L,bi,bj)*vref(i,j,L,bi,bj)) ) + 1. ) |
83 |
molod |
1.3 |
qdyn2(i,j,L,bi,bj) = sqrt(sqrtarg) |
84 |
|
|
qdyn1(i,j,L,bi,bj) = qdyn2(i,j,L,bi,bj) * |
85 |
molod |
1.2 |
. (uref(i,j,L,bi,bj)/vref(i,j,L,bi,bj)) |
86 |
molod |
1.3 |
else |
87 |
|
|
qdyn1(i,j,L,bi,bj) = qplumeav(i,j,L) |
88 |
|
|
qdyn2(i,j,L,bi,bj) = 0. |
89 |
molod |
1.2 |
enddo |
90 |
|
|
enddo |
91 |
molod |
1.1 |
enddo |
92 |
molod |
1.2 |
endif |
93 |
molod |
1.1 |
|
94 |
|
|
return |
95 |
|
|
end |