Changeset 97 for palm/trunk/SOURCE/diffusion_e.f90

Timestamp:

Jun 21, 2007 8:23:15 AM (17 years ago)

Author:

raasch

Message:

New:
---
ocean version including prognostic equation for salinity and equation of state for seawater. Routine buoyancy can be used with both temperature and density.
+ inipar-parameters bc_sa_t, bottom_salinityflux, ocean, sa_surface, sa_vertical_gradient, sa_vertical_gradient_level, top_salinityflux

advec_s_bc, average_3d_data, boundary_conds, buoyancy, check_parameters, data_output_2d, data_output_3d, diffusion_e, flow_statistics, header, init_grid, init_3d_model, modules, netcdf, parin, production_e, prognostic_equations, read_var_list, sum_up_3d_data, swap_timelevel, time_integration, user_interface, write_var_list, write_3d_binary

New:
eqn_state_seawater, init_ocean

Changed:

---

inipar-parameter use_pt_reference renamed use_reference

hydro_press renamed hyp, routine calc_mean_pt_profile renamed calc_mean_profile

format adjustments for the ocean version (run_control)

advec_particles, buoyancy, calc_liquid_water_content, check_parameters, diffusion_e, diffusivities, header, init_cloud_physics, modules, production_e, prognostic_equations, run_control

Errors:

---

Bugfix: height above topography instead of height above level k=0 is used for calculating the mixing length (diffusion_e and diffusivities).

Bugfix: error in boundary condition for TKE removed (advec_s_bc)

advec_s_bc, diffusion_e, prognostic_equations

File:

• palm/trunk/SOURCE/diffusion_e.f90 (modified) (9 diffs)

palm/trunk/SOURCE/diffusion_e.f90

@@ -8 +8 @@
 ! This is also a bugfix, because the height above the topography is now
 ! used instead of the height above level k=0.
+! theta renamed var, dpt_dz renamed dvar_dz, +new argument var_reference
+! use_pt_reference renamed use_reference
 !
 ! Former revisions:
@@ -50 +52 @@
 ! Call for all grid points
 !------------------------------------------------------------------------------!
-    SUBROUTINE diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, theta, &
-                            rif, tend, zu, zw )
+    SUBROUTINE diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, var, &
+                            var_reference, rif, tend, zu, zw )
 
        USE control_parameters
@@ -61 +63 @@
 
        INTEGER ::  i, j, k
-       REAL            ::  dpt_dz, l_stable, phi_m
+       REAL            ::  dvar_dz, l_stable, phi_m, var_reference
        REAL            ::  ddzu(1:nzt+1), dd2zu(1:nzt), ddzw(1:nzt+1), &
                            l_grid(1:nzt), zu(0:nzt+1), zw(0:nzt+1)
        REAL, DIMENSION(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) :: diss, tend
        REAL, DIMENSION(:,:), POINTER   ::  rif
-       REAL, DIMENSION(:,:,:), POINTER ::  e, km, theta
+       REAL, DIMENSION(:,:,:), POINTER ::  e, km, var
        REAL, DIMENSION(nzb+1:nzt,nys:nyn) ::  dissipation, l, ll
  
@@ -73 +75 @@
 !--    This if clause must be outside the k-loop because otherwise
 !--    runtime errors occur with -C hopt on NEC
-       IF ( use_pt_reference )  THEN
+       IF ( use_reference )  THEN
 
           DO  i = nxl, nxr
@@ -91 +93 @@
 !
 !--                Calculate the mixing length (for dissipation)
-                   dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k)
-                   IF ( dpt_dz > 0.0 ) THEN
+                   dvar_dz = atmos_ocean_sign * &
+                             ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
+                   IF ( dvar_dz > 0.0 ) THEN
                       l_stable = 0.76 * SQRT( e(k,j,i) ) / &
-                                        SQRT( g / pt_reference * dpt_dz ) + 1E-5
+                                 SQRT( g / var_reference * dvar_dz ) + 1E-5
                    ELSE
                       l_stable = l_grid(k)
@@ -180 +183 @@
 !
 !--                Calculate the mixing length (for dissipation)
-                   dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k)
-                   IF ( dpt_dz > 0.0 ) THEN
+                   dvar_dz = atmos_ocean_sign * &
+                             ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
+                   IF ( dvar_dz > 0.0 ) THEN
                       l_stable = 0.76 * SQRT( e(k,j,i) ) / &
-                                        SQRT( g / theta(k,j,i) * dpt_dz ) + 1E-5
+                                        SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5
                    ELSE
                       l_stable = l_grid(k)
@@ -270 +274 @@
 !------------------------------------------------------------------------------!
     SUBROUTINE diffusion_e_ij( i, j, ddzu, dd2zu, ddzw, diss, e, km, l_grid, &
-                               theta, rif, tend, zu, zw )
+                               var, var_reference, rif, tend, zu, zw )
 
        USE control_parameters
@@ -280 +284 @@
 
        INTEGER         ::  i, j, k
-       REAL            ::  dpt_dz, l_stable, phi_m
+       REAL            ::  dvar_dz, l_stable, phi_m, var_reference
        REAL            ::  ddzu(1:nzt+1), dd2zu(1:nzt), ddzw(1:nzt+1), &
                            l_grid(1:nzt), zu(0:nzt+1), zw(0:nzt+1)
        REAL, DIMENSION(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) ::  diss, tend
        REAL, DIMENSION(:,:), POINTER   ::  rif
-       REAL, DIMENSION(:,:,:), POINTER ::  e, km, theta
+       REAL, DIMENSION(:,:,:), POINTER ::  e, km, var
        REAL, DIMENSION(nzb+1:nzt)    ::  dissipation, l, ll
 
@@ -303 +307 @@
 !--    Calculate the mixing length (for dissipation)
        DO  k = nzb_s_inner(j,i)+1, nzt
-          dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k)
-          IF ( dpt_dz > 0.0 ) THEN
-             IF ( use_pt_reference )  THEN
+          dvar_dz = atmos_ocean_sign * &
+                    ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
+          IF ( dvar_dz > 0.0 ) THEN
+             IF ( use_reference )  THEN
                 l_stable = 0.76 * SQRT( e(k,j,i) ) / &
-                                  SQRT( g / pt_reference * dpt_dz ) + 1E-5
+                                  SQRT( g / var_reference * dvar_dz ) + 1E-5
              ELSE
                 l_stable = 0.76 * SQRT( e(k,j,i) ) / &
-                                  SQRT( g / theta(k,j,i) * dpt_dz ) + 1E-5
+                                  SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5
              ENDIF
           ELSE