source: palm/trunk/SOURCE/disturb_field.f90 @ 2000

Last change on this file since 2000 was 2000, checked in by knoop, 5 years ago

Forced header and separation lines into 80 columns

  • Property svn:keywords set to Id
File size: 8.9 KB
RevLine 
[1682]1!> @file disturb_field.f90
[2000]2!------------------------------------------------------------------------------!
[1036]3! This file is part of PALM.
4!
[2000]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.
[1036]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!
[1818]17! Copyright 1997-2016 Leibniz Universitaet Hannover
[2000]18!------------------------------------------------------------------------------!
[1036]19!
[484]20! Current revisions:
[1318]21! ------------------
[2000]22! Forced header and separation lines into 80 columns
[1354]23!
[1321]24! Former revisions:
25! -----------------
26! $Id: disturb_field.f90 2000 2016-08-20 18:09:15Z knoop $
27!
[1683]28! 1682 2015-10-07 23:56:08Z knoop
29! Code annotations made doxygen readable
30!
[1426]31! 1425 2014-07-05 10:57:53Z knoop
32! bugfix: Parallel random number generator loop: print-statement no needed
33!
[1401]34! 1400 2014-05-09 14:03:54Z knoop
35! Parallel random number generator added
36!
[1354]37! 1353 2014-04-08 15:21:23Z heinze
38! REAL constants provided with KIND-attribute
39!
[1323]40! 1322 2014-03-20 16:38:49Z raasch
41! REAL constants defined as wp-kind
42!
[1321]43! 1320 2014-03-20 08:40:49Z raasch
[1320]44! ONLY-attribute added to USE-statements,
45! kind-parameters added to all INTEGER and REAL declaration statements,
46! kinds are defined in new module kinds,
47! revision history before 2012 removed,
48! comment fields (!:) to be used for variable explanations added to
49! all variable declaration statements
[1]50!
[1037]51! 1036 2012-10-22 13:43:42Z raasch
52! code put under GPL (PALM 3.9)
53!
[1]54! Revision 1.1  1998/02/04 15:40:45  raasch
55! Initial revision
56!
57!
58! Description:
59! ------------
[1682]60!> Imposing a random perturbation on a 3D-array.
61!> On parallel computers, the random number generator is as well called for all
62!> gridpoints of the total domain to ensure, regardless of the number of PEs
63!> used, that the elements of the array have the same values in the same
64!> order in every case. The perturbation range is steered by dist_range.
[1]65!------------------------------------------------------------------------------!
[1682]66 SUBROUTINE disturb_field( nzb_uv_inner, dist1, field )
67 
[1]68
[1320]69    USE control_parameters,   &
70        ONLY:  dist_nxl, dist_nxr, dist_nyn, dist_nys, dist_range,             &
71               disturbance_amplitude, disturbance_created,                     &
72               disturbance_level_ind_b, disturbance_level_ind_t, iran,         &
73               random_generator, topography
74               
75    USE cpulog,                                                                &
76        ONLY:  cpu_log, log_point
77       
78    USE indices,                                                               &
79        ONLY:  nbgp, nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt
80       
81    USE kinds
82   
83    USE random_function_mod,                                                   &
84        ONLY: random_function
[1400]85       
86    USE random_generator_parallel,                                             &
87        ONLY:  random_number_parallel, random_seed_parallel, random_dummy,     &
88               id_random_array, seq_random_array
[1]89
90    IMPLICIT NONE
91
[1682]92    INTEGER(iwp) ::  i  !<
93    INTEGER(iwp) ::  j  !<
94    INTEGER(iwp) ::  k  !<
[1320]95   
[1682]96    INTEGER(iwp) ::  nzb_uv_inner(nysg:nyng,nxlg:nxrg) !<
[1]97
[1682]98    REAL(wp) ::  randomnumber  !<
[1320]99   
[1682]100    REAL(wp) ::  dist1(nzb:nzt+1,nysg:nyng,nxlg:nxrg)  !<
101    REAL(wp) ::  field(nzb:nzt+1,nysg:nyng,nxlg:nxrg)  !<
[1320]102   
[1682]103    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  dist2  !<
[1]104
105
106    CALL cpu_log( log_point(20), 'disturb_field', 'start' )
107
108!
109!-- Create an additional temporary array and initialize the arrays needed
110!-- to store the disturbance
[667]111    ALLOCATE( dist2(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
[1353]112    dist1 = 0.0_wp
113    dist2 = 0.0_wp
[1]114
115!
116!-- Create the random perturbation and store it on temporary array
117    IF ( random_generator == 'numerical-recipes' )  THEN
118       DO  i = dist_nxl(dist_range), dist_nxr(dist_range)
119          DO  j = dist_nys(dist_range), dist_nyn(dist_range)
120             DO  k = disturbance_level_ind_b, disturbance_level_ind_t
[1322]121                randomnumber = 3.0_wp * disturbance_amplitude *                &
122                               ( random_function( iran ) - 0.5_wp )
[1320]123                IF ( nxl <= i  .AND.  nxr >= i  .AND.  nys <= j  .AND.         &
124                     nyn >= j )                                                &
[1]125                THEN
126                   dist1(k,j,i) = randomnumber
127                ENDIF
128             ENDDO
129          ENDDO
130       ENDDO
[1400]131    ELSEIF ( random_generator == 'random-parallel' )  THEN
132       DO  i = dist_nxl(dist_range), dist_nxr(dist_range)
133          DO  j = dist_nys(dist_range), dist_nyn(dist_range)
134             CALL random_seed_parallel( put=seq_random_array(:, j, i) )
135             DO  k = disturbance_level_ind_b, disturbance_level_ind_t
136                CALL random_number_parallel( random_dummy )
137                randomnumber = 3.0_wp * disturbance_amplitude *                &
138                               ( random_dummy - 0.5_wp )
139                IF ( nxl <= i  .AND.  nxr >= i  .AND.  nys <= j  .AND.         &
140                     nyn >= j )                                                &
141                THEN
142                   dist1(k,j,i) = randomnumber
143                ENDIF
144             ENDDO
145             CALL random_seed_parallel( get=seq_random_array(:, j, i) )
146          ENDDO
147       ENDDO
[1]148    ELSEIF ( random_generator == 'system-specific' )  THEN
149       DO  i = dist_nxl(dist_range), dist_nxr(dist_range)
150          DO  j = dist_nys(dist_range), dist_nyn(dist_range)
151             DO  k = disturbance_level_ind_b, disturbance_level_ind_t
152#if defined( __nec )
[1322]153                randomnumber = 3.0_wp * disturbance_amplitude *                &
154                               ( RANDOM( 0 ) - 0.5_wp )
[1]155#else
156                CALL RANDOM_NUMBER( randomnumber )
[1322]157                randomnumber = 3.0_wp * disturbance_amplitude *                &
158                                ( randomnumber - 0.5_wp )
[1]159#endif
[1320]160                IF ( nxl <= i .AND. nxr >= i .AND. nys <= j .AND. nyn >= j )   &
[1]161                THEN
162                   dist1(k,j,i) = randomnumber
163                ENDIF
164             ENDDO
165          ENDDO
166       ENDDO
167
168    ENDIF
169
170!
171!-- Exchange of ghost points for the random perturbation
[420]172
[667]173    CALL exchange_horiz( dist1, nbgp )
[1]174!
175!-- Applying the Shuman filter in order to smooth the perturbations.
176!-- Neighboured grid points in all three directions are used for the
177!-- filter operation.
[420]178!-- Loop has been splitted to make runs reproducible on HLRN systems using
179!-- compiler option -O3
180     DO  i = nxl, nxr
181        DO  j = nys, nyn
[1]182          DO  k = disturbance_level_ind_b-1, disturbance_level_ind_t+1
[1320]183             dist2(k,j,i) = ( dist1(k,j,i-1) + dist1(k,j,i+1)                  &
184                            + dist1(k,j+1,i) + dist1(k+1,j,i)                  &
[1322]185                            ) / 12.0_wp
[420]186          ENDDO
187          DO  k = disturbance_level_ind_b-1, disturbance_level_ind_t+1
188              dist2(k,j,i) = dist2(k,j,i) + ( dist1(k,j-1,i) + dist1(k-1,j,i)  &
[1322]189                            + 6.0_wp * dist1(k,j,i)                            &
190                            ) / 12.0_wp
[1]191          ENDDO
[420]192        ENDDO
193     ENDDO
[1]194
195!
196!-- Exchange of ghost points for the filtered perturbation.
197!-- Afterwards, filter operation and exchange of ghost points are repeated.
[667]198    CALL exchange_horiz( dist2, nbgp )
[420]199
[1]200    DO  i = nxl, nxr
201       DO  j = nys, nyn
202          DO  k = disturbance_level_ind_b-2, disturbance_level_ind_t+2
203             dist1(k,j,i) = ( dist2(k,j,i-1) + dist2(k,j,i+1) + dist2(k,j-1,i) &
204                            + dist2(k,j+1,i) + dist2(k+1,j,i) + dist2(k-1,j,i) &
[1322]205                            + 6.0_wp * dist2(k,j,i)                            &
206                            ) / 12.0_wp
[1]207          ENDDO
208       ENDDO
209    ENDDO
[420]210
[667]211    CALL exchange_horiz( dist1, nbgp )
[1]212
213!
214!-- Remove perturbations below topography (including one gridpoint above it
215!-- in order to allow for larger timesteps at the beginning of the simulation
216!-- (diffusion criterion))
217    IF ( TRIM( topography ) /= 'flat' )  THEN
[667]218       DO  i = nxlg, nxrg
219          DO  j = nysg, nyng
[1353]220             dist1(nzb:nzb_uv_inner(j,i)+1,j,i) = 0.0_wp
[1]221          ENDDO
222       ENDDO
223    ENDIF
224
225!
226!-- Random perturbation is added to the array to be disturbed.
[667]227    DO  i = nxlg, nxrg
228       DO  j = nysg, nyng
[1]229          DO  k = disturbance_level_ind_b-2, disturbance_level_ind_t+2
230             field(k,j,i) = field(k,j,i) + dist1(k,j,i)
231          ENDDO
232       ENDDO
233    ENDDO
234
235!
236!-- Deallocate the temporary array
237    DEALLOCATE( dist2 )
238
239!
240!-- Set a flag, which indicates that a random perturbation is imposed
241    disturbance_created = .TRUE.
242
243
244    CALL cpu_log( log_point(20), 'disturb_field', 'stop' )
245
246
247 END SUBROUTINE disturb_field
Note: See TracBrowser for help on using the repository browser.