Home


About us
PALM
PALM-4U
Installation
Release notes
Documentation
Getting Help
Gallery
Developers area

source: palm/trunk/SOURCE/diffusion_s.f90 @ 161

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

New:
---

Plant canopy model of Watanabe (2004,BLM 112,307-341) added.
It can be switched on by the inipar parameter plant_canopy.
The inipar parameter canopy_mode can be used to prescribe a
plant canopy type. The default case is a homogeneous plant
canopy. Heterogeneous distributions of the leaf area
density and the canopy drag coefficient can be defined in the
new routine user_init_plant_canopy (user_interface).
The inipar parameters lad_surface, lad_vertical_gradient and
lad_vertical_gradient_level can be used in order to
prescribe the vertical profile of leaf area density. The
inipar parameter drag_coefficient determines the canopy
drag coefficient.
Finally, the inipar parameter pch_index determines the
index of the upper boundary of the plant canopy.

Allow new case bc_uv_t = 'dirichlet_0' for channel flow.

For unknown variables (CASE DEFAULT) call new subroutine user_data_output_dvrp

Pressure boundary conditions for vertical walls added to the multigrid solver.
They are applied using new wall flag arrays (wall_flags_..) which are defined
for each grid level. New argument gls added to routine user_init_grid
(user_interface).

Frequence of sorting particles can be controlled with new particles_par
parameter dt_sort_particles. Sorting is moved from the SGS timestep loop in
advec_particles after the end of this loop.

advec_particles, check_parameters, data_output_dvrp, header, init_3d_model, init_grid, init_particles, init_pegrid, modules, package_parin, parin, plant_canopy_model, read_var_list, read_3d_binary, user_interface, write_var_list, write_3d_binary

Changed:

Redefine initial nzb_local as the actual total size of topography (later the
extent of topography in nzb_local is reduced by 1dx at the E topography walls
and by 1dy at the N topography walls to form the basis for nzb_s_inner);
for consistency redefine 'single_building' case.

Vertical profiles now based on nzb_s_inner; they are divided by
ngp_2dh_s_inner (scalars, procucts of scalars) and ngp_2dh (staggered velocity
components and their products, procucts of scalars and velocity components),
respectively.

Allow two instead of one digit to specify isosurface and slicer variables.

Status of 3D-volume NetCDF data file only depends on switch netcdf_64bit_3d (check_open)

prognostic_equations include the respective wall_*flux in the parameter list of
calls of diffusion_s. Same as before, only the values of wall_heatflux(0:4)
can be assigned. At present, wall_humidityflux, wall_qflux, wall_salinityflux,
and wall_scalarflux are kept zero. diffusion_s uses the respective wall_*flux
instead of wall_heatflux. This update serves two purposes:

  • it avoids errors in calculations with humidity/scalar/salinity and prescribed

non-zero wall_heatflux,

  • it prepares PALM for a possible assignment of wall fluxes of

humidity/scalar/salinity in a future release.

buoyancy, check_open, data_output_dvrp, diffusion_s, diffusivities, flow_statistics, header, init_3d_model, init_dvrp, init_grid, modules, prognostic_equations

Errors:

Bugfix: summation of sums_l_l in diffusivities.

Several bugfixes in the ocean part: Initial density rho is calculated
(init_ocean). Error in initializing u_init and v_init removed
(check_parameters). Calculation of density flux now starts from
nzb+1 (production_e).

Bugfix: pleft/pright changed to pnorth/psouth in sendrecv of particle tail
numbers along y, small bugfixes in the SGS part (advec_particles)

Bugfix: model_string needed a default value (combine_plot_fields)

Bugfix: wavenumber calculation for even nx in routines maketri (poisfft)

Bugfix: assignment of fluxes at walls

Bugfix: absolute value of f must be used when calculating the Blackadar mixing length (init_1d_model)

advec_particles, check_parameters, combine_plot_fields, diffusion_s, diffusivities, init_ocean, init_1d_model, poisfft, production_e
  • Property svn:keywords set to Id
File size: 10.8 KB
Line 
1 MODULE diffusion_s_mod
2
3!------------------------------------------------------------------------------!
4! Actual revisions:
5! -----------------
6!
7!
8! Former revisions:
9! -----------------
10! $Id: diffusion_s.f90 139 2007-11-29 09:37:41Z letzel $
11!
12! 129 2007-10-30 12:12:24Z letzel
13! replace wall_heatflux by wall_s_flux that is now included in the parameter
14! list, bugfix for assignment of fluxes at walls
15!
16! 20 2007-02-26 00:12:32Z raasch
17! Bugfix: ddzw dimensioned 1:nzt"+1"
18! Calculation extended for gridpoint nzt, fluxes can be given at top,
19! +s_flux_t in parameter list, s_flux renamed s_flux_b
20!
21! RCS Log replace by Id keyword, revision history cleaned up
22!
23! Revision 1.8  2006/02/23 10:34:17  raasch
24! nzb_2d replaced by nzb_s_outer in horizontal diffusion and by nzb_s_inner
25! or nzb_diff_s_inner, respectively, in vertical diffusion, prescribed surface
26! fluxes at vertically oriented topography
27!
28! Revision 1.1  2000/04/13 14:54:02  schroeter
29! Initial revision
30!
31!
32! Description:
33! ------------
34! Diffusion term of scalar quantities (temperature and water content)
35!------------------------------------------------------------------------------!
36
37    PRIVATE
38    PUBLIC diffusion_s
39
40    INTERFACE diffusion_s
41       MODULE PROCEDURE diffusion_s
42       MODULE PROCEDURE diffusion_s_ij
43    END INTERFACE diffusion_s
44
45 CONTAINS
46
47
48!------------------------------------------------------------------------------!
49! Call for all grid points
50!------------------------------------------------------------------------------!
51    SUBROUTINE diffusion_s( ddzu, ddzw, kh, s, s_flux_b, s_flux_t, &
52                            wall_s_flux, tend )
53
54       USE control_parameters
55       USE grid_variables
56       USE indices
57
58       IMPLICIT NONE
59
60       INTEGER ::  i, j, k
61       REAL    ::  vertical_gridspace
62       REAL    ::  ddzu(1:nzt+1), ddzw(1:nzt+1)
63       REAL    ::  tend(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1)
64       REAL    ::  wall_s_flux(0:4)
65       REAL, DIMENSION(:,:),   POINTER ::  s_flux_b, s_flux_t
66       REAL, DIMENSION(:,:,:), POINTER ::  kh, s
67
68       DO  i = nxl, nxr
69          DO  j = nys,nyn
70!
71!--          Compute horizontal diffusion
72             DO  k = nzb_s_outer(j,i)+1, nzt
73
74                tend(k,j,i) = tend(k,j,i)                                     &
75                                          + 0.5 * (                           &
76                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
77                      - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
78                                                  ) * ddx2                    &
79                                          + 0.5 * (                           &
80                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
81                      - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
82                                                  ) * ddy2
83             ENDDO
84
85!
86!--          Apply prescribed horizontal wall heatflux where necessary
87             IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0 ) ) &
88             THEN
89                DO  k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
90
91                   tend(k,j,i) = tend(k,j,i)                                  &
92                                          + 0.5 * ( fwxp(j,i) *               &
93                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
94                        + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)                &
95                                                   -fwxm(j,i) *               &
96                        ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
97                        + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)                &
98                                                  ) * ddx2                    &
99                                          + 0.5 * ( fwyp(j,i) *               &
100                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
101                        + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)                &
102                                                   -fwym(j,i) *               &
103                        ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
104                        + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)                &
105                                                  ) * ddy2
106                ENDDO
107             ENDIF
108
109!
110!--          Compute vertical diffusion. In case that surface fluxes have been
111!--          prescribed or computed at bottom and/or top, index k starts/ends at
112!--          nzb+2 or nzt-1, respectively.
113             DO  k = nzb_diff_s_inner(j,i), nzt_diff
114
115                tend(k,j,i) = tend(k,j,i)                                     &
116                                       + 0.5 * (                              &
117            ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
118          - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k)   &
119                                               ) * ddzw(k)
120             ENDDO
121
122!
123!--          Vertical diffusion at the first computational gridpoint along
124!--          z-direction
125             IF ( use_surface_fluxes )  THEN
126
127                k = nzb_s_inner(j,i)+1
128
129                tend(k,j,i) = tend(k,j,i)                                     &
130                                       + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )    &
131                                               * ( s(k+1,j,i)-s(k,j,i) )      &
132                                               * ddzu(k+1)                    &
133                                           + s_flux_b(j,i)                    &
134                                         ) * ddzw(k)
135
136             ENDIF
137
138!
139!--          Vertical diffusion at the last computational gridpoint along
140!--          z-direction
141             IF ( use_top_fluxes )  THEN
142
143                k = nzt
144
145                tend(k,j,i) = tend(k,j,i)                                     &
146                                       + ( - s_flux_t(j,i)                    &
147                                           - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )  &
148                                                 * ( s(k,j,i)-s(k-1,j,i) )    &
149                                                 * ddzu(k)                    &
150                                         ) * ddzw(k)
151
152             ENDIF
153
154          ENDDO
155       ENDDO
156
157    END SUBROUTINE diffusion_s
158
159
160!------------------------------------------------------------------------------!
161! Call for grid point i,j
162!------------------------------------------------------------------------------!
163    SUBROUTINE diffusion_s_ij( i, j, ddzu, ddzw, kh, s, s_flux_b, s_flux_t, &
164                               wall_s_flux, tend )
165
166       USE control_parameters
167       USE grid_variables
168       USE indices
169
170       IMPLICIT NONE
171
172       INTEGER ::  i, j, k
173       REAL    ::  vertical_gridspace
174       REAL    ::  ddzu(1:nzt+1), ddzw(1:nzt+1)
175       REAL    ::  tend(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1)
176       REAL    ::  wall_s_flux(0:4)
177       REAL, DIMENSION(:,:),   POINTER ::  s_flux_b, s_flux_t
178       REAL, DIMENSION(:,:,:), POINTER ::  kh, s
179
180!
181!--    Compute horizontal diffusion
182       DO  k = nzb_s_outer(j,i)+1, nzt
183
184          tend(k,j,i) = tend(k,j,i)                                           &
185                                          + 0.5 * (                           &
186                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
187                      - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
188                                                  ) * ddx2                    &
189                                          + 0.5 * (                           &
190                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
191                      - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
192                                                  ) * ddy2
193       ENDDO
194
195!
196!--    Apply prescribed horizontal wall heatflux where necessary
197       IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0 ) ) &
198       THEN
199          DO  k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
200
201             tend(k,j,i) = tend(k,j,i)                                        &
202                                          + 0.5 * ( fwxp(j,i) *               &
203                        ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
204                        + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)                &
205                                                   -fwxm(j,i) *               &
206                        ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
207                        + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)                &
208                                                  ) * ddx2                    &
209                                          + 0.5 * ( fwyp(j,i) *               &
210                        ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
211                        + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)                &
212                                                   -fwym(j,i) *               &
213                        ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
214                        + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)                &
215                                                  ) * ddy2
216          ENDDO
217       ENDIF
218
219!
220!--    Compute vertical diffusion. In case that surface fluxes have been
221!--    prescribed or computed at bottom and/or top, index k starts/ends at
222!--    nzb+2 or nzt-1, respectively.
223       DO  k = nzb_diff_s_inner(j,i), nzt_diff
224
225          tend(k,j,i) = tend(k,j,i)                                           &
226                                       + 0.5 * (                              &
227            ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
228          - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k)   &
229                                               ) * ddzw(k)
230       ENDDO
231
232!
233!--    Vertical diffusion at the first computational gridpoint along z-direction
234       IF ( use_surface_fluxes )  THEN
235
236          k = nzb_s_inner(j,i)+1
237
238          tend(k,j,i) = tend(k,j,i) + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )  &
239                                            * ( s(k+1,j,i)-s(k,j,i) )    &
240                                            * ddzu(k+1)                  &
241                                        + s_flux_b(j,i)                  &
242                                      ) * ddzw(k)
243
244       ENDIF
245
246!
247!--    Vertical diffusion at the last computational gridpoint along z-direction
248       IF ( use_top_fluxes )  THEN
249
250          k = nzt
251
252          tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(j,i)                  &
253                                      - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )  &
254                                            * ( s(k,j,i)-s(k-1,j,i) )    &
255                                            * ddzu(k)                    &
256                                      ) * ddzw(k)
257
258       ENDIF
259
260    END SUBROUTINE diffusion_s_ij
261
262 END MODULE diffusion_s_mod
Note: See TracBrowser for help on using the repository browser.
                                                                                                                                                                                                                                                                                                                                                                               
  | Impressum | ©Leibniz Universität Hannover |