1 | ! |
---|
2 | ! $Id: friction_p.F 1910 2013-11-29 08:40:25Z idelkadi $ |
---|
3 | ! |
---|
4 | c======================================================================= |
---|
5 | SUBROUTINE friction_p(ucov,vcov,pdt) |
---|
6 | USE parallel_lmdz |
---|
7 | USE control_mod |
---|
8 | #ifdef CPP_IOIPSL |
---|
9 | USE IOIPSL |
---|
10 | #else |
---|
11 | ! if not using IOIPSL, we still need to use (a local version of) getin |
---|
12 | USE ioipsl_getincom |
---|
13 | #endif |
---|
14 | IMPLICIT NONE |
---|
15 | |
---|
16 | !======================================================================= |
---|
17 | ! |
---|
18 | ! Friction for the Newtonian case: |
---|
19 | ! -------------------------------- |
---|
20 | ! 2 possibilities (depending on flag 'friction_type' |
---|
21 | ! friction_type=0 : A friction that is only applied to the lowermost |
---|
22 | ! atmospheric layer |
---|
23 | ! friction_type=1 : Friction applied on all atmospheric layer (but |
---|
24 | ! (default) with stronger magnitude near the surface; see |
---|
25 | ! iniacademic.F) |
---|
26 | !======================================================================= |
---|
27 | |
---|
28 | #include "dimensions.h" |
---|
29 | #include "paramet.h" |
---|
30 | #include "comgeom2.h" |
---|
31 | #include "comconst.h" |
---|
32 | #include "iniprint.h" |
---|
33 | #include "academic.h" |
---|
34 | |
---|
35 | ! arguments: |
---|
36 | REAL,INTENT(inout) :: ucov( iip1,jjp1,llm ) |
---|
37 | REAL,INTENT(inout) :: vcov( iip1,jjm,llm ) |
---|
38 | REAL,INTENT(in) :: pdt ! time step |
---|
39 | |
---|
40 | ! local variables: |
---|
41 | REAL modv(iip1,jjp1),zco,zsi |
---|
42 | REAL vpn,vps,upoln,upols,vpols,vpoln |
---|
43 | REAL u2(iip1,jjp1),v2(iip1,jjm) |
---|
44 | INTEGER i,j,l |
---|
45 | REAL,PARAMETER :: cfric=1.e-5 |
---|
46 | LOGICAL,SAVE :: firstcall=.true. |
---|
47 | INTEGER,SAVE :: friction_type=1 |
---|
48 | CHARACTER(len=20) :: modname="friction_p" |
---|
49 | CHARACTER(len=80) :: abort_message |
---|
50 | !$OMP THREADPRIVATE(firstcall,friction_type) |
---|
51 | integer :: jjb,jje |
---|
52 | |
---|
53 | !$OMP SINGLE |
---|
54 | IF (firstcall) THEN |
---|
55 | ! set friction type |
---|
56 | call getin("friction_type",friction_type) |
---|
57 | if ((friction_type.lt.0).or.(friction_type.gt.1)) then |
---|
58 | abort_message="wrong friction type" |
---|
59 | write(lunout,*)'Friction: wrong friction type',friction_type |
---|
60 | call abort_gcm(modname,abort_message,42) |
---|
61 | endif |
---|
62 | firstcall=.false. |
---|
63 | ENDIF |
---|
64 | !$OMP END SINGLE COPYPRIVATE(friction_type,firstcall) |
---|
65 | |
---|
66 | if (friction_type.eq.0) then ! friction on first layer only |
---|
67 | !$OMP SINGLE |
---|
68 | c calcul des composantes au carre du vent naturel |
---|
69 | jjb=jj_begin |
---|
70 | jje=jj_end+1 |
---|
71 | if (pole_sud) jje=jj_end |
---|
72 | |
---|
73 | do j=jjb,jje |
---|
74 | do i=1,iip1 |
---|
75 | u2(i,j)=ucov(i,j,1)*ucov(i,j,1)*unscu2(i,j) |
---|
76 | enddo |
---|
77 | enddo |
---|
78 | |
---|
79 | jjb=jj_begin-1 |
---|
80 | jje=jj_end+1 |
---|
81 | if (pole_nord) jjb=jj_begin |
---|
82 | if (pole_sud) jje=jj_end-1 |
---|
83 | |
---|
84 | do j=jjb,jje |
---|
85 | do i=1,iip1 |
---|
86 | v2(i,j)=vcov(i,j,1)*vcov(i,j,1)*unscv2(i,j) |
---|
87 | enddo |
---|
88 | enddo |
---|
89 | |
---|
90 | c calcul du module de V en dehors des poles |
---|
91 | jjb=jj_begin |
---|
92 | jje=jj_end+1 |
---|
93 | if (pole_nord) jjb=jj_begin+1 |
---|
94 | if (pole_sud) jje=jj_end-1 |
---|
95 | |
---|
96 | do j=jjb,jje |
---|
97 | do i=2,iip1 |
---|
98 | modv(i,j)=sqrt(0.5*(u2(i-1,j)+u2(i,j)+v2(i,j-1)+v2(i,j))) |
---|
99 | enddo |
---|
100 | modv(1,j)=modv(iip1,j) |
---|
101 | enddo |
---|
102 | |
---|
103 | c les deux composantes du vent au pole sont obtenues comme |
---|
104 | c premiers modes de fourier de v pres du pole |
---|
105 | if (pole_nord) then |
---|
106 | |
---|
107 | upoln=0. |
---|
108 | vpoln=0. |
---|
109 | |
---|
110 | do i=2,iip1 |
---|
111 | zco=cos(rlonv(i))*(rlonu(i)-rlonu(i-1)) |
---|
112 | zsi=sin(rlonv(i))*(rlonu(i)-rlonu(i-1)) |
---|
113 | vpn=vcov(i,1,1)/cv(i,1) |
---|
114 | upoln=upoln+zco*vpn |
---|
115 | vpoln=vpoln+zsi*vpn |
---|
116 | enddo |
---|
117 | vpn=sqrt(upoln*upoln+vpoln*vpoln)/pi |
---|
118 | do i=1,iip1 |
---|
119 | c modv(i,1)=vpn |
---|
120 | modv(i,1)=modv(i,2) |
---|
121 | enddo |
---|
122 | |
---|
123 | endif |
---|
124 | |
---|
125 | if (pole_sud) then |
---|
126 | |
---|
127 | upols=0. |
---|
128 | vpols=0. |
---|
129 | do i=2,iip1 |
---|
130 | zco=cos(rlonv(i))*(rlonu(i)-rlonu(i-1)) |
---|
131 | zsi=sin(rlonv(i))*(rlonu(i)-rlonu(i-1)) |
---|
132 | vps=vcov(i,jjm,1)/cv(i,jjm) |
---|
133 | upols=upols+zco*vps |
---|
134 | vpols=vpols+zsi*vps |
---|
135 | enddo |
---|
136 | vps=sqrt(upols*upols+vpols*vpols)/pi |
---|
137 | do i=1,iip1 |
---|
138 | c modv(i,jjp1)=vps |
---|
139 | modv(i,jjp1)=modv(i,jjm) |
---|
140 | enddo |
---|
141 | |
---|
142 | endif |
---|
143 | |
---|
144 | c calcul du frottement au sol. |
---|
145 | |
---|
146 | jjb=jj_begin |
---|
147 | jje=jj_end |
---|
148 | if (pole_nord) jjb=jj_begin+1 |
---|
149 | if (pole_sud) jje=jj_end-1 |
---|
150 | |
---|
151 | do j=jjb,jje |
---|
152 | do i=1,iim |
---|
153 | ucov(i,j,1)=ucov(i,j,1) |
---|
154 | s -cfric*pdt*0.5*(modv(i+1,j)+modv(i,j))*ucov(i,j,1) |
---|
155 | enddo |
---|
156 | ucov(iip1,j,1)=ucov(1,j,1) |
---|
157 | enddo |
---|
158 | |
---|
159 | jjb=jj_begin |
---|
160 | jje=jj_end |
---|
161 | if (pole_sud) jje=jj_end-1 |
---|
162 | |
---|
163 | do j=jjb,jje |
---|
164 | do i=1,iip1 |
---|
165 | vcov(i,j,1)=vcov(i,j,1) |
---|
166 | s -cfric*pdt*0.5*(modv(i,j+1)+modv(i,j))*vcov(i,j,1) |
---|
167 | enddo |
---|
168 | vcov(iip1,j,1)=vcov(1,j,1) |
---|
169 | enddo |
---|
170 | !$OMP END SINGLE |
---|
171 | endif ! of if (friction_type.eq.0) |
---|
172 | |
---|
173 | if (friction_type.eq.1) then |
---|
174 | ! for ucov() |
---|
175 | jjb=jj_begin |
---|
176 | jje=jj_end |
---|
177 | if (pole_nord) jjb=jj_begin+1 |
---|
178 | if (pole_sud) jje=jj_end-1 |
---|
179 | |
---|
180 | !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) |
---|
181 | do l=1,llm |
---|
182 | ucov(1:iip1,jjb:jje,l)=ucov(1:iip1,jjb:jje,l)* |
---|
183 | & (1.-pdt*kfrict(l)) |
---|
184 | enddo |
---|
185 | !$OMP END DO NOWAIT |
---|
186 | |
---|
187 | ! for vcoc() |
---|
188 | jjb=jj_begin |
---|
189 | jje=jj_end |
---|
190 | if (pole_sud) jje=jj_end-1 |
---|
191 | |
---|
192 | !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) |
---|
193 | do l=1,llm |
---|
194 | vcov(1:iip1,jjb:jje,l)=vcov(1:iip1,jjb:jje,l)* |
---|
195 | & (1.-pdt*kfrict(l)) |
---|
196 | enddo |
---|
197 | !$OMP END DO |
---|
198 | endif ! of if (friction_type.eq.1) |
---|
199 | |
---|
200 | RETURN |
---|
201 | END |
---|
202 | |
---|