Changeset 95 for palm/trunk/SOURCE


Ignore:
Timestamp:
Jun 2, 2007 4:48:38 PM (18 years ago)
Author:
raasch
Message:

further preliminary uncomplete changes for ocean version

Location:
palm/trunk/SOURCE
Files:
1 added
18 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/CURRENT_MODIFICATIONS

    r94 r95  
    22---
    33ocean version including prognostic equation for salinity and equation of state for seawater
    4 + inipar-parameters ocean, sa_surface, sa_vertical_gradient, sa_vertical_gradient_level, top_salinity_flux
     4+ inipar-parameters bc_sa_t, bottom_salinityflux, ocean, sa_surface, sa_vertical_gradient, sa_vertical_gradient_level, top_salinityflux
    55
    6 check_parameters, diffusion_e, flow_statistics, header, init_grid, init_3d_model, modules, parin, prognostic_equations, read_var_list, write_var_list, write_3d_binary
     6boundary_conds, check_parameters, diffusion_e, flow_statistics, header, init_grid, init_3d_model, modules, parin, prognostic_equations, read_var_list, swap_timelevel, time_integration, user_interface, write_var_list, write_3d_binary
    77
     8New:
     9eqn_state_seawater, init_ocean
    810
    911Changed:
    1012-------
     13hydro_press renamed hyp
    1114
    12 
    13 
     15advec_particles
    1416
    1517
     
    3638
    3739default disturbance level
     40
     41Bott-Chlond scheme to be extended for salinity
     42
     43Upstream-spline scheme to be extended for salinity?????
  • palm/trunk/SOURCE/Makefile

    r82 r95  
    44# Actual revisions:
    55# -----------------
    6 # +local_flush
     6# +eqn_state_seawater, init_ocean
    77#
    88# Former revisions:
    99# -----------------
    1010# $Id$
     11#
     12# 82 2007-04-16 15:40:52Z raasch
     13# +local_flush
    1114#
    1215# 58 2007-03-09 14:27:38Z raasch
     
    4548        data_output_3d.f90 diffusion_e.f90 diffusion_s.f90 diffusion_u.f90 \
    4649        diffusion_v.f90 diffusion_w.f90 diffusivities.f90 disturb_field.f90 \
    47         disturb_heatflux.f90 \exchange_horiz.f90 exchange_horiz_2d.f90 \
     50        disturb_heatflux.f90 eqn_state_seawater.f90 exchange_horiz.f90 exchange_horiz_2d.f90 \
    4851        fft_xy.f90 flow_statistics.f90 global_min_max.f90 header.f90 \
    4952        impact_of_latent_heat.f90 init_1d_model.f90 init_3d_model.f90 \
    5053        init_advec.f90 init_cloud_physics.f90 init_dvrp.f90 init_grid.f90 \
    51         init_particles.f90 init_pegrid.f90 init_pt_anomaly.f90 \
     54        init_ocean.f90 init_particles.f90 init_pegrid.f90 init_pt_anomaly.f90 \
    5255        init_rankine.f90 init_slope.f90 interaction_droplets_ptq.f90 \
    5356        local_flush.f90 local_getenv.f90 local_stop.f90 local_system.f90 local_tremain.f90 \
     
    7679        data_output_3d.o diffusion_e.o diffusion_s.o diffusion_u.o \
    7780        diffusion_v.o diffusion_w.o diffusivities.o disturb_field.o \
    78         disturb_heatflux.o exchange_horiz.o exchange_horiz_2d.o fft_xy.o \
     81        disturb_heatflux.o eqn_state_seawater.o exchange_horiz.o exchange_horiz_2d.o fft_xy.o \
    7982        flow_statistics.o global_min_max.o header.o impact_of_latent_heat.o \
    8083        init_1d_model.o init_3d_model.o init_advec.o init_cloud_physics.o \
    81         init_dvrp.o init_grid.o init_particles.o init_pegrid.o \
     84        init_dvrp.o init_grid.o init_ocean.o init_particles.o init_pegrid.o \
    8285        init_pt_anomaly.o init_rankine.o init_slope.o \
    8386        interaction_droplets_ptq.o local_flush.o local_getenv.o local_stop.o \
     
    164167disturb_field.o: modules.o random_function.o
    165168disturb_heatflux.o: modules.o
     169eqn_state_seawater.o: modules.o
    166170exchange_horiz.o: modules.o
    167171exchange_horiz_2d.o: modules.o
     
    177181init_dvrp.o: modules.o
    178182init_grid.o: modules.o
     183init_ocean.o: modules.o eqn_state_seawater.o
    179184init_particles.o: modules.o random_function.o
    180185init_pegrid.o: modules.o fft_xy.o poisfft.o poisfft_hybrid.o
  • palm/trunk/SOURCE/advec_particles.f90

    r82 r95  
    44! Actual revisions:
    55! -----------------
     6! hydro_press renamed hyp
    67! TEST: PRINT statements on unit 9 (commented out)
    78!
     
    238239!
    239240!--       Calculate real temperature and saturation vapor pressure
    240           p_int = hydro_press(k) + ( particles(n)%z - zu(k) ) / dz * &
    241                                    ( hydro_press(k+1) - hydro_press(k) )
     241          p_int = hyp(k) + ( particles(n)%z - zu(k) ) / dz * ( hyp(k+1)-hyp(k) )
    242242          t_int = pt_int * ( p_int / 100000.0 )**0.286
    243243
  • palm/trunk/SOURCE/boundary_conds.f90

    r77 r95  
    44! Actual revisions:
    55! -----------------
    6 !
     6! Boundary conditions for salinity added
    77!
    88! Former revisions:
     
    128128
    129129!
     130!--    Boundary conditions for salinity
     131       IF ( ocean )  THEN
     132!
     133!--       Bottom boundary: Neumann condition because salinity flux is always
     134!--       given
     135          DO  i = nxl-1, nxr+1
     136             DO  j = nys-1, nyn+1
     137                sa_p(nzb_s_inner(j,i),j,i) = sa_p(nzb_s_inner(j,i)+1,j,i)
     138             ENDDO
     139          ENDDO
     140
     141!
     142!--       Top boundary: Dirichlet or Neumann
     143          IF ( ibc_sa_t == 0 )  THEN
     144             sa_p(nzt+1,:,:) = sa(nzt+1,:,:)
     145          ELSEIF ( ibc_sa_t == 1 )  THEN
     146             sa_p(nzt+1,:,:) = sa_p(nzt,:,:)
     147          ENDIF
     148
     149       ENDIF
     150
     151!
    130152!--    Boundary conditions for total water content or scalar,
    131 !--    bottom and surface boundary (see also temperature)
     153!--    bottom and top boundary (see also temperature)
    132154       IF ( humidity  .OR.  passive_scalar )  THEN
    133155!
  • palm/trunk/SOURCE/calc_liquid_water_content.f90

    r39 r95  
    44! Actual revisions:
    55! -----------------
    6 !
     6! hydro_press renamed hyp
    77!
    88! Former revisions:
     
    5959!
    6060!--          Compute approximation of saturation humidity
    61              q_s = 0.622 * e_s / &
    62                    ( hydro_press(k) - 0.378 * e_s )
     61             q_s = 0.622 * e_s / ( hyp(k) - 0.378 * e_s )
    6362
    6463!
  • palm/trunk/SOURCE/check_parameters.f90

    r94 r95  
    10051005
    10061006!
     1007!-- Boundary conditions for salinity
     1008    IF ( ocean )  THEN
     1009       IF ( bc_sa_t == 'dirichlet' )  THEN
     1010          ibc_sa_t = 0
     1011       ELSEIF ( bc_sa_t == 'neumann' )  THEN
     1012          ibc_sa_t = 1
     1013       ELSE
     1014          IF ( myid == 0 )  THEN
     1015             PRINT*, '+++ check_parameters:'
     1016             PRINT*, '    unknown boundary condition: bc_sa_t = ', bc_sa_t
     1017          ENDIF
     1018          CALL local_stop
     1019       ENDIF
     1020
     1021       IF ( top_salinityflux == 9999999.9 )  constant_top_salinityflux = .FALSE.
     1022
     1023!
     1024!--    A fixed salinity at the top implies Dirichlet boundary condition for
     1025!--    salinity. In this case specification of a constant salinity flux is
     1026!--    forbidden.
     1027       IF ( ibc_sa_t == 0  .AND.   constant_top_salinityflux  .AND. &
     1028            top_salinityflux /= 0.0 )  THEN
     1029          IF ( myid == 0 )  THEN
     1030             PRINT*, '+++ check_parameters:'
     1031             PRINT*, '    boundary_condition: bc_sa_t = ', bc_sa_t
     1032             PRINT*, '    is not allowed with constant_top_salinityflux = ', &
     1033                          '.TRUE.'
     1034          ENDIF
     1035          CALL local_stop
     1036       ENDIF
     1037
     1038    ENDIF
     1039
     1040!
    10071041!-- In case of humidity or passive scalar, set boundary conditions for total
    10081042!-- water content / scalar
  • palm/trunk/SOURCE/init_3d_model.f90

    r94 r95  
    77! Actual revisions:
    88! -----------------
    9 ! Initialization of salinity
     9! Initialization of salinity, call of init_ocean
    1010!
    1111! Former revisions:
     
    188188
    189189    IF ( ocean )  THEN
    190        ALLOCATE( saswst_1(nys-1:nyn+1,nxl-1:nxr+1) )
    191        ALLOCATE( r(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1),    &
     190       ALLOCATE( saswsb_1(nys-1:nyn+1,nxl-1:nxr+1), &
     191                 saswst_1(nys-1:nyn+1,nxl-1:nxr+1) )
     192       ALLOCATE( rho(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1),  &
    192193                 sa_1(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1), &
    193194                 sa_2(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1), &
     
    303304
    304305       IF ( ocean )  THEN
     306          saswsb => saswsb_1
    305307          saswst => saswst_1
    306308          sa     => sa_1;    sa_p => sa_2;    tsa_m => sa_3
     
    602604
    603605             IF ( ocean )  THEN
     606                saswsb = bottom_salinityflux
    604607                saswst = top_salinityflux
    605608             ENDIF
     
    821824       CALL local_stop
    822825    ENDIF
     826
     827!
     828!-- If required, initialize quantities needed for the ocean model
     829    IF ( ocean )  CALL init_ocean
    823830
    824831!
  • palm/trunk/SOURCE/init_cloud_physics.f90

    r4 r95  
    44! Actual revisions:
    55! -----------------
    6 !
     6! hydro_press renamed hyp
    77!
    88! Former revisions:
     
    3535    REAL    ::  t_surface
    3636
    37     ALLOCATE( hydro_press(nzb:nzt+1), pt_d_t(nzb:nzt+1), t_d_pt(nzb:nzt+1) )
     37    ALLOCATE( hyp(nzb:nzt+1), pt_d_t(nzb:nzt+1), t_d_pt(nzb:nzt+1) )
    3838
    3939!
     
    5454!-- pt / t : ratio of potential and actual temperature (pt_d_t)
    5555!-- t / pt : ratio of actual and potential temperature (t_d_pt)
    56 !-- p_0(z) : vertical profile of the hydrostatic pressure (hydro_press)
     56!-- p_0(z) : vertical profile of the hydrostatic pressure (hyp)
    5757    t_surface = pt_surface * ( surface_pressure / 1000.0 )**0.286
    5858    DO  k = nzb, nzt+1
    59        hydro_press(k) = surface_pressure * 100.0 * &
    60                         ( (t_surface - g/cp * zu(k)) / t_surface )**(1.0/0.286)
    61        pt_d_t(k)      = ( 100000.0 / hydro_press(k) )**0.286
    62        t_d_pt(k)      = 1.0 / pt_d_t(k)       
     59       hyp(k)    = surface_pressure * 100.0 * &
     60                   ( (t_surface - g/cp * zu(k)) / t_surface )**(1.0/0.286)
     61       pt_d_t(k) = ( 100000.0 / hyp(k) )**0.286
     62       t_d_pt(k) = 1.0 / pt_d_t(k)       
    6363    ENDDO
    6464
  • palm/trunk/SOURCE/modules.f90

    r94 r95  
    77! +ocean, r, + salinity variables
    88! defaults of .._vertical_gradient_levels changed from -1.0 to -9999999.9
     9! hydro_press renamed hyp
    910!
    1011! Former revisions:
     
    9495
    9596    REAL, DIMENSION(:), ALLOCATABLE ::                                         &
    96           ddzu, dd2zu, dzu, ddzw, dzw, km_damp_x, km_damp_y, l_grid, pt_init,  &
    97           q_init, rdf, sa_init, ug, u_init, u_nzb_p1_for_vfc, vg, v_init,      &
    98           v_nzb_p1_for_vfc, zu, zw
     97          ddzu, dd2zu, dzu, ddzw, dzw, hyp, km_damp_x, km_damp_y, l_grid,      &
     98          pt_init, q_init, rdf, sa_init, ug, u_init, u_nzb_p1_for_vfc, vg,     &
     99          v_init, v_nzb_p1_for_vfc, zu, zw
    99100
    100101    REAL, DIMENSION(:,:), ALLOCATABLE ::                                       &
     
    102103
    103104    REAL, DIMENSION(:,:), ALLOCATABLE, TARGET ::                               &
    104           qsws_1, qsws_2, qswst_1, qswst_2, rif_1, rif_2, saswst_1, shf_1,     &
    105           shf_2, tswst_1, tswst_2, usws_1, usws_2, vsws_1, vsws_2
     105          qsws_1, qsws_2, qswst_1, qswst_2, rif_1, rif_2, saswsb_1, saswst_1,  &
     106          shf_1, shf_2, tswst_1, tswst_2, usws_1, usws_2, vsws_1, vsws_2
    106107
    107108    REAL, DIMENSION(:,:), POINTER ::                                           &
    108           qsws, qsws_m, qswst, qswst_m, rif, rif_m, saswst, shf, shf_m, tswst, &
    109           tswst_m, usws, usws_m, vsws, vsws_m
     109          qsws, qsws_m, qswst, qswst_m, rif, rif_m, saswsb, saswst, shf,      &
     110          shf_m, tswst, tswst_m, usws, usws_m, vsws, vsws_m
    110111
    111112    REAL, DIMENSION(:,:,:), ALLOCATABLE ::                                     &
     
    118119    REAL, DIMENSION(:,:,:), ALLOCATABLE, TARGET ::                             &
    119120          e_1, e_2, e_3, kh_1, kh_2, km_1, km_2, p, pt_1, pt_2, pt_3, q_1,     &
    120           q_2, q_3, ql_1, ql_2, r, sa_1, sa_2, sa_3, u_1, u_2, u_3, v_1, v_2,  &
    121           v_3, vpt_1, vpt_2, w_1, w_2, w_3
     121          q_2, q_3, ql_1, ql_2, rho, sa_1, sa_2, sa_3, u_1, u_2, u_3, v_1,     &
     122          v_2, v_3, vpt_1, vpt_2, w_1, w_2, w_3
    122123
    123124    REAL, DIMENSION(:,:,:), POINTER ::                                         &
     
    167168              mass_of_solute, molecular_weight_of_solute,                      &
    168169              prec_time_const = 0.001, ql_crit = 0.0005, rho_l = 1.0E3,        &
    169               r_d = 287.0, r_v = 461.51, rho_surface,                          &
    170               thermal_conductivity_l = 2.43E-2
    171 
    172     REAL, DIMENSION(:), ALLOCATABLE :: hydro_press, pt_d_t, t_d_pt
     170              r_d = 287.0, r_v = 461.51, thermal_conductivity_l = 2.43E-2
     171
     172    REAL, DIMENSION(:), ALLOCATABLE   ::  pt_d_t, t_d_pt
    173173
    174174    REAL, DIMENSION(:,:), ALLOCATABLE ::  precipitation_amount, &
     
    241241                            bc_q_b = 'dirichlet', bc_q_t = 'neumann', &
    242242                            bc_s_b = 'dirichlet', bc_s_t = 'neumann', &
     243                            bc_sa_t = 'neumann', &
    243244                            bc_uv_b = 'dirichlet', bc_uv_t = 'dirichlet', &
    244245                            dissipation_1d = 'as_in_3d_model', &
     
    272273                dvrp_filecount = 0, dz_stretch_level_index, gamma_mg, &
    273274                grid_level, ibc_e_b, ibc_p_b, ibc_p_t, ibc_pt_b, ibc_pt_t, &
    274                 ibc_q_b, ibc_q_t, ibc_uv_b, ibc_uv_t, &
     275                ibc_q_b, ibc_q_t, ibc_sa_t, ibc_uv_b, ibc_uv_t, &
    275276                inflow_disturbance_begin = -1, inflow_disturbance_end = -1, &
    276277                intermediate_timestep_count, intermediate_timestep_count_max, &
     
    302303                conserve_volume_flow = .FALSE., constant_diffusion = .FALSE., &
    303304                constant_heatflux = .TRUE., constant_top_heatflux = .TRUE., &
     305                constant_top_salinityflux = .TRUE., &
    304306                constant_waterflux = .TRUE., create_disturbances = .TRUE., &
    305307                cut_spline_overshoot = .TRUE., &
     
    333335             averaging_interval = 0.0, averaging_interval_pr = 9999999.9, &
    334336             averaging_interval_sp = 9999999.9, bc_pt_t_val, bc_q_t_val, &
     337             bottom_salinityflux = 0.0, &
    335338             building_height = 50.0, building_length_x = 50.0, &
    336339             building_length_y = 50.0, building_wall_left = 9999999.9, &
     
    363366             q_surface_initial_change = 0.0, rayleigh_damping_factor = -1.0, &
    364367             rayleigh_damping_height = -1.0, residual_limit = 1.0E-4, &
    365              restart_time = 9999999.9, rif_max = 1.0, rif_min = -5.0, &
    366              roughness_length = 0.1, sa_surface = 35.0, simulated_time = 0.0, &
    367              simulated_time_at_begin, sin_alpha_surface, &
     368             restart_time = 9999999.9, rho_ref, rho_surface, rif_max = 1.0, &
     369             rif_min = -5.0, roughness_length = 0.1, sa_surface = 35.0, &
     370             simulated_time = 0.0, simulated_time_at_begin, sin_alpha_surface, &
    368371             skip_time_data_output = 0.0, skip_time_data_output_av = 9999999.9,&
    369372             skip_time_dopr = 9999999.9, skip_time_dosp = 9999999.9, &
     
    380383             time_do_sla = 0.0, time_dvrp = 0.0, time_prel = 0.0, &
    381384             time_restart = 9999999.9, time_run_control = 0.0, &
    382              top_heatflux = 9999999.9, top_salinityflux = 0.0, &
     385             top_heatflux = 9999999.9, top_salinityflux = 9999999.9, &
    383386             ug_surface = 0.0, u_gtrans = 0.0, &
    384387             ups_limit_e = 0.0, ups_limit_pt = 0.0, ups_limit_u = 0.0, &
     
    393396             q_vertical_gradient_level(10) = -1.0, &
    394397             s_vertical_gradient(10) = 0.0, &
    395              s_vertical_gradient_level(10) = -1.0,
     398             s_vertical_gradient_level(10) = -1.0, &
    396399             sa_vertical_gradient(10) = 0.0, &
    397              sa_vertical_gradient_level(10) = -1.0, threshold(20) = 0.0, &
     400             sa_vertical_gradient_level(10) = -9999999.9, threshold(20) = 0.0, &
    398401             tsc(10) = (/ 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 /), &
    399402             ug_vertical_gradient(10) = 0.0, &
  • palm/trunk/SOURCE/parin.f90

    r94 r95  
    44! Actual revisions:
    55! -----------------
    6 ! +ocean, sa_surface, sa_vertical_gradient, sa_vertical_gradient_level,
    7 ! top_salinity_flux in inipar,
     6! +bc_sa_t, bottom_salinityflux, ocean, sa_surface, sa_vertical_gradient,
     7! sa_vertical_gradient_level, top_salinityflux in inipar,
    88! sa_init is allocated
    99!
     
    6363    NAMELIST /inipar/  adjust_mixing_length, alpha_surface, bc_e_b, bc_lr, &
    6464                       bc_ns, bc_p_b, bc_p_t, bc_pt_b, bc_pt_t, bc_q_b, &
    65              bc_q_t,bc_s_b, bc_s_t, bc_uv_b, bc_uv_t, building_height, &
     65             bc_q_t,bc_s_b, bc_s_t, bc_sa_t, bc_uv_b, bc_uv_t, &
     66             bottom_salinityflux, building_height, &
    6667             building_length_x, building_length_y, building_wall_left, &
    6768             building_wall_south, cloud_droplets, cloud_physics, &
     
    8788             surface_waterflux, s_surface, s_surface_initial_change, &
    8889             s_vertical_gradient, s_vertical_gradient_level, top_heatflux, &
    89              top_salinity_flux, timestep_scheme, topography, ug_surface, &
     90             top_salinityflux, timestep_scheme, topography, ug_surface, &
    9091             ug_vertical_gradient, ug_vertical_gradient_level, ups_limit_e, &
    9192             ups_limit_pt, ups_limit_u, ups_limit_v, ups_limit_w, &
  • palm/trunk/SOURCE/prognostic_equations.f90

    r94 r95  
    44! Actual revisions:
    55! -----------------
     6! prognostic equation for salinity
    67! new argument zw in calls of diffusion_e
    78!
     
    414415
    415416!
     417!-- If required, compute prognostic equation for salinity
     418    IF ( ocean )  THEN
     419
     420       CALL cpu_log( log_point(37), 'sa-equation', 'start' )
     421
     422!
     423!--    sa-tendency terms with communication
     424       sat = tsc(1)
     425       sbt = tsc(2)
     426       IF ( scalar_advec == 'bc-scheme' )  THEN
     427
     428          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
     429!
     430!--          Bott-Chlond scheme always uses Euler time step when leapfrog is
     431!--          switched on. Thus:
     432             sat = 1.0
     433             sbt = 1.0
     434          ENDIF
     435          tend = 0.0
     436          CALL advec_s_bc( sa, 'sa' )
     437       ELSE
     438          IF ( tsc(2) /= 2.0 )  THEN
     439             IF ( scalar_advec == 'ups-scheme' )  THEN
     440                tend = 0.0
     441                CALL advec_s_ups( sa, 'sa' )
     442             ENDIF
     443          ENDIF
     444       ENDIF
     445
     446!
     447!--    sa terms with no communication
     448       DO  i = nxl, nxr
     449          DO  j = nys, nyn
     450!
     451!--          Tendency-terms
     452             IF ( scalar_advec == 'bc-scheme' )  THEN
     453                CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
     454                                  tend )
     455             ELSE
     456                IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) THEN
     457                   tend(:,j,i) = 0.0
     458                   CALL advec_s_pw( i, j, sa )
     459                ELSE
     460                   IF ( scalar_advec /= 'ups-scheme' )  THEN
     461                      tend(:,j,i) = 0.0
     462                      CALL advec_s_up( i, j, sa )
     463                   ENDIF
     464                ENDIF
     465                CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
     466                                  tend )
     467             ENDIF
     468       
     469             CALL user_actions( i, j, 'sa-tendency' )
     470
     471!
     472!--          Prognostic equation for salinity
     473             DO  k = nzb_s_inner(j,i)+1, nzt
     474                sa_p(k,j,i) = sat * sa(k,j,i) +                                &
     475                              dt_3d * (                                        &
     476                                     sbt * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
     477                                      ) -                                      &
     478                              tsc(5) * rdf(k) * ( sa(k,j,i) - sa_init(k) )
     479                IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
     480             ENDDO
     481
     482!
     483!--          Calculate tendencies for the next Runge-Kutta step
     484             IF ( timestep_scheme(1:5) == 'runge' )  THEN
     485                IF ( intermediate_timestep_count == 1 )  THEN
     486                   DO  k = nzb_s_inner(j,i)+1, nzt
     487                      tsa_m(k,j,i) = tend(k,j,i)
     488                   ENDDO
     489                ELSEIF ( intermediate_timestep_count < &
     490                         intermediate_timestep_count_max )  THEN
     491                   DO  k = nzb_s_inner(j,i)+1, nzt
     492                      tsa_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     493                                      5.3125 * tsa_m(k,j,i)
     494                   ENDDO
     495                ENDIF
     496             ENDIF
     497
     498          ENDDO
     499       ENDDO
     500
     501       CALL cpu_log( log_point(37), 'sa-equation', 'stop' )
     502
     503    ENDIF
     504
     505!
    416506!-- If required, compute prognostic equation for total water content / scalar
    417507    IF ( humidity  .OR.  passive_scalar )  THEN
     
    647737! be called for the standard Piascek-Williams advection scheme.
    648738!
    649 ! The call of this subroutine is embedded in two DO loops over i and j, thus
    650 ! communication between CPUs is not allowed in this subroutine.
     739! Here the calls of most subroutines are embedded in two DO loops over i and j,
     740! so communication between CPUs is not allowed (does not make sense) within
     741! these loops.
    651742!
    652743! (Optimized to avoid cache missings, i.e. for Power4/5-architectures.)
     
    883974                   ENDDO
    884975                ENDIF
     976             ENDIF
     977
     978!
     979!--          If required, compute prognostic equation for salinity
     980             IF ( ocean )  THEN
     981
     982!
     983!--             Tendency-terms for salinity
     984                tend(:,j,i) = 0.0
     985                IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
     986                THEN
     987                   CALL advec_s_pw( i, j, sa )
     988                ELSE
     989                   CALL advec_s_up( i, j, sa )
     990                ENDIF
     991                CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
     992                                  tend )
     993       
     994                CALL user_actions( i, j, 'sa-tendency' )
     995
     996!
     997!--             Prognostic equation for salinity
     998                DO  k = nzb_s_inner(j,i)+1, nzt
     999                   sa_p(k,j,i) = tsc(1) * sa(k,j,i) +                          &
     1000                                 dt_3d * (                                     &
     1001                                  tsc(2) * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
     1002                                         ) -                                   &
     1003                                tsc(5) * rdf(k) * ( sa(k,j,i) - sa_init(k) )
     1004                   IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
     1005                ENDDO
     1006
     1007!
     1008!--             Calculate tendencies for the next Runge-Kutta step
     1009                IF ( timestep_scheme(1:5) == 'runge' )  THEN
     1010                   IF ( intermediate_timestep_count == 1 )  THEN
     1011                      DO  k = nzb_s_inner(j,i)+1, nzt
     1012                         tsa_m(k,j,i) = tend(k,j,i)
     1013                      ENDDO
     1014                   ELSEIF ( intermediate_timestep_count < &
     1015                            intermediate_timestep_count_max )  THEN
     1016                      DO  k = nzb_s_inner(j,i)+1, nzt
     1017                         tsa_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     1018                                         5.3125 * tsa_m(k,j,i)
     1019                      ENDDO
     1020                   ENDIF
     1021                ENDIF
     1022
    8851023             ENDIF
    8861024
     
    13631501
    13641502!
     1503!-- If required, compute prognostic equation for salinity
     1504    IF ( ocean )  THEN
     1505
     1506       CALL cpu_log( log_point(37), 'sa-equation', 'start' )
     1507
     1508!
     1509!--    sa-tendency terms with communication
     1510       sat = tsc(1)
     1511       sbt = tsc(2)
     1512       IF ( scalar_advec == 'bc-scheme' )  THEN
     1513
     1514          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
     1515!
     1516!--          Bott-Chlond scheme always uses Euler time step when leapfrog is
     1517!--          switched on. Thus:
     1518             sat = 1.0
     1519             sbt = 1.0
     1520          ENDIF
     1521          tend = 0.0
     1522          CALL advec_s_bc( sa, 'sa' )
     1523       ELSE
     1524          IF ( tsc(2) /= 2.0 )  THEN
     1525             IF ( scalar_advec == 'ups-scheme' )  THEN
     1526                tend = 0.0
     1527                CALL advec_s_ups( sa, 'sa' )
     1528             ENDIF
     1529          ENDIF
     1530       ENDIF
     1531
     1532!
     1533!--    Scalar/q-tendency terms with no communication
     1534       IF ( scalar_advec == 'bc-scheme' )  THEN
     1535          CALL diffusion_s( ddzu, ddzw, kh, sa, saswsb, saswst, tend )
     1536       ELSE
     1537          IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     1538             tend = 0.0
     1539             CALL advec_s_pw( sa )
     1540          ELSE
     1541             IF ( scalar_advec /= 'ups-scheme' )  THEN
     1542                tend = 0.0
     1543                CALL advec_s_up( sa )
     1544             ENDIF
     1545          ENDIF
     1546          CALL diffusion_s( ddzu, ddzw, kh, sa, saswsb, saswst, tend )
     1547       ENDIF
     1548       
     1549       CALL user_actions( 'sa-tendency' )
     1550
     1551!
     1552!--    Prognostic equation for salinity
     1553       DO  i = nxl, nxr
     1554          DO  j = nys, nyn
     1555             DO  k = nzb_s_inner(j,i)+1, nzt
     1556                sa_p(k,j,i) = sat * sa(k,j,i) +                                &
     1557                              dt_3d * (                                        &
     1558                                     sbt * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
     1559                                      ) -                                      &
     1560                              tsc(5) * rdf(k) * ( sa(k,j,i) - sa_init(k) )
     1561                IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
     1562             ENDDO
     1563          ENDDO
     1564       ENDDO
     1565
     1566!
     1567!--    Calculate tendencies for the next Runge-Kutta step
     1568       IF ( timestep_scheme(1:5) == 'runge' )  THEN
     1569          IF ( intermediate_timestep_count == 1 )  THEN
     1570             DO  i = nxl, nxr
     1571                DO  j = nys, nyn
     1572                   DO  k = nzb_s_inner(j,i)+1, nzt
     1573                      tsa_m(k,j,i) = tend(k,j,i)
     1574                   ENDDO
     1575                ENDDO
     1576             ENDDO
     1577          ELSEIF ( intermediate_timestep_count < &
     1578                   intermediate_timestep_count_max )  THEN
     1579             DO  i = nxl, nxr
     1580                DO  j = nys, nyn
     1581                   DO  k = nzb_s_inner(j,i)+1, nzt
     1582                      tsa_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     1583                                      5.3125 * tsa_m(k,j,i)
     1584                   ENDDO
     1585                ENDDO
     1586             ENDDO
     1587          ENDIF
     1588       ENDIF
     1589
     1590       CALL cpu_log( log_point(37), 'sa-equation', 'stop' )
     1591
     1592    ENDIF
     1593
     1594!
    13651595!-- If required, compute prognostic equation for total water content / scalar
    13661596    IF ( humidity  .OR.  passive_scalar )  THEN
  • palm/trunk/SOURCE/read_3d_binary.f90

    r94 r95  
    44! Actual revisions:
    55! -----------------
    6 ! +sa, saswst
     6! +sa, saswsb, saswst
    77!
    88! Former revisions:
     
    284284          CASE ( 'sa' )
    285285             READ ( 13 )  sa
     286          CASE ( 'saswsb' )
     287             READ ( 13 )  saswsb
    286288          CASE ( 'saswst' )
    287289             READ ( 13 )  saswst
  • palm/trunk/SOURCE/read_var_list.f90

    r94 r95  
    44! Actual revisions:
    55! -----------------
    6 ! +ocean, sa_init, sa_surface, sa_vertical_gradient,
    7 ! sa_vertical_gradient_level, top_salinity_flux
     6! +bc_sa_t, ocean, sa_init, sa_surface, sa_vertical_gradient,
     7! sa_vertical_gradient_level, bottom/top_salinity_flux
    88!
    99! Former revisions:
     
    161161          CASE ( 'bc_s_t' )
    162162             READ ( 13 )  bc_s_t
     163          CASE ( 'bc_sa_t' )
     164             READ ( 13 )  bc_sa_t
    163165          CASE ( 'bc_uv_b' )
    164166             READ ( 13 )  bc_uv_b
    165167          CASE ( 'bc_uv_t' )
    166168             READ ( 13 )  bc_uv_t
     169          CASE ( 'bottom_salinityflux' )
     170             READ ( 13 )  bottom_salinityflux
    167171          CASE ( 'building_height' )
    168172             READ ( 13 )  building_height
  • palm/trunk/SOURCE/swap_timelevel.f90

    r77 r95  
    44! Actual revisions:
    55! -----------------
    6 !
     6! Swaping of salinity
    77!
    88! Former revisions:
     
    106106             pt => pt_1;  pt_p => pt_2
    107107             IF ( .NOT. constant_diffusion )  THEN
    108                 e => e_1;  e_p => e_2
    109              ENDIF
    110              IF ( humidity  .OR.  passive_scalar )  THEN
    111                 q => q_1;  q_p => q_2
     108                e => e_1;    e_p => e_2
     109             ENDIF
     110             IF ( ocean )  THEN
     111                sa => sa_1;  sa_p => sa_2
     112             ENDIF
     113             IF ( humidity  .OR.  passive_scalar )  THEN
     114                q => q_1;    q_p => q_2
    112115             ENDIF
    113116
     
    152155             pt => pt_2;  pt_p => pt_1
    153156             IF ( .NOT. constant_diffusion )  THEN
    154                 e => e_2;  e_p => e_1
    155              ENDIF
    156              IF ( humidity  .OR.  passive_scalar )  THEN
    157                 q => q_2;  q_p => q_1
     157                e => e_2;    e_p => e_1
     158             ENDIF
     159             IF ( ocean )  THEN
     160                sa => sa_2;  sa_p => sa_1
     161             ENDIF
     162             IF ( humidity  .OR.  passive_scalar )  THEN
     163                q => q_2;    q_p => q_1
    158164             ENDIF
    159165
  • palm/trunk/SOURCE/time_integration.f90

    r90 r95  
    44! Actual revisions:
    55! -----------------
    6 !
     6! Ghostpoint exchange for salinity and density
    77!
    88! Former revisions:
     
    150150          CALL exchange_horiz( pt_p )
    151151          IF ( .NOT. constant_diffusion       )  CALL exchange_horiz( e_p )
     152          IF ( ocean )  THEN
     153             CALL exchange_horiz( sa_p )
     154             CALL exchange_horiz( rho )
     155          ENDIF
    152156          IF ( humidity  .OR.  passive_scalar )  CALL exchange_horiz( q_p )
    153157          IF ( cloud_droplets )  THEN
  • palm/trunk/SOURCE/user_interface.f90

    r90 r95  
    44! Actual revisions:
    55! -----------------
    6 !
     6! user action for salinity added
    77!
    88! Former revisions:
     
    393393
    394394
     395          CASE ( 'sa-tendency' )
     396
     397
    395398          CASE ( 'e-tendency' )
    396399
     
    443446
    444447          CASE ( 'pt-tendency' )
     448
     449
     450          CASE ( 'sa-tendency' )
    445451
    446452
  • palm/trunk/SOURCE/write_3d_binary.f90

    r94 r95  
    44! Actual revisions:
    55! -----------------
    6 ! +sa, saswst
     6! +sa, saswsb, saswst
    77!
    88! Former revisions:
     
    134134    IF ( ocean )  THEN
    135135       WRITE ( 14 )  'sa                  ';  WRITE ( 14 )  sa
     136       WRITE ( 14 )  'saswsb              ';  WRITE ( 14 )  saswsb
    136137       WRITE ( 14 )  'saswst              ';  WRITE ( 14 )  saswst
    137138    ENDIF
  • palm/trunk/SOURCE/write_var_list.f90

    r94 r95  
    44! Actual revisions:
    55! -----------------
    6 ! +ocean, sa_init, sa_surface, sa_vertical_gradient,
    7 ! sa_vertical_gradient_level, top_salinity_flux
     6! +bc_sa_t, ocean, sa_init, sa_surface, sa_vertical_gradient,
     7! sa_vertical_gradient_level, bottom/top_salinity_flux
    88!
    99! Former revisions:
     
    107107    WRITE ( 14 )  'bc_s_t                        '
    108108    WRITE ( 14 )  bc_s_t
     109    WRITE ( 14 )  'bc_sa_t                       '
     110    WRITE ( 14 )  bc_sa_t
    109111    WRITE ( 14 )  'bc_uv_b                       '
    110112    WRITE ( 14 )  bc_uv_b
    111113    WRITE ( 14 )  'bc_uv_t                       '
    112114    WRITE ( 14 )  bc_uv_t
     115    WRITE ( 14 )  'bottom_salinityflux           '
     116    WRITE ( 14 )  bottom_salinityflux
    113117    WRITE ( 14 )  'building_height               '
    114118    WRITE ( 14 )  building_height
Note: See TracChangeset for help on using the changeset viewer.