source: palm/trunk/SOURCE/init_slope.f90 @ 623

Last change on this file since 623 was 623, checked in by raasch, 14 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 3.5 KB
Line 
1 SUBROUTINE init_slope
2
3!------------------------------------------------------------------------------!
4! Current revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: init_slope.f90 623 2010-12-10 08:52:17Z raasch $
11!
12! 622 2010-12-10 08:08:13Z raasch
13! optional barriers included in order to speed up collective operations
14!
15! Feb. 2007
16! RCS Log replace by Id keyword, revision history cleaned up
17!
18! Revision 1.5  2006/02/23 12:35:34  raasch
19! nanz_2dh renamed ngp_2dh
20!
21! Revision 1.1  2000/04/27 07:06:24  raasch
22! Initial revision
23!
24!
25! Description:
26! ------------
27! Initialization of the temperature field and other variables used in case
28! of a sloping surface.
29! Remember: when a sloping surface is used, only one constant temperature
30!           gradient is allowed!
31!------------------------------------------------------------------------------!
32
33    USE arrays_3d
34    USE constants
35    USE grid_variables
36    USE indices
37    USE pegrid
38    USE control_parameters
39
40    IMPLICIT NONE
41
42    INTEGER ::  i, j, k
43    REAL    ::  alpha, height, pt_value, radius
44    REAL, DIMENSION(:), ALLOCATABLE ::  pt_init_local
45
46!
47!-- Calculate reference temperature field needed for computing buoyancy
48    ALLOCATE( pt_slope_ref(nzb:nzt+1,nxl-1:nxr+1) )
49
50    DO  i = nxl-1, nxr+1
51       DO  k = nzb, nzt+1
52
53!
54!--       Compute height of grid-point relative to lower left corner of
55!--       the total domain.
56!--       First compute the distance between the actual grid point and the
57!--       lower left corner as well as the angle between the line connecting
58!--       these points and the bottom of the model.
59          IF ( k /= nzb )  THEN
60             radius = SQRT( ( i * dx )**2 + zu(k)**2 )
61             height = zu(k)
62          ELSE
63             radius = SQRT( ( i * dx )**2 )
64             height = 0.0
65          ENDIF
66          IF ( radius /= 0.0 )  THEN
67             alpha = ASIN( height / radius )
68          ELSE
69             alpha = 0.0
70          ENDIF
71!
72!--       Compute temperatures in the rotated coordinate system
73          alpha    = alpha + alpha_surface / 180.0 * pi
74          pt_value = pt_surface + radius * SIN( alpha ) * &
75                                  pt_vertical_gradient(1) / 100.0
76          pt_slope_ref(k,i) = pt_value
77       ENDDO               
78    ENDDO
79
80!
81!-- Temperature difference between left and right boundary of the total domain,
82!-- used for the cyclic boundary in x-direction
83    pt_slope_offset = (nx+1) * dx * sin_alpha_surface * &
84                      pt_vertical_gradient(1) / 100.0
85
86
87!
88!-- Following action must only be executed for initial runs
89    IF ( TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
90!
91!--    Set initial temperature equal to the reference temperature field
92       DO  j = nys-1, nyn+1
93          pt(:,j,:) = pt_slope_ref
94       ENDDO
95
96!
97!--    Recompute the mean initial temperature profile (mean along x-direction of
98!--    the rotated coordinate system)
99       ALLOCATE( pt_init_local(nzb:nzt+1) )
100       pt_init_local = 0.0
101       DO  i = nxl, nxr
102          DO  j =  nys, nyn
103             DO  k = nzb, nzt+1
104                pt_init_local(k) = pt_init_local(k) + pt(k,j,i)
105             ENDDO
106          ENDDO
107       ENDDO
108
109#if defined( __parallel )
110       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
111       CALL MPI_ALLREDUCE( pt_init_local, pt_init, nzt+2-nzb, MPI_REAL, &
112                            MPI_SUM, comm2d, ierr )
113#else
114       pt_init = pt_init_local
115#endif
116
117       pt_init = pt_init / ngp_2dh(0)
118       DEALLOCATE( pt_init_local )
119
120    ENDIF
121
122 END SUBROUTINE init_slope
Note: See TracBrowser for help on using the repository browser.