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

Last change on this file since 3 was 3, checked in by raasch, 17 years ago

RCS Log replace by Id keyword, revision history cleaned up

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