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

Last change on this file since 791 was 668, checked in by suehring, 14 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 5.8 KB
Line 
1 SUBROUTINE disturb_field( nzb_uv_inner, dist1, field )
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6!
7! Former revisions:
8! -----------------
9! $Id: disturb_field.f90 668 2010-12-23 13:22:58Z raasch $
10!
11! 667 2010-12-23 12:06:00Z suehring/gryschka
12! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
13! Calls of exchange_horiz are modified.
14!
15! 420 2010-01-13 15:10:53Z franke
16! Loop was split to make runs reproducible when using ifort compiler
17!
18! 75 2007-03-22 09:54:05Z raasch
19! xrp, ynp eliminated, 2nd+3rd argument removed from exchange horiz
20!
21! RCS Log replace by Id keyword, revision history cleaned up
22!
23! Revision 1.11  2006/08/04 14:31:59  raasch
24! izuf renamed iran
25!
26! Revision 1.1  1998/02/04 15:40:45  raasch
27! Initial revision
28!
29!
30! Description:
31! ------------
32! Imposing a random perturbation on a 3D-array.
33! On parallel computers, the random number generator is as well called for all
34! gridpoints of the total domain to ensure, regardless of the number of PEs
35! used, that the elements of the array have the same values in the same
36! order in every case. The perturbation range is steered by dist_range.
37!------------------------------------------------------------------------------!
38
39    USE control_parameters
40    USE cpulog
41    USE grid_variables
42    USE indices
43    USE interfaces
44    USE random_function_mod
45
46    IMPLICIT NONE
47
48    INTEGER ::  i, j, k
49    INTEGER ::  nzb_uv_inner(nysg:nyng,nxlg:nxrg)
50
51    REAL    ::  randomnumber,                             &
52                dist1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &
53                field(nzb:nzt+1,nysg:nyng,nxlg:nxrg)
54    REAL, DIMENSION(:,:,:), ALLOCATABLE ::  dist2
55
56
57    CALL cpu_log( log_point(20), 'disturb_field', 'start' )
58
59!
60!-- Create an additional temporary array and initialize the arrays needed
61!-- to store the disturbance
62    ALLOCATE( dist2(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
63    dist1 = 0.0
64    dist2 = 0.0
65
66!
67!-- Create the random perturbation and store it on temporary array
68    IF ( random_generator == 'numerical-recipes' )  THEN
69       DO  i = dist_nxl(dist_range), dist_nxr(dist_range)
70          DO  j = dist_nys(dist_range), dist_nyn(dist_range)
71             DO  k = disturbance_level_ind_b, disturbance_level_ind_t
72                randomnumber = 3.0 * disturbance_amplitude * &
73                               ( random_function( iran ) - 0.5 )
74                IF ( nxl <= i  .AND.  nxr >= i  .AND.  nys <= j  .AND.  &
75                     nyn >= j ) &
76                THEN
77                   dist1(k,j,i) = randomnumber
78                ENDIF
79             ENDDO
80          ENDDO
81       ENDDO
82    ELSEIF ( random_generator == 'system-specific' )  THEN
83       DO  i = dist_nxl(dist_range), dist_nxr(dist_range)
84          DO  j = dist_nys(dist_range), dist_nyn(dist_range)
85             DO  k = disturbance_level_ind_b, disturbance_level_ind_t
86#if defined( __nec )
87                randomnumber = 3.0 * disturbance_amplitude * &
88                               ( RANDOM( 0 ) - 0.5 )
89#else
90                CALL RANDOM_NUMBER( randomnumber )
91                randomnumber = 3.0 * disturbance_amplitude * &
92                                ( randomnumber - 0.5 )
93#endif
94                IF ( nxl <= i .AND. nxr >= i .AND. nys <= j .AND. nyn >= j ) &
95                THEN
96                   dist1(k,j,i) = randomnumber
97                ENDIF
98             ENDDO
99          ENDDO
100       ENDDO
101
102    ENDIF
103
104!
105!-- Exchange of ghost points for the random perturbation
106
107    CALL exchange_horiz( dist1, nbgp )
108!
109!-- Applying the Shuman filter in order to smooth the perturbations.
110!-- Neighboured grid points in all three directions are used for the
111!-- filter operation.
112!-- Loop has been splitted to make runs reproducible on HLRN systems using
113!-- compiler option -O3
114     DO  i = nxl, nxr
115        DO  j = nys, nyn
116          DO  k = disturbance_level_ind_b-1, disturbance_level_ind_t+1
117             dist2(k,j,i) = ( dist1(k,j,i-1) + dist1(k,j,i+1) &
118                            + dist1(k,j+1,i) + dist1(k+1,j,i) &
119                            ) / 12.0
120          ENDDO
121          DO  k = disturbance_level_ind_b-1, disturbance_level_ind_t+1
122              dist2(k,j,i) = dist2(k,j,i) + ( dist1(k,j-1,i) + dist1(k-1,j,i)  &
123                            + 6.0 * dist1(k,j,i)                               &
124                            ) / 12.0
125          ENDDO
126        ENDDO
127     ENDDO
128
129!
130!-- Exchange of ghost points for the filtered perturbation.
131!-- Afterwards, filter operation and exchange of ghost points are repeated.
132    CALL exchange_horiz( dist2, nbgp )
133
134    DO  i = nxl, nxr
135       DO  j = nys, nyn
136          DO  k = disturbance_level_ind_b-2, disturbance_level_ind_t+2
137             dist1(k,j,i) = ( dist2(k,j,i-1) + dist2(k,j,i+1) + dist2(k,j-1,i) &
138                            + dist2(k,j+1,i) + dist2(k+1,j,i) + dist2(k-1,j,i) &
139                            + 6.0 * dist2(k,j,i)                               &
140                            ) / 12.0
141          ENDDO
142       ENDDO
143    ENDDO
144
145    CALL exchange_horiz( dist1, nbgp )
146
147!
148!-- Remove perturbations below topography (including one gridpoint above it
149!-- in order to allow for larger timesteps at the beginning of the simulation
150!-- (diffusion criterion))
151    IF ( TRIM( topography ) /= 'flat' )  THEN
152       DO  i = nxlg, nxrg
153          DO  j = nysg, nyng
154             dist1(nzb:nzb_uv_inner(j,i)+1,j,i) = 0.0
155          ENDDO
156       ENDDO
157    ENDIF
158
159!
160!-- Random perturbation is added to the array to be disturbed.
161    DO  i = nxlg, nxrg
162       DO  j = nysg, nyng
163          DO  k = disturbance_level_ind_b-2, disturbance_level_ind_t+2
164             field(k,j,i) = field(k,j,i) + dist1(k,j,i)
165          ENDDO
166       ENDDO
167    ENDDO
168
169!
170!-- Deallocate the temporary array
171    DEALLOCATE( dist2 )
172
173!
174!-- Set a flag, which indicates that a random perturbation is imposed
175    disturbance_created = .TRUE.
176
177
178    CALL cpu_log( log_point(20), 'disturb_field', 'stop' )
179
180
181 END SUBROUTINE disturb_field
Note: See TracBrowser for help on using the repository browser.