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

Last change on this file since 3948 was 3725, checked in by raasch, 6 years ago

modifications to avoid compiler warnings about unused variables, temperton-fft: GOTO statements replaced, file re-formatted corresponding to coding standards, ssh-calls for compilations on remote systems modified to avoid output of login messages on specific systems changed again (palmbuild, reverted as before r3549), error messages for failed restarts extended (palmrun)

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