source: palm/trunk/SOURCE/buoyancy.f90 @ 3253

Last change on this file since 3253 was 3241, checked in by raasch, 6 years ago

various changes to avoid compiler warnings (mainly removal of unused variables)

  • Property svn:keywords set to Id
File size: 11.4 KB
Line 
1!> @file buoyancy.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2018 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: buoyancy.f90 3241 2018-09-12 15:02:00Z suehring $
27! unused Variable removed
28!
29! 2718 2018-01-02 08:49:38Z maronga
30! Corrected "Former revisions" section
31!
32! 2696 2017-12-14 17:12:51Z kanani
33! Change in file header (GPL part)
34!
35! 2233 2017-05-30 18:08:54Z suehring
36!
37! 2232 2017-05-30 17:47:52Z suehring
38! Adjustments to new topography concept
39!
40! 2118 2017-01-17 16:38:49Z raasch
41! OpenACC version of subroutine removed
42!
43! 2000 2016-08-20 18:09:15Z knoop
44! Forced header and separation lines into 80 columns
45!
46! 1873 2016-04-18 14:50:06Z maronga
47! Module renamed (removed _mod)
48!
49!
50! 1850 2016-04-08 13:29:27Z maronga
51! Module renamed
52!
53!
54! 1682 2015-10-07 23:56:08Z knoop
55! Code annotations made doxygen readable
56!
57! 1374 2014-04-25 12:55:07Z raasch
58! missing variables added to ONLY list
59!
60! 1365 2014-04-22 15:03:56Z boeske
61! Calculation of reference state in subroutine calc_mean_profile moved to
62! subroutine time_integration,
63! subroutine calc_mean_profile moved to new file calc_mean_profile.f90
64!
65! 1353 2014-04-08 15:21:23Z heinze
66! REAL constants provided with KIND-attribute
67!
68! 1320 2014-03-20 08:40:49Z raasch
69! ONLY-attribute added to USE-statements,
70! kind-parameters added to all INTEGER and REAL declaration statements,
71! kinds are defined in new module kinds,
72! revision history before 2012 removed,
73! comment fields (!:) to be used for variable explanations added to
74! all variable declaration statements
75!
76! 1257 2013-11-08 15:18:40Z raasch
77! vector length (32) removed from openacc clause
78!
79! 1241 2013-10-30 11:36:58Z heinze
80! Generalize calc_mean_profile for wider use: use additional steering
81! character loc
82!
83! 1179 2013-06-14 05:57:58Z raasch
84! steering of reference state revised (var_reference and pr removed from
85! parameter list), use_reference renamed use_single_reference_value
86!
87! 1171 2013-05-30 11:27:45Z raasch
88! openacc statements added to use_reference-case in accelerator version
89!
90! 1153 2013-05-10 14:33:08Z raasch
91! code adjustments of accelerator version required by PGI 12.3 / CUDA 5.0
92!
93! 1128 2013-04-12 06:19:32Z raasch
94! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
95! j_north
96!
97! 1036 2012-10-22 13:43:42Z raasch
98! code put under GPL (PALM 3.9)
99!
100! 1015 2012-09-27 09:23:24Z raasch
101! accelerator version (*_acc) added
102!
103! 1010 2012-09-20 07:59:54Z raasch
104! cpp switch __nopointer added for pointer free version
105!
106! Revision 1.1  1997/08/29 08:56:48  raasch
107! Initial revision
108!
109!
110! Description:
111! ------------
112!> Buoyancy term of the third component of the equation of motion.
113!> @attention Humidity is not regarded when using a sloping surface!
114!------------------------------------------------------------------------------!
115 MODULE buoyancy_mod
116 
117
118    PRIVATE
119    PUBLIC buoyancy
120
121    INTERFACE buoyancy
122       MODULE PROCEDURE buoyancy
123       MODULE PROCEDURE buoyancy_ij
124    END INTERFACE buoyancy
125
126 CONTAINS
127
128
129!------------------------------------------------------------------------------!
130! Description:
131! ------------
132!> Call for all grid points
133!------------------------------------------------------------------------------!
134    SUBROUTINE buoyancy( var, wind_component )
135
136       USE arrays_3d,                                                          &
137           ONLY:  pt, pt_slope_ref, ref_state, tend
138
139       USE control_parameters,                                                 &
140           ONLY:  atmos_ocean_sign, cos_alpha_surface, g, message_string,      &
141                  pt_surface, sin_alpha_surface, sloping_surface
142
143       USE indices,                                                            &
144           ONLY:  nxl, nxlg, nxlu, nxr, nxrg, nyn, nyng, nys, nysg, nzb,       &
145                  nzt, wall_flags_0
146
147       USE kinds
148
149       USE pegrid
150
151
152       IMPLICIT NONE
153
154       INTEGER(iwp) ::  i              !<
155       INTEGER(iwp) ::  j              !<
156       INTEGER(iwp) ::  k              !<
157       INTEGER(iwp) ::  wind_component !<
158       
159#if defined( __nopointer )
160       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  var !<
161#else
162       REAL(wp), DIMENSION(:,:,:), POINTER ::  var
163#endif
164
165
166       IF ( .NOT. sloping_surface )  THEN
167!
168!--       Normal case: horizontal surface
169          DO  i = nxl, nxr
170             DO  j = nys, nyn
171                DO  k = nzb+1, nzt-1
172                   tend(k,j,i) = tend(k,j,i) + atmos_ocean_sign * g * 0.5_wp *  &
173                          (                                                     &
174                             ( var(k,j,i)   - ref_state(k) )   / ref_state(k) + &
175                             ( var(k+1,j,i) - ref_state(k+1) ) / ref_state(k+1) &
176                          ) * MERGE( 1.0_wp, 0.0_wp,                            &
177                                     BTEST( wall_flags_0(k,j,i), 0 ) )
178                ENDDO
179             ENDDO
180          ENDDO
181
182       ELSE
183!
184!--       Buoyancy term for a surface with a slope in x-direction. The equations
185!--       for both the u and w velocity-component contain proportionate terms.
186!--       Temperature field at time t=0 serves as environmental temperature.
187!--       Reference temperature (pt_surface) is the one at the lower left corner
188!--       of the total domain.
189          IF ( wind_component == 1 )  THEN
190
191             DO  i = nxlu, nxr
192                DO  j = nys, nyn
193                   DO  k = nzb+1, nzt-1
194                      tend(k,j,i) = tend(k,j,i) + g * sin_alpha_surface *         &
195                           0.5_wp * ( ( pt(k,j,i-1)         + pt(k,j,i)         ) &
196                                    - ( pt_slope_ref(k,i-1) + pt_slope_ref(k,i) ) &
197                                    ) / pt_surface                                &
198                                      * MERGE( 1.0_wp, 0.0_wp,                    &
199                                               BTEST( wall_flags_0(k,j,i), 0 ) )
200                   ENDDO
201                ENDDO
202             ENDDO
203
204          ELSEIF ( wind_component == 3 )  THEN
205
206             DO  i = nxl, nxr
207                DO  j = nys, nyn
208                   DO  k = nzb+1, nzt-1
209                      tend(k,j,i) = tend(k,j,i) + g * cos_alpha_surface *         &
210                           0.5_wp * ( ( pt(k,j,i)         + pt(k+1,j,i)         ) &
211                                    - ( pt_slope_ref(k,i) + pt_slope_ref(k+1,i) ) &
212                                    ) / pt_surface                                &
213                                      * MERGE( 1.0_wp, 0.0_wp,                    &
214                                               BTEST( wall_flags_0(k,j,i), 0 ) )
215                   ENDDO
216                ENDDO
217            ENDDO
218
219          ELSE
220             
221             WRITE( message_string, * ) 'no term for component "',             &
222                                       wind_component,'"'
223             CALL message( 'buoyancy', 'PA0159', 1, 2, 0, 6, 0 )
224
225          ENDIF
226
227       ENDIF
228
229    END SUBROUTINE buoyancy
230
231
232!------------------------------------------------------------------------------!
233! Description:
234! ------------
235!> Call for grid point i,j
236!> @attention PGI-compiler creates SIGFPE if opt>1 is used! Therefore, opt=1 is
237!>            forced by compiler-directive.
238!------------------------------------------------------------------------------!
239!pgi$r opt=1
240    SUBROUTINE buoyancy_ij( i, j, var, wind_component )
241
242       USE arrays_3d,                                                          &
243           ONLY:  pt, pt_slope_ref, ref_state, tend
244
245       USE control_parameters,                                                 &
246           ONLY:  atmos_ocean_sign, cos_alpha_surface, g, message_string,      &
247                  pt_surface, sin_alpha_surface, sloping_surface
248
249       USE indices,                                                            &
250           ONLY:  nxlg, nxrg, nyng, nysg, nzb, nzt, wall_flags_0
251
252       USE kinds
253
254       USE pegrid
255
256
257       IMPLICIT NONE
258
259       INTEGER(iwp) ::  i              !<
260       INTEGER(iwp) ::  j              !<
261       INTEGER(iwp) ::  k              !<
262       INTEGER(iwp) ::  wind_component !<
263       
264#if defined( __nopointer )
265       REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ::  var !<
266#else
267       REAL(wp), DIMENSION(:,:,:), POINTER ::  var
268#endif
269
270
271       IF ( .NOT. sloping_surface )  THEN
272!
273!--       Normal case: horizontal surface
274          DO  k = nzb+1, nzt-1
275              tend(k,j,i) = tend(k,j,i) + atmos_ocean_sign * g * 0.5_wp * (    &
276                        ( var(k,j,i)   - ref_state(k)   ) / ref_state(k)   +   &
277                        ( var(k+1,j,i) - ref_state(k+1) ) / ref_state(k+1)     &
278                                                                          )    &
279                                      * MERGE( 1.0_wp, 0.0_wp,                 &
280                                               BTEST( wall_flags_0(k,j,i), 0 ) )
281          ENDDO
282
283       ELSE
284!
285!--       Buoyancy term for a surface with a slope in x-direction. The equations
286!--       for both the u and w velocity-component contain proportionate terms.
287!--       Temperature field at time t=0 serves as environmental temperature.
288!--       Reference temperature (pt_surface) is the one at the lower left corner
289!--       of the total domain.
290          IF ( wind_component == 1 )  THEN
291
292             DO  k = nzb+1, nzt-1
293                tend(k,j,i) = tend(k,j,i) + g * sin_alpha_surface *               &
294                           0.5_wp * ( ( pt(k,j,i-1)         + pt(k,j,i)         ) &
295                                    - ( pt_slope_ref(k,i-1) + pt_slope_ref(k,i) ) &
296                                    ) / pt_surface                                &
297                                      * MERGE( 1.0_wp, 0.0_wp,                    &
298                                               BTEST( wall_flags_0(k,j,i), 0 ) )
299             ENDDO
300
301          ELSEIF ( wind_component == 3 )  THEN
302
303             DO  k = nzb+1, nzt-1
304                tend(k,j,i) = tend(k,j,i) + g * cos_alpha_surface *               &
305                           0.5_wp * ( ( pt(k,j,i)         + pt(k+1,j,i)         ) &
306                                    - ( pt_slope_ref(k,i) + pt_slope_ref(k+1,i) ) &
307                                    ) / pt_surface                                &
308                                      * MERGE( 1.0_wp, 0.0_wp,                    &
309                                               BTEST( wall_flags_0(k,j,i), 0 ) )
310             ENDDO
311
312          ELSE
313
314             WRITE( message_string, * ) 'no term for component "',             &
315                                       wind_component,'"'
316             CALL message( 'buoyancy', 'PA0159', 1, 2, 0, 6, 0 )
317
318          ENDIF
319
320       ENDIF
321
322    END SUBROUTINE buoyancy_ij
323
324 END MODULE buoyancy_mod
Note: See TracBrowser for help on using the repository browser.