Changeset 97 for palm/trunk/SOURCE
- Timestamp:
- Jun 21, 2007 8:23:15 AM (17 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 18 edited
-
CURRENT_MODIFICATIONS (modified) (2 diffs)
-
advec_s_bc.f90 (modified) (3 diffs)
-
buoyancy.f90 (modified) (8 diffs)
-
check_parameters.f90 (modified) (10 diffs)
-
diffusion_e.f90 (modified) (9 diffs)
-
diffusivities.f90 (modified) (5 diffs)
-
eqn_state_seawater.f90 (modified) (6 diffs)
-
flow_statistics.f90 (modified) (6 diffs)
-
header.f90 (modified) (11 diffs)
-
init_3d_model.f90 (modified) (1 diff)
-
init_ocean.f90 (modified) (2 diffs)
-
modules.f90 (modified) (6 diffs)
-
netcdf.f90 (modified) (5 diffs)
-
palm.f90 (modified) (1 diff)
-
production_e.f90 (modified) (3 diffs)
-
prognostic_equations.f90 (modified) (12 diffs)
-
run_control.f90 (modified) (2 diffs)
-
time_integration.f90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/CURRENT_MODIFICATIONS
r96 r97 1 1 New: 2 2 --- 3 ocean version including prognostic equation for salinity and equation of state for seawater4 + inipar-parameters bc_sa_t, bottom_salinityflux, ocean, sa_surface, sa_vertical_gradient, sa_vertical_gradient_level, top_salinityflux5 3 6 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, parin, prognostic_equations, read_var_list, sum_up_3d_data, swap_timelevel, time_integration, user_interface, write_var_list, write_3d_binary7 8 New:9 eqn_state_seawater, init_ocean10 4 11 5 Changed: 12 6 ------- 13 hydro_press renamed hyp, routine calc_mean_pt_profile renamed calc_mean_profile14 15 advec_particles, buoyancy, calc_liquid_water_content, init_cloud_physics, modules, prognostic_equations16 7 17 8 … … 19 10 ------ 20 11 21 Bugfix: height above topography instead of height above level k=0 is used for calculating the mixing length (diffusion_e and diffusivities).22 12 23 diffusion_e, prognostic_equations24 25 26 27 28 To be completed (adjusted to ocean version):29 ----------------30 31 boundary conditions for salinity (check_parameters)32 output of salinity profile (header)33 34 initialization of ug,vg profiles35 36 calculation of w*37 output of theta* and z_i (e.g. run_control file)38 39 default disturbance level40 41 Bott-Chlond scheme to be extended for salinity42 43 Upstream-spline scheme to be extended for salinity????? -
palm/trunk/SOURCE/advec_s_bc.f90
r77 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! Advection of salinity included 7 ! Bugfix: Error in boundary condition for TKE removed 7 8 ! 8 9 ! Former revisions: … … 755 756 ENDIF 756 757 758 ELSEIF ( sk_char == 'sa' ) THEN 759 760 ! 761 !-- Salinity boundary condition at the bottom boundary. 762 !-- So far, always Neumann (i.e. here zero gradient) is used 763 DO i = nxl, nxr 764 DO j = nys, nyn 765 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 766 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 767 ENDDO 768 ENDDO 769 770 ! 771 !-- Salinity boundary condition at the top boundary. 772 !-- Dirichlet or Neumann (zero gradient) 773 DO i = nxl, nxr 774 DO j = nys, nyn 775 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) 776 sk_p(nzt+3,j,i) = sk_p(nzt+1,j,i) 777 ENDDO 778 ENDDO 779 757 780 ELSEIF ( sk_char == 'q' ) THEN 758 781 759 782 ! 760 !-- Specific humidity boundary condition at the bottom boundary 761 IF ( ibc_q_b == 0 ) THEN 762 ! 763 !-- Dirichlet (fixed surface humidity) 764 DO i = nxl, nxr 765 DO j = nys, nyn 766 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 767 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 768 ENDDO 769 ENDDO 770 771 ELSE 772 ! 773 !-- Neumann (i.e. here zero gradient) 774 DO i = nxl, nxr 775 DO j = nys, nyn 776 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 777 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 778 ENDDO 779 ENDDO 780 781 ENDIF 783 !-- Specific humidity boundary condition at the bottom boundary. 784 !-- Dirichlet (fixed surface humidity) or Neumann (i.e. zero gradient) 785 DO i = nxl, nxr 786 DO j = nys, nyn 787 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 788 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 789 ENDDO 790 ENDDO 782 791 783 792 ! … … 809 818 ! 810 819 !-- TKE boundary condition at bottom and top boundary (generally Neumann) 811 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 812 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 813 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) 814 sk_p(nzt+3,j,i) = sk_p(nzt+1,j,i) 820 DO i = nxl, nxr 821 DO j = nys, nyn 822 sk_p(nzb-1,j,i) = sk_p(nzb,j,i) 823 sk_p(nzb-2,j,i) = sk_p(nzb,j,i) 824 sk_p(nzt+2,j,i) = sk_p(nzt+1,j,i) 825 sk_p(nzt+3,j,i) = sk_p(nzt+1,j,i) 826 ENDDO 827 ENDDO 815 828 816 829 ELSE -
palm/trunk/SOURCE/buoyancy.f90
r96 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Routine reneralized to be used with temperature AND density: 7 ! argument theta renamed var, new argument var_reference, 8 ! use_pt_reference renamed use_reference, 6 9 ! calc_mean_pt_profile renamed calc_mean_profile 7 10 ! … … 46 49 ! Call for all grid points 47 50 !------------------------------------------------------------------------------! 48 SUBROUTINE buoyancy( theta, wind_component, pr )51 SUBROUTINE buoyancy( var, var_reference, wind_component, pr ) 49 52 50 53 USE arrays_3d … … 57 60 58 61 INTEGER :: i, j, k, pr, wind_component 59 REAL, DIMENSION(:,:,:), POINTER :: theta 62 REAL :: var_reference 63 REAL, DIMENSION(:,:,:), POINTER :: var 60 64 61 65 … … 63 67 ! 64 68 !-- Normal case: horizontal surface 65 IF ( use_pt_reference ) THEN 66 DO i = nxl, nxr 67 DO j = nys, nyn 68 DO k = nzb_s_inner(j,i)+1, nzt-1 69 tend(k,j,i) = tend(k,j,i) + g * 0.5 * ( & 70 ( theta(k,j,i) - hom(k,1,pr,0) ) / pt_reference + & 71 ( theta(k+1,j,i) - hom(k+1,1,pr,0) ) / pt_reference & 69 IF ( use_reference ) THEN 70 DO i = nxl, nxr 71 DO j = nys, nyn 72 DO k = nzb_s_inner(j,i)+1, nzt-1 73 tend(k,j,i) = tend(k,j,i) + atmos_ocean_sign * g * 0.5 * & 74 ( & 75 ( var(k,j,i) - hom(k,1,pr,0) ) / var_reference + & 76 ( var(k+1,j,i) - hom(k+1,1,pr,0) ) / var_reference & 72 77 ) 73 78 ENDDO … … 78 83 DO j = nys, nyn 79 84 DO k = nzb_s_inner(j,i)+1, nzt-1 80 tend(k,j,i) = tend(k,j,i) + g * 0.5 * ( & 81 ( theta(k,j,i) - hom(k,1,pr,0) ) / hom(k,1,pr,0) + & 82 ( theta(k+1,j,i) - hom(k+1,1,pr,0) ) / hom(k+1,1,pr,0) & 85 tend(k,j,i) = tend(k,j,i) + atmos_ocean_sign * g * 0.5 * & 86 ( & 87 ( var(k,j,i) - hom(k,1,pr,0) ) / hom(k,1,pr,0) + & 88 ( var(k+1,j,i) - hom(k+1,1,pr,0) ) / hom(k+1,1,pr,0) & 83 89 ) 84 90 ENDDO … … 136 142 ! Call for grid point i,j 137 143 !------------------------------------------------------------------------------! 138 SUBROUTINE buoyancy_ij( i, j, theta, wind_component, pr )144 SUBROUTINE buoyancy_ij( i, j, var, var_reference, wind_component, pr ) 139 145 140 146 USE arrays_3d … … 147 153 148 154 INTEGER :: i, j, k, pr, wind_component 149 REAL, DIMENSION(:,:,:), POINTER :: theta 155 REAL :: var_reference 156 REAL, DIMENSION(:,:,:), POINTER :: var 150 157 151 158 … … 153 160 ! 154 161 !-- Normal case: horizontal surface 155 IF ( use_ pt_reference ) THEN156 DO k = nzb_s_inner(j,i)+1, nzt-1 157 tend(k,j,i) = tend(k,j,i) + g * 0.5 * (&158 ( theta(k,j,i) - hom(k,1,pr,0) ) / pt_reference + &159 ( theta(k+1,j,i) - hom(k+1,1,pr,0) ) / pt_reference &160 )161 ENDDO 162 ELSE 163 DO k = nzb_s_inner(j,i)+1, nzt-1 164 tend(k,j,i) = tend(k,j,i) + g * 0.5 * (&165 ( theta(k,j,i) - hom(k,1,pr,0) ) / hom(k,1,pr,0) + &166 ( theta(k+1,j,i) - hom(k+1,1,pr,0) ) / hom(k+1,1,pr,0) &167 )162 IF ( use_reference ) THEN 163 DO k = nzb_s_inner(j,i)+1, nzt-1 164 tend(k,j,i) = tend(k,j,i) + atmos_ocean_sign * g * 0.5 * ( & 165 ( var(k,j,i) - hom(k,1,pr,0) ) / var_reference + & 166 ( var(k+1,j,i) - hom(k+1,1,pr,0) ) / var_reference & 167 ) 168 ENDDO 169 ELSE 170 DO k = nzb_s_inner(j,i)+1, nzt-1 171 tend(k,j,i) = tend(k,j,i) + atmos_ocean_sign * g * 0.5 * ( & 172 ( var(k,j,i) - hom(k,1,pr,0) ) / hom(k,1,pr,0) + & 173 ( var(k+1,j,i) - hom(k+1,1,pr,0) ) / hom(k+1,1,pr,0) & 174 ) 168 175 ENDDO 169 176 ENDIF -
palm/trunk/SOURCE/check_parameters.f90
r96 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Initial salinity profile is calculated. 6 ! Initial salinity profile is calculated, salinity boundary conditions are 7 ! checked, 7 8 ! z_max_do1d is checked only in case of ocean = .f., 8 ! +initial temperature profile for the ocean version, 9 ! +initial temperature and geostrophic velocity profiles for the ocean version, 10 ! use_pt_reference renamed use_reference 9 11 ! 10 12 ! Former revisions: … … 163 165 WRITE( action, '(A,A)' ) 'timestep_scheme = ', timestep_scheme 164 166 ENDIF 167 IF ( momentum_advec == 'ups-scheme' ) THEN 168 WRITE( action, '(A,A)' ) 'momentum_advec = ', momentum_advec 169 ENDIF 165 170 IF ( action /= ' ' ) THEN 166 171 IF ( myid == 0 ) THEN … … 447 452 i = 1 448 453 gradient = 0.0 449 ug_vertical_gradient_level_ind(1) = 0 450 ug(0) = ug_surface 451 DO k = 1, nzt+1 452 IF ( ug_vertical_gradient_level(i) < zu(k) .AND. & 453 ug_vertical_gradient_level(i) >= 0.0 ) THEN 454 gradient = ug_vertical_gradient(i) / 100.0 455 ug_vertical_gradient_level_ind(i) = k - 1 456 i = i + 1 457 IF ( i > 10 ) THEN 458 IF ( myid == 0 ) THEN 459 PRINT*, '+++ check_parameters: upper bound 10 of array', & 460 ' "ug_vertical_gradient_level_ind" exceeded' 461 ENDIF 462 CALL local_stop 463 ENDIF 464 ENDIF 465 IF ( gradient /= 0.0 ) THEN 466 IF ( k /= 1 ) THEN 467 ug(k) = ug(k-1) + dzu(k) * gradient 454 455 IF ( .NOT. ocean ) THEN 456 457 ug_vertical_gradient_level_ind(1) = 0 458 ug(0) = ug_surface 459 DO k = 1, nzt+1 460 IF ( ug_vertical_gradient_level(i) < zu(k) .AND. & 461 ug_vertical_gradient_level(i) >= 0.0 ) THEN 462 gradient = ug_vertical_gradient(i) / 100.0 463 ug_vertical_gradient_level_ind(i) = k - 1 464 i = i + 1 465 IF ( i > 10 ) THEN 466 IF ( myid == 0 ) THEN 467 PRINT*, '+++ check_parameters: upper bound 10 of array', & 468 ' "ug_vertical_gradient_level_ind" exceeded' 469 ENDIF 470 CALL local_stop 471 ENDIF 472 ENDIF 473 IF ( gradient /= 0.0 ) THEN 474 IF ( k /= 1 ) THEN 475 ug(k) = ug(k-1) + dzu(k) * gradient 476 ELSE 477 ug(k) = ug_surface + 0.5 * dzu(k) * gradient 478 ENDIF 468 479 ELSE 469 ug(k) = ug_surface + 0.5 * dzu(k) * gradient 470 ENDIF 471 ELSE 472 ug(k) = ug(k-1) 473 ENDIF 474 ENDDO 480 ug(k) = ug(k-1) 481 ENDIF 482 ENDDO 483 484 ELSE 485 486 ug_vertical_gradient_level_ind(1) = nzt+1 487 DO k = nzt, 0, -1 488 IF ( ug_vertical_gradient_level(i) > zu(k) .AND. & 489 ug_vertical_gradient_level(i) <= 0.0 ) THEN 490 gradient = ug_vertical_gradient(i) / 100.0 491 ug_vertical_gradient_level_ind(i) = k + 1 492 i = i + 1 493 IF ( i > 10 ) THEN 494 IF ( myid == 0 ) THEN 495 PRINT*, '+++ check_parameters: upper bound 10 of array', & 496 ' "ug_vertical_gradient_level_ind" exceeded' 497 ENDIF 498 CALL local_stop 499 ENDIF 500 ENDIF 501 IF ( gradient /= 0.0 ) THEN 502 IF ( k /= nzt ) THEN 503 ug(k) = ug(k+1) - dzu(k+1) * gradient 504 ELSE 505 ug(k) = ug_surface - 0.5 * dzu(k+1) * gradient 506 ug(k+1) = ug_surface + 0.5 * dzu(k+1) * gradient 507 ENDIF 508 ELSE 509 ug(k) = ug(k+1) 510 ENDIF 511 ENDDO 512 513 ENDIF 475 514 476 515 u_init = ug … … 478 517 ! 479 518 !-- In case of no given gradients for ug, choose a vanishing gradient 480 IF ( ug_vertical_gradient_level(1) == - 1.0) THEN519 IF ( ug_vertical_gradient_level(1) == -9999999.9 ) THEN 481 520 ug_vertical_gradient_level(1) = 0.0 482 521 ENDIF … … 488 527 i = 1 489 528 gradient = 0.0 490 vg_vertical_gradient_level_ind(1) = 0 491 vg(0) = vg_surface 492 DO k = 1, nzt+1 493 IF ( vg_vertical_gradient_level(i) < zu(k) .AND. & 494 vg_vertical_gradient_level(i) >= 0.0 ) THEN 495 gradient = vg_vertical_gradient(i) / 100.0 496 vg_vertical_gradient_level_ind(i) = k - 1 497 i = i + 1 498 IF ( i > 10 ) THEN 499 IF ( myid == 0 ) THEN 500 PRINT*, '+++ check_parameters: upper bound 10 of array', & 501 ' "vg_vertical_gradient_level_ind" exceeded' 502 ENDIF 503 CALL local_stop 504 ENDIF 505 ENDIF 506 IF ( gradient /= 0.0 ) THEN 507 IF ( k /= 1 ) THEN 508 vg(k) = vg(k-1) + dzu(k) * gradient 529 530 IF ( .NOT. ocean ) THEN 531 532 vg_vertical_gradient_level_ind(1) = 0 533 vg(0) = vg_surface 534 DO k = 1, nzt+1 535 IF ( vg_vertical_gradient_level(i) < zu(k) .AND. & 536 vg_vertical_gradient_level(i) >= 0.0 ) THEN 537 gradient = vg_vertical_gradient(i) / 100.0 538 vg_vertical_gradient_level_ind(i) = k - 1 539 i = i + 1 540 IF ( i > 10 ) THEN 541 IF ( myid == 0 ) THEN 542 PRINT*, '+++ check_parameters: upper bound 10 of array', & 543 ' "vg_vertical_gradient_level_ind" exceeded' 544 ENDIF 545 CALL local_stop 546 ENDIF 547 ENDIF 548 IF ( gradient /= 0.0 ) THEN 549 IF ( k /= 1 ) THEN 550 vg(k) = vg(k-1) + dzu(k) * gradient 551 ELSE 552 vg(k) = vg_surface + 0.5 * dzu(k) * gradient 553 ENDIF 509 554 ELSE 510 vg(k) = vg_surface + 0.5 * dzu(k) * gradient 511 ENDIF 512 ELSE 513 vg(k) = vg(k-1) 514 ENDIF 515 ENDDO 555 vg(k) = vg(k-1) 556 ENDIF 557 ENDDO 558 559 ELSE 560 561 vg_vertical_gradient_level_ind(1) = 0 562 DO k = nzt, 0, -1 563 IF ( vg_vertical_gradient_level(i) > zu(k) .AND. & 564 vg_vertical_gradient_level(i) <= 0.0 ) THEN 565 gradient = vg_vertical_gradient(i) / 100.0 566 vg_vertical_gradient_level_ind(i) = k + 1 567 i = i + 1 568 IF ( i > 10 ) THEN 569 IF ( myid == 0 ) THEN 570 PRINT*, '+++ check_parameters: upper bound 10 of array', & 571 ' "vg_vertical_gradient_level_ind" exceeded' 572 ENDIF 573 CALL local_stop 574 ENDIF 575 ENDIF 576 IF ( gradient /= 0.0 ) THEN 577 IF ( k /= nzt ) THEN 578 vg(k) = vg(k+1) - dzu(k+1) * gradient 579 ELSE 580 vg(k) = vg_surface - 0.5 * dzu(k+1) * gradient 581 vg(k+1) = vg_surface + 0.5 * dzu(k+1) * gradient 582 ENDIF 583 ELSE 584 vg(k) = vg(k+1) 585 ENDIF 586 ENDDO 587 588 ENDIF 516 589 517 590 v_init = vg … … 519 592 ! 520 593 !-- In case of no given gradients for vg, choose a vanishing gradient 521 IF ( vg_vertical_gradient_level(1) == - 1.0) THEN594 IF ( vg_vertical_gradient_level(1) == -9999999.9 ) THEN 522 595 vg_vertical_gradient_level(1) = 0.0 523 596 ENDIF … … 707 780 708 781 ! 709 !-- Reference temperature to be used in buoyancy terms 710 IF ( pt_reference /= 9999999.9 ) use_pt_reference = .TRUE. 782 !-- Ocean runs always use reference values in the buoyancy term. Therefore 783 !-- set the reference temperature equal to the surface temperature. 784 IF ( ocean .AND. pt_reference == 9999999.9 ) pt_reference = pt_surface 785 786 ! 787 !-- Reference value has to be used in buoyancy terms 788 IF ( pt_reference /= 9999999.9 ) use_reference = .TRUE. 789 790 ! 791 !-- Sign of buoyancy/stability terms 792 IF ( ocean ) atmos_ocean_sign = -1.0 793 794 ! 795 !-- Ocean version is using flux boundary conditions at the top 796 IF ( ocean ) use_top_fluxes = .TRUE. 711 797 712 798 ! … … 1020 1106 1021 1107 IF ( top_salinityflux == 9999999.9 ) constant_top_salinityflux = .FALSE. 1108 IF ( ibc_sa_t == 1 .AND. top_salinityflux == 9999999.9 ) THEN 1109 IF ( myid == 0 ) THEN 1110 PRINT*, '+++ check_parameters:' 1111 PRINT*, ' boundary_condition: bc_sa_t = ', bc_sa_t 1112 PRINT*, ' requires to set top_salinityflux ' 1113 ENDIF 1114 CALL local_stop 1115 ENDIF 1022 1116 1023 1117 ! … … 2474 2568 ! 2475 2569 !-- Determine upper and lower hight level indices for random perturbations 2476 IF ( disturbance_level_b == -1.0 ) THEN 2477 disturbance_level_b = zu(nzb+3) 2478 disturbance_level_ind_b = nzb + 3 2570 IF ( disturbance_level_b == -9999999.9 ) THEN 2571 IF ( ocean ) THEN 2572 disturbance_level_b = zu((nzt*2)/3) 2573 disturbance_level_ind_b = ( nzt * 2 ) / 3 2574 ELSE 2575 disturbance_level_b = zu(nzb+3) 2576 disturbance_level_ind_b = nzb + 3 2577 ENDIF 2479 2578 ELSEIF ( disturbance_level_b < zu(3) ) THEN 2480 2579 IF ( myid == 0 ) PRINT*, '+++ check_parameters: disturbance_level_b=',& … … 2496 2595 ENDIF 2497 2596 2498 IF ( disturbance_level_t == -1.0 ) THEN 2499 disturbance_level_t = zu(nzt/3) 2500 disturbance_level_ind_t = nzt / 3 2597 IF ( disturbance_level_t == -9999999.9 ) THEN 2598 IF ( ocean ) THEN 2599 disturbance_level_t = zu(nzt-3) 2600 disturbance_level_ind_t = nzt - 3 2601 ELSE 2602 disturbance_level_t = zu(nzt/3) 2603 disturbance_level_ind_t = nzt / 3 2604 ENDIF 2501 2605 ELSEIF ( disturbance_level_t > zu(nzt-2) ) THEN 2502 2606 IF ( myid == 0 ) PRINT*, '+++ check_parameters: disturbance_level_t=',& -
palm/trunk/SOURCE/diffusion_e.f90
r94 r97 8 8 ! This is also a bugfix, because the height above the topography is now 9 9 ! used instead of the height above level k=0. 10 ! theta renamed var, dpt_dz renamed dvar_dz, +new argument var_reference 11 ! use_pt_reference renamed use_reference 10 12 ! 11 13 ! Former revisions: … … 50 52 ! Call for all grid points 51 53 !------------------------------------------------------------------------------! 52 SUBROUTINE diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, theta, &53 rif, tend, zu, zw )54 SUBROUTINE diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, var, & 55 var_reference, rif, tend, zu, zw ) 54 56 55 57 USE control_parameters … … 61 63 62 64 INTEGER :: i, j, k 63 REAL :: d pt_dz, l_stable, phi_m65 REAL :: dvar_dz, l_stable, phi_m, var_reference 64 66 REAL :: ddzu(1:nzt+1), dd2zu(1:nzt), ddzw(1:nzt+1), & 65 67 l_grid(1:nzt), zu(0:nzt+1), zw(0:nzt+1) 66 68 REAL, DIMENSION(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) :: diss, tend 67 69 REAL, DIMENSION(:,:), POINTER :: rif 68 REAL, DIMENSION(:,:,:), POINTER :: e, km, theta70 REAL, DIMENSION(:,:,:), POINTER :: e, km, var 69 71 REAL, DIMENSION(nzb+1:nzt,nys:nyn) :: dissipation, l, ll 70 72 … … 73 75 !-- This if clause must be outside the k-loop because otherwise 74 76 !-- runtime errors occur with -C hopt on NEC 75 IF ( use_ pt_reference ) THEN77 IF ( use_reference ) THEN 76 78 77 79 DO i = nxl, nxr … … 91 93 ! 92 94 !-- Calculate the mixing length (for dissipation) 93 dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k) 94 IF ( dpt_dz > 0.0 ) THEN 95 dvar_dz = atmos_ocean_sign * & 96 ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k) 97 IF ( dvar_dz > 0.0 ) THEN 95 98 l_stable = 0.76 * SQRT( e(k,j,i) ) / & 96 SQRT( g / pt_reference * dpt_dz ) + 1E-599 SQRT( g / var_reference * dvar_dz ) + 1E-5 97 100 ELSE 98 101 l_stable = l_grid(k) … … 180 183 ! 181 184 !-- Calculate the mixing length (for dissipation) 182 dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k) 183 IF ( dpt_dz > 0.0 ) THEN 185 dvar_dz = atmos_ocean_sign * & 186 ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k) 187 IF ( dvar_dz > 0.0 ) THEN 184 188 l_stable = 0.76 * SQRT( e(k,j,i) ) / & 185 SQRT( g / theta(k,j,i) * dpt_dz ) + 1E-5189 SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5 186 190 ELSE 187 191 l_stable = l_grid(k) … … 270 274 !------------------------------------------------------------------------------! 271 275 SUBROUTINE diffusion_e_ij( i, j, ddzu, dd2zu, ddzw, diss, e, km, l_grid, & 272 theta, rif, tend, zu, zw )276 var, var_reference, rif, tend, zu, zw ) 273 277 274 278 USE control_parameters … … 280 284 281 285 INTEGER :: i, j, k 282 REAL :: d pt_dz, l_stable, phi_m286 REAL :: dvar_dz, l_stable, phi_m, var_reference 283 287 REAL :: ddzu(1:nzt+1), dd2zu(1:nzt), ddzw(1:nzt+1), & 284 288 l_grid(1:nzt), zu(0:nzt+1), zw(0:nzt+1) 285 289 REAL, DIMENSION(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) :: diss, tend 286 290 REAL, DIMENSION(:,:), POINTER :: rif 287 REAL, DIMENSION(:,:,:), POINTER :: e, km, theta291 REAL, DIMENSION(:,:,:), POINTER :: e, km, var 288 292 REAL, DIMENSION(nzb+1:nzt) :: dissipation, l, ll 289 293 … … 303 307 !-- Calculate the mixing length (for dissipation) 304 308 DO k = nzb_s_inner(j,i)+1, nzt 305 dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k) 306 IF ( dpt_dz > 0.0 ) THEN 307 IF ( use_pt_reference ) THEN 309 dvar_dz = atmos_ocean_sign * & 310 ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k) 311 IF ( dvar_dz > 0.0 ) THEN 312 IF ( use_reference ) THEN 308 313 l_stable = 0.76 * SQRT( e(k,j,i) ) / & 309 SQRT( g / pt_reference * dpt_dz ) + 1E-5314 SQRT( g / var_reference * dvar_dz ) + 1E-5 310 315 ELSE 311 316 l_stable = 0.76 * SQRT( e(k,j,i) ) / & 312 SQRT( g / theta(k,j,i) * dpt_dz ) + 1E-5317 SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5 313 318 ENDIF 314 319 ELSE -
palm/trunk/SOURCE/diffusivities.f90
r94 r97 1 SUBROUTINE diffusivities( theta)1 SUBROUTINE diffusivities( var, var_reference ) 2 2 3 3 !------------------------------------------------------------------------------! … … 7 7 ! This is also a bugfix, because the height above the topography is now 8 8 ! used instead of the height above level k=0. 9 ! theta renamed var, dpt_dz renamed dvar_dz, +new argument var_reference 10 ! use_pt_reference renamed use_reference 9 11 ! 10 12 ! Former revisions: … … 41 43 INTEGER :: i, j, k, omp_get_thread_num, sr, tn 42 44 43 REAL :: d pt_dz, l_stable, phi_m = 1.045 REAL :: dvar_dz, l_stable, var_reference 44 46 45 REAL :: theta(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) 47 REAL, SAVE :: phi_m = 1.0 48 49 REAL :: var(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1) 46 50 47 51 REAL, DIMENSION(1:nzt) :: l, ll, sqrt_e … … 58 62 ! 59 63 !-- Compute the turbulent diffusion coefficient for momentum 60 !$OMP PARALLEL PRIVATE (d pt_dz,i,j,k,l,ll,l_stable,phi_m,sqrt_e,sr,tn)64 !$OMP PARALLEL PRIVATE (dvar_dz,i,j,k,l,ll,l_stable,phi_m,sqrt_e,sr,tn) 61 65 !$ tn = omp_get_thread_num() 62 66 … … 94 98 !-- Determine the mixing length 95 99 DO k = nzb_s_inner(j,i)+1, nzt 96 dpt_dz = ( theta(k+1,j,i) - theta(k-1,j,i) ) * dd2zu(k) 97 IF ( dpt_dz > 0.0 ) THEN 98 IF ( use_pt_reference ) THEN 100 dvar_dz = atmos_ocean_sign * & ! inverse effect of pt/rho gradient 101 ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k) 102 IF ( dvar_dz > 0.0 ) THEN 103 IF ( use_reference ) THEN 99 104 l_stable = 0.76 * sqrt_e(k) / & 100 SQRT( g / pt_reference * dpt_dz ) + 1E-5105 SQRT( g / var_reference * dvar_dz ) + 1E-5 101 106 ELSE 102 107 l_stable = 0.76 * sqrt_e(k) / & 103 SQRT( g / theta(k,j,i) * dpt_dz ) + 1E-5108 SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5 104 109 ENDIF 105 110 ELSE -
palm/trunk/SOURCE/eqn_state_seawater.f90
r96 r97 8 8 ! Former revisions: 9 9 ! ----------------- 10 ! $Id :$10 ! $Id$ 11 11 ! 12 12 ! Initial revision … … 18 18 ! salinity, and pressure. 19 19 ! For coefficients see Jackett et al., 2006: J. Atm. Ocean Tech. 20 ! eqn_state_seawater calculates the potential density referred at hyp(0). 21 ! eqn_state_seawater_func calculates density. 20 22 !------------------------------------------------------------------------------! 21 23 … … 70 72 DO i = nxl, nxr 71 73 DO j = nys, nyn 72 DO k = nzb_ u_inner(j,i)+1, nzt74 DO k = nzb_s_inner(j,i)+1, nzt 73 75 ! 74 76 !-- Pressure is needed in dbar 75 p1 = hyp( k) * 1E-477 p1 = hyp(0) * 1E-4 76 78 p2 = p1 * p1 77 79 p3 = p2 * p1 … … 102 104 103 105 ENDDO 106 ! 107 !-- Neumann conditions are assumed at bottom and top boundary 108 rho(nzt+1,j,i) = rho(nzt,j,i) 109 rho(nzb_s_inner(j,i),j,i) = rho(nzb_s_inner(j,i)+1,j,i) 104 110 ENDDO 105 111 ENDDO … … 122 128 REAL :: p1, p2, p3, pt1, pt2, pt3, pt4, sa1, sa15, sa2 123 129 124 DO k = nzb_ u_inner(j,i)+1, nzt130 DO k = nzb_s_inner(j,i)+1, nzt 125 131 ! 126 132 !-- Pressure is needed in dbar 127 p1 = hyp( k) * 1E-4133 p1 = hyp(0) * 1E-4 128 134 p2 = p1 * p1 129 135 p3 = p2 * p1 … … 152 158 ) 153 159 ENDDO 160 ! 161 !-- Neumann conditions are assumed at bottom and top boundary 162 rho(nzt+1,j,i) = rho(nzt,j,i) 163 rho(nzb_s_inner(j,i),j,i) = rho(nzb_s_inner(j,i)+1,j,i) 154 164 155 165 END SUBROUTINE eqn_state_seawater_ij -
palm/trunk/SOURCE/flow_statistics.f90
r96 r97 5 5 ! ----------------- 6 6 ! Statistics for ocean version (salinity, density) added, 7 ! calculation of Deardorff velocity scale adjusted to be used with the ocean8 ! version (HAS STILL TO BE COMPLETED!!!)7 ! calculation of z_i and Deardorff velocity scale adjusted to be used with 8 ! the ocean version 9 9 ! 10 10 ! Former revisions: … … 558 558 ! 559 559 !-- Salinity flux and density (density does not belong to here, 560 !-- but so far there is no suitable place to calculate)560 !-- but so far there is no other suitable place to calculate) 561 561 IF ( ocean ) THEN 562 562 pts = 0.5 * ( sa(k,j,i) - hom(k,1,23,sr) + & … … 611 611 612 612 ! 613 !-- Density at top follows Neumann condition 614 IF ( ocean ) sums_l(nzt+1,64,tn) = sums_l(nzt,64,tn) 615 616 ! 613 617 !-- Divergence of vertical flux of resolved scale energy and pressure 614 618 !-- fluctuations. First calculate the products, then the divergence. … … 709 713 ! 710 714 !-- Fluxes at the surface must be zero (e.g. due to the Prandtl-layer) 711 sums (nzb,58) = 0.0712 sums (nzb,59) = 0.0713 sums (nzb,60) = 0.0714 sums (nzb,61) = 0.0715 sums (nzb,62) = 0.0716 sums (nzb,63) = 0.0715 sums_l(nzb,58,tn) = 0.0 716 sums_l(nzb,59,tn) = 0.0 717 sums_l(nzb,60,tn) = 0.0 718 sums_l(nzb,61,tn) = 0.0 719 sums_l(nzb,62,tn) = 0.0 720 sums_l(nzb,63,tn) = 0.0 717 721 718 722 ENDIF … … 857 861 z_i(1) = 0.0 858 862 first = .TRUE. 859 ! IF ( .NOT. ocean ) THEN 863 IF ( ocean ) THEN 864 DO k = nzt, nzb+1, -1 865 IF ( first .AND. hom(k,1,18,sr) < 0.0 ) THEN 866 first = .FALSE. 867 height = zw(k) 868 ENDIF 869 IF ( hom(k,1,18,sr) < 0.0 .AND. & 870 hom(k-1,1,18,sr) > hom(k,1,18,sr) ) THEN 871 IF ( zw(k) < 1.5 * height ) THEN 872 z_i(1) = zw(k) 873 ELSE 874 z_i(1) = height 875 ENDIF 876 EXIT 877 ENDIF 878 ENDDO 879 ELSE 860 880 DO k = nzb, nzt-1 861 881 IF ( first .AND. hom(k,1,18,sr) < 0.0 ) THEN … … 873 893 ENDIF 874 894 ENDDO 875 ! ELSE 876 ! ENDIF 877 878 ! 879 !-- Second scheme: Starting from the top model boundary, look for the first 880 !-- characteristic kink in the temperature profile, where the originally 881 !-- stable stratification notably weakens. 895 ENDIF 896 897 ! 898 !-- Second scheme: Starting from the top/bottom model boundary, look for 899 !-- the first characteristic kink in the temperature profile, where the 900 !-- originally stable stratification notably weakens. 882 901 z_i(2) = 0.0 883 DO k = nzt-1, nzb+1, -1 884 IF ( ( hom(k+1,1,4,sr) - hom(k,1,4,sr) ) > & 885 2.0 * ( hom(k,1,4,sr) - hom(k-1,1,4,sr) ) ) THEN 886 z_i(2) = zu(k) 887 EXIT 888 ENDIF 889 ENDDO 902 IF ( ocean ) THEN 903 DO k = nzb+1, nzt-1 904 IF ( ( hom(k,1,4,sr) - hom(k-1,1,4,sr) ) > & 905 2.0 * ( hom(k+1,1,4,sr) - hom(k,1,4,sr) ) ) THEN 906 z_i(2) = zu(k) 907 EXIT 908 ENDIF 909 ENDDO 910 ELSE 911 DO k = nzt-1, nzb+1, -1 912 IF ( ( hom(k+1,1,4,sr) - hom(k,1,4,sr) ) > & 913 2.0 * ( hom(k,1,4,sr) - hom(k-1,1,4,sr) ) ) THEN 914 z_i(2) = zu(k) 915 EXIT 916 ENDIF 917 ENDDO 918 ENDIF 890 919 891 920 hom(nzb+6,1,pr_palm,sr) = z_i(1) -
palm/trunk/SOURCE/header.f90
r94 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Output of model height, stretch level, Prandtl-layer height and initial7 ! temperature profile adjusted to be used also with the ocean version.6 ! Adjustments for the ocean version. 7 ! use_pt_reference renamed use_reference 8 8 ! 9 9 ! Former revisions: … … 103 103 PRINT*,'+++ header: unknown action(s): ',initializing_actions 104 104 ENDIF 105 ENDIF 106 IF ( ocean ) THEN 107 run_classification = 'ocean - ' // run_classification 108 ELSE 109 run_classification = 'atmosphere - ' // run_classification 105 110 ENDIF 106 111 … … 391 396 roben = 'e(nzt+1) = e(nzt) = e(nzt-1)' 392 397 393 WRITE ( io, 301 ) runten, roben 394 395 ENDIF 396 397 IF ( humidity .OR. passive_scalar ) THEN 398 IF ( humidity ) THEN 399 IF ( ibc_q_b == 0 ) THEN 400 runten = 'q(0) = q_surface' 401 ELSE 402 runten = 'q(0) = q(1)' 403 ENDIF 404 IF ( ibc_q_t == 0 ) THEN 405 roben = 'q(nzt) = q_top' 406 ELSE 407 roben = 'q(nzt) = q(nzt-1) + dq/dz' 408 ENDIF 409 ELSE 410 IF ( ibc_q_b == 0 ) THEN 411 runten = 's(0) = s_surface' 412 ELSE 413 runten = 's(0) = s(1)' 414 ENDIF 415 IF ( ibc_q_t == 0 ) THEN 416 roben = 's(nzt) = s_top' 417 ELSE 418 roben = 's(nzt) = s(nzt-1) + ds/dz' 419 ENDIF 420 ENDIF 421 422 WRITE ( io, 302 ) runten, roben 423 398 WRITE ( io, 301 ) 'e', runten, roben 399 400 ENDIF 401 402 IF ( ocean ) THEN 403 runten = 'sa(0) = sa(1)' 404 IF ( ibc_sa_t == 0 ) THEN 405 roben = 'sa(nzt+1) = sa_surface' 406 ELSE 407 roben = 'sa(nzt+1) = sa(nzt)' 408 ENDIF 409 WRITE ( io, 301 ) 'sa', runten, roben 410 ENDIF 411 412 IF ( humidity ) THEN 413 IF ( ibc_q_b == 0 ) THEN 414 runten = 'q(0) = q_surface' 415 ELSE 416 runten = 'q(0) = q(1)' 417 ENDIF 418 IF ( ibc_q_t == 0 ) THEN 419 roben = 'q(nzt) = q_top' 420 ELSE 421 roben = 'q(nzt) = q(nzt-1) + dq/dz' 422 ENDIF 423 WRITE ( io, 301 ) 'q', runten, roben 424 ENDIF 425 426 IF ( passive_scalar ) THEN 427 IF ( ibc_q_b == 0 ) THEN 428 runten = 's(0) = s_surface' 429 ELSE 430 runten = 's(0) = s(1)' 431 ENDIF 432 IF ( ibc_q_t == 0 ) THEN 433 roben = 's(nzt) = s_top' 434 ELSE 435 roben = 's(nzt) = s(nzt-1) + ds/dz' 436 ENDIF 437 WRITE ( io, 301 ) 's', runten, roben 424 438 ENDIF 425 439 … … 442 456 IF ( constant_top_heatflux ) THEN 443 457 WRITE ( io, 306 ) top_heatflux 458 ENDIF 459 IF ( ocean .AND. constant_top_salinityflux ) THEN 460 WRITE ( io, 309 ) top_salinityflux 444 461 ENDIF 445 462 IF ( humidity .OR. passive_scalar ) THEN … … 887 904 !-- Other quantities 888 905 WRITE ( io, 411 ) g 889 IF ( use_pt_reference ) WRITE ( io, 412 ) pt_reference 906 IF ( use_reference ) THEN 907 IF ( ocean ) THEN 908 WRITE ( io, 412 ) prho_reference 909 ELSE 910 WRITE ( io, 413 ) pt_reference 911 ENDIF 912 ENDIF 890 913 891 914 ! … … 1011 1034 TRIM( gradients ), TRIM( slices ) 1012 1035 ENDIF 1036 ENDIF 1037 1038 ! 1039 !-- Initial salinity profile 1040 !-- Building output strings, starting with surface salinity 1041 IF ( ocean ) THEN 1042 WRITE ( temperatures, '(F6.2)' ) sa_surface 1043 gradients = '------' 1044 slices = ' 0' 1045 coordinates = ' 0.0' 1046 i = 1 1047 DO WHILE ( sa_vertical_gradient_level_ind(i) /= -9999 ) 1048 1049 WRITE (coor_chr,'(F7.2)') sa_init(sa_vertical_gradient_level_ind(i)) 1050 temperatures = TRIM( temperatures ) // ' ' // TRIM( coor_chr ) 1051 1052 WRITE (coor_chr,'(F7.2)') sa_vertical_gradient(i) 1053 gradients = TRIM( gradients ) // ' ' // TRIM( coor_chr ) 1054 1055 WRITE (coor_chr,'(I7)') sa_vertical_gradient_level_ind(i) 1056 slices = TRIM( slices ) // ' ' // TRIM( coor_chr ) 1057 1058 WRITE (coor_chr,'(F7.1)') sa_vertical_gradient_level(i) 1059 coordinates = TRIM( coordinates ) // ' ' // TRIM( coor_chr ) 1060 1061 i = i + 1 1062 ENDDO 1063 1064 WRITE ( io, 425 ) TRIM( coordinates ), TRIM( temperatures ), & 1065 TRIM( gradients ), TRIM( slices ) 1013 1066 ENDIF 1014 1067 … … 1135 1188 1136 1189 99 FORMAT (1X,78('-')) 1137 100 FORMAT (/1X,'*************************',11X, 28('-')/ &1190 100 FORMAT (/1X,'*************************',11X,42('-')/ & 1138 1191 1X,'* ',A,' *',11X,A/ & 1139 1X,'*************************',11X, 28('-')// &1192 1X,'*************************',11X,42('-')// & 1140 1193 ' Date: ',A8,11X,'Run: ',A20/ & 1141 1194 ' Time: ',A8,11X,'Run-No.: ',I2.2/ & … … 1256 1309 ' B. bound.: ',A/ & 1257 1310 ' T. bound.: ',A) 1258 301 FORMAT (/' e'// & 1259 ' B. bound.: ',A/ & 1260 ' T. bound.: ',A) 1261 302 FORMAT (/' q'// & 1311 301 FORMAT (/' ',A// & 1262 1312 ' B. bound.: ',A/ & 1263 1313 ' T. bound.: ',A) … … 1268 1318 ' zp = ',F6.2,' m z0 = ',F6.4,' m kappa = ',F4.2/ & 1269 1319 ' Rif value range: ',F6.2,' <= rif <=',F6.2) 1270 306 FORMAT (' Predefined constant heatflux: ',F 6.3,' K m/s')1320 306 FORMAT (' Predefined constant heatflux: ',F9.6,' K m/s') 1271 1321 307 FORMAT (' Heatflux has a random normal distribution') 1272 1322 308 FORMAT (' Predefined surface temperature') 1323 309 FORMAT (' Predefined constant salinityflux: ',F9.6,' psu m/s') 1273 1324 310 FORMAT (//' 1D-Model:'// & 1274 1325 ' Rif value range: ',F6.2,' <= rif <=',F6.2) … … 1352 1403 ' f* = ',F9.6,' 1/s') 1353 1404 411 FORMAT (/' Gravity : g = ',F4.1,' m/s**2') 1354 412 FORMAT (/' Reference temperature in buoyancy terms: ',F8.4,' K') 1405 412 FORMAT (/' Reference density in buoyancy terms: ',F8.3,' kg/m**3') 1406 413 FORMAT (/' Reference temperature in buoyancy terms: ',F8.4,' K') 1355 1407 415 FORMAT (/' Cloud physics parameters:'/ & 1356 1408 ' ------------------------'/) … … 1382 1434 424 FORMAT (/' Characteristic levels of the geo. wind component vg:'// & 1383 1435 ' Height: ',A,' m'/ & 1384 ' vg: ',A,' m/ S'/ &1436 ' vg: ',A,' m/s'/ & 1385 1437 ' Gradient: ',A,' 1/100s'/ & 1386 1438 ' Gridpoint: ',A) 1439 425 FORMAT (/' Characteristic levels of the initial salinity profile:'// & 1440 ' Height: ',A,' m'/ & 1441 ' Salinity: ',A,' psu'/ & 1442 ' Gradient: ',A,' psu/100m'/ & 1443 ' Gridpoint: ',A) 1387 1444 450 FORMAT (//' LES / Turbulence quantities:'/ & 1388 1445 ' ---------------------------'/) -
palm/trunk/SOURCE/init_3d_model.f90
r96 r97 529 529 hom(:,1,24,:) = SPREAD( kh(:,nys,nxl), 2, statistic_regions+1 ) 530 530 531 IF ( ocean ) THEN 532 ! 533 !-- Store initial salinity profile 534 hom(:,1,26,:) = SPREAD( sa(:,nys,nxl), 2, statistic_regions+1 ) 535 ENDIF 531 536 532 537 IF ( humidity ) THEN -
palm/trunk/SOURCE/init_ocean.f90
r96 r97 40 40 hyp(nzt+1) = surface_pressure * 100.0 41 41 42 hyp(nzt) = hyp(nzt+1) + rho_surface * g * 0.5 * dzu(nzt+1)43 rho_ref 42 hyp(nzt) = hyp(nzt+1) + rho_surface * g * 0.5 * dzu(nzt+1) 43 rho_reference = rho_surface * 0.5 * dzu(nzt+1) 44 44 45 45 DO k = nzt-1, 0, -1 … … 48 48 pt_l = 0.5 * ( pt_init(k) + pt_init(k+1) ) 49 49 50 rho_l 50 rho_l = eqn_state_seawater_func( hyp(k+1), pt_l, sa_l ) 51 51 52 hyp(k) = hyp(k+1) + rho_l * g * dzu(k+1)53 rho_ref = rho_ref+ rho_l * dzu(k+1)52 hyp(k) = hyp(k+1) + rho_l * g * dzu(k+1) 53 rho_reference = rho_reference + rho_l * dzu(k+1) 54 54 55 55 ENDDO 56 56 57 rho_ref = rho_ref / ( zw(nzt) - zu(nzb) ) 58 print*, '*** rho_ref = ', rho_ref 57 rho_reference = rho_reference / ( zw(nzt) - zu(nzb) ) 58 59 ! 60 !-- Calculate the reference potential density 61 prho_reference = 0.0 62 DO k = 0, nzt 63 64 sa_l = 0.5 * ( sa_init(k) + sa_init(k+1) ) 65 pt_l = 0.5 * ( pt_init(k) + pt_init(k+1) ) 66 67 prho_reference = prho_reference + dzu(k+1) * & 68 eqn_state_seawater_func( hyp(0), pt_l, sa_l ) 69 70 ENDDO 71 72 prho_reference = prho_reference / ( zu(nzt) - zu(nzb) ) 59 73 60 74 -
palm/trunk/SOURCE/modules.f90
r96 r97 5 5 ! Actual revisions: 6 6 ! ----------------- 7 ! + ocean, r, + salinity variables7 ! +atmos_ocean_sign, ocean, r, + salinity variables 8 8 ! defaults of .._vertical_gradient_levels changed from -1.0 to -9999999.9 9 ! hydro_press renamed hyp 9 ! hydro_press renamed hyp, use_pt_reference renamed use_reference 10 10 ! 11 11 ! Former revisions: … … 324 324 sloping_surface = .FALSE., stop_dt = .FALSE., & 325 325 terminate_run = .FALSE., use_prior_plot1d_parameters = .FALSE.,& 326 use_ pt_reference = .FALSE., use_surface_fluxes = .FALSE., &326 use_reference = .FALSE., use_surface_fluxes = .FALSE., & 327 327 use_top_fluxes = .FALSE., use_ug_for_galilei_tr = .TRUE., & 328 328 use_upstream_for_tke = .FALSE., wall_adjustment = .TRUE. … … 333 333 REAL :: advected_distance_x = 0.0, advected_distance_y = 0.0, & 334 334 alpha_surface = 0.0, asselin_filter_factor = 0.1, & 335 atmos_ocean_sign = 1.0, & 335 336 averaging_interval = 0.0, averaging_interval_pr = 9999999.9, & 336 337 averaging_interval_sp = 9999999.9, bc_pt_t_val, bc_q_t_val, & … … 340 341 building_wall_south = 9999999.9, cfl_factor = -1.0, & 341 342 cos_alpha_surface, disturbance_amplitude = 0.25, & 342 disturbance_energy_limit = 0.01, disturbance_level_b = -1.0, & 343 disturbance_level_t = -1.0, dt = -1.0, dt_averaging_input = 0.0, & 343 disturbance_energy_limit = 0.01, & 344 disturbance_level_b = -9999999.9, & 345 disturbance_level_t = -9999999.9, & 346 dt = -1.0, dt_averaging_input = 0.0, & 344 347 dt_averaging_input_pr = 9999999.9, dt_data_output = 9999999.9, & 345 348 dt_data_output_av = 9999999.9, dt_disturb = 9999999.9, & … … 360 363 overshoot_limit_w = 0.0, particle_maximum_age = 9999999.9, & 361 364 phi = 55.0, prandtl_number = 1.0, & 362 precipitation_amount_interval = 9999999.9, &363 pt_reference = 9999999.9, &364 pt_s lope_offset = 0.0, pt_surface = 300.0, &365 pt_surface_initial_change = 0.0, q_surface = 0.0, &366 q_surface_initial_change = 0.0, rayleigh_damping_factor= -1.0, &367 r ayleigh_damping_height = -1.0, residual_limit = 1.0E-4, &368 r estart_time = 9999999.9, rho_ref, rho_surface, rif_max = 1.0, &365 precipitation_amount_interval = 9999999.9, prho_reference, & 366 pt_reference = 9999999.9, pt_slope_offset = 0.0, & 367 pt_surface = 300.0, pt_surface_initial_change = 0.0, & 368 q_surface = 0.0, q_surface_initial_change = 0.0, & 369 rayleigh_damping_factor = -1.0, rayleigh_damping_height = -1.0, & 370 residual_limit = 1.0E-4, restart_time = 9999999.9, rho_reference, & 371 rho_surface, rif_max = 1.0, & 369 372 rif_min = -5.0, roughness_length = 0.1, sa_surface = 35.0, & 370 373 simulated_time = 0.0, simulated_time_at_begin, sin_alpha_surface, & … … 401 404 tsc(10) = (/ 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 /), & 402 405 ug_vertical_gradient(10) = 0.0, & 403 ug_vertical_gradient_level(10) = - 1.0, &406 ug_vertical_gradient_level(10) = -9999999.9, & 404 407 vg_vertical_gradient(10) = 0.0, & 405 vg_vertical_gradient_level(10) = - 1.0, &408 vg_vertical_gradient_level(10) = -9999999.9, & 406 409 volume_flow(1:2) = 0.0, volume_flow_area(1:2) = 0.0, & 407 410 volume_flow_initial(1:2) = 0.0, wall_heatflux(0:4) = 0.0 -
palm/trunk/SOURCE/netcdf.f90
r77 r97 7 7 ! Current revisions: 8 8 ! ------------------ 9 ! 9 ! Grids defined for rho and sa 10 10 ! 11 11 ! Former revisions: … … 345 345 !-- Most variables are defined on the scalar grid 346 346 CASE ( 'e', 'p', 'pc', 'pr', 'pt', 'q', 'ql', 'ql_c', 'ql_v', & 347 'ql_vp', 'qv', ' s', 'vpt' )347 'ql_vp', 'qv', 'rho', 's', 'sa', 'vpt' ) 348 348 349 349 grid_x = 'x' … … 870 870 CASE ( 'e_xy', 'p_xy', 'pc_xy', 'pr_xy', 'pt_xy', 'q_xy',& 871 871 'ql_xy', 'ql_c_xy', 'ql_v_xy', 'ql_vp_xy', & 872 'qv_xy', ' s_xy', 'vpt_xy' )872 'qv_xy', 'rho_xy', 's_xy', 'sa_xy', 'vpt_xy' ) 873 873 874 874 grid_x = 'x' … … 1414 1414 CASE ( 'e_xz', 'p_xz', 'pc_xz', 'pr_xz', 'pt_xz', 'q_xz', & 1415 1415 'ql_xz', 'ql_c_xz', 'ql_v_xz', 'ql_vp_xz', 'qv_xz', & 1416 ' s_xz', 'vpt_xz' )1416 'rho_xz', 's_xz', 'sa_xz', 'vpt_xz' ) 1417 1417 1418 1418 grid_x = 'x' … … 1925 1925 CASE ( 'e_yz', 'p_yz', 'pc_yz', 'pr_yz', 'pt_yz', 'q_yz', & 1926 1926 'ql_yz', 'ql_c_yz', 'ql_v_yz', 'ql_vp_yz', 'qv_yz', & 1927 ' s_yz', 'vpt_yz' )1927 'rho_yz', 's_yz', 'sa_yz', 'vpt_yz' ) 1928 1928 1929 1929 grid_x = 'x' -
palm/trunk/SOURCE/palm.f90
r90 r97 63 63 INTEGER :: i, run_description_header_i(80) 64 64 65 version = 'PALM 3. 2b'65 version = 'PALM 3.3' 66 66 67 67 #if defined( __parallel ) -
palm/trunk/SOURCE/production_e.f90
r77 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! energy production by density flux (in ocean) added 7 ! use_pt_reference renamed use_reference 7 8 ! 8 9 ! Former revisions: … … 364 365 IF ( .NOT. humidity ) THEN 365 366 366 IF ( use_pt_reference ) THEN 367 368 DO j = nys, nyn 369 DO k = nzb_diff_s_inner(j,i), nzt_diff 370 tend(k,j,i) = tend(k,j,i) - kh(k,j,i) * g/pt_reference * & 367 IF ( use_reference ) THEN 368 369 IF ( ocean ) THEN 370 ! 371 !-- So far in the ocean no special treatment of density flux in 372 !-- the bottom and top surface layer 373 DO j = nys, nyn 374 DO k = nzb_s_inner(j,i), nzt 375 tend(k,j,i) = tend(k,j,i) + & 376 kh(k,j,i) * g / prho_reference * & 377 ( rho(k+1,j,i)-rho(k-1,j,i) ) * dd2zu(k) 378 ENDDO 379 ENDDO 380 381 ELSE 382 383 DO j = nys, nyn 384 DO k = nzb_diff_s_inner(j,i), nzt_diff 385 tend(k,j,i) = tend(k,j,i) + & 386 kh(k,j,i) * g / pt_reference * & 371 387 ( pt(k+1,j,i) - pt(k-1,j,i) ) * dd2zu(k) 388 ENDDO 389 390 IF ( use_surface_fluxes ) THEN 391 k = nzb_diff_s_inner(j,i)-1 392 tend(k,j,i) = tend(k,j,i) + g / pt_reference * shf(j,i) 393 ENDIF 394 395 IF ( use_top_fluxes ) THEN 396 k = nzt 397 tend(k,j,i) = tend(k,j,i) + g / pt_reference * & 398 tswst(j,i) 399 ENDIF 372 400 ENDDO 373 401 374 IF ( use_surface_fluxes ) THEN 375 k = nzb_diff_s_inner(j,i)-1 376 tend(k,j,i) = tend(k,j,i) + g / pt_reference * shf(j,i) 377 ENDIF 378 379 IF ( use_top_fluxes ) THEN 380 k = nzt 381 tend(k,j,i) = tend(k,j,i) + g / pt_reference * tswst(j,i) 382 ENDIF 383 ENDDO 402 ENDIF 384 403 385 404 ELSE 386 405 387 DO j = nys, nyn 388 DO k = nzb_diff_s_inner(j,i), nzt_diff 389 tend(k,j,i) = tend(k,j,i) - kh(k,j,i) * g / pt(k,j,i) * & 406 IF ( ocean ) THEN 407 ! 408 !-- So far in the ocean no special treatment of density flux in 409 !-- the bottom and top surface layer 410 DO j = nys, nyn 411 DO k = nzb_s_inner(j,i), nzt 412 tend(k,j,i) = tend(k,j,i) - & 413 kh(k,j,i) * g / rho(k,j,i) * & 414 ( rho(k+1,j,i)-rho(k-1,j,i) ) * dd2zu(k) 415 ENDDO 416 ENDDO 417 418 ELSE 419 420 DO j = nys, nyn 421 DO k = nzb_diff_s_inner(j,i), nzt_diff 422 tend(k,j,i) = tend(k,j,i) - & 423 kh(k,j,i) * g / pt(k,j,i) * & 390 424 ( pt(k+1,j,i) - pt(k-1,j,i) ) * dd2zu(k) 425 ENDDO 426 427 IF ( use_surface_fluxes ) THEN 428 k = nzb_diff_s_inner(j,i)-1 429 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * shf(j,i) 430 ENDIF 431 432 IF ( use_top_fluxes ) THEN 433 k = nzt 434 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * tswst(j,i) 435 ENDIF 391 436 ENDDO 392 437 393 IF ( use_surface_fluxes ) THEN 394 k = nzb_diff_s_inner(j,i)-1 395 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * shf(j,i) 396 ENDIF 397 398 IF ( use_top_fluxes ) THEN 399 k = nzt 400 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * tswst(j,i) 401 ENDIF 402 ENDDO 438 ENDIF 403 439 404 440 ENDIF … … 762 798 IF ( .NOT. humidity ) THEN 763 799 764 IF ( use_pt_reference ) THEN 765 766 DO k = nzb_diff_s_inner(j,i), nzt_diff 767 tend(k,j,i) = tend(k,j,i) - kh(k,j,i) * g / pt_reference * & 800 IF ( use_reference ) THEN 801 802 IF ( ocean ) THEN 803 ! 804 !-- So far in the ocean no special treatment of density flux in the 805 !-- bottom and top surface layer 806 DO k = nzb_s_inner(j,i), nzt 807 tend(k,j,i) = tend(k,j,i) + kh(k,j,i) * g / prho_reference * & 808 ( rho(k+1,j,i) - rho(k-1,j,i) ) * dd2zu(k) 809 ENDDO 810 811 ELSE 812 813 DO k = nzb_diff_s_inner(j,i), nzt_diff 814 tend(k,j,i) = tend(k,j,i) - kh(k,j,i) * g / pt_reference * & 768 815 ( pt(k+1,j,i) - pt(k-1,j,i) ) * dd2zu(k) 769 ENDDO 770 771 IF ( use_surface_fluxes ) THEN 772 k = nzb_diff_s_inner(j,i)-1 773 tend(k,j,i) = tend(k,j,i) + g / pt_reference * shf(j,i) 816 ENDDO 817 818 IF ( use_surface_fluxes ) THEN 819 k = nzb_diff_s_inner(j,i)-1 820 tend(k,j,i) = tend(k,j,i) + g / pt_reference * shf(j,i) 821 ENDIF 822 823 IF ( use_top_fluxes ) THEN 824 k = nzt 825 tend(k,j,i) = tend(k,j,i) + g / pt_reference * tswst(j,i) 826 ENDIF 827 774 828 ENDIF 775 829 776 IF ( use_top_fluxes ) THEN 777 k = nzt 778 tend(k,j,i) = tend(k,j,i) + g / pt_reference * tswst(j,i) 830 ELSE 831 832 IF ( ocean ) THEN 833 ! 834 !-- So far in the ocean no special treatment of density flux in the 835 !-- bottom and top surface layer 836 DO k = nzb_s_inner(j,i), nzt 837 tend(k,j,i) = tend(k,j,i) + kh(k,j,i) * g / rho(k,j,i) * & 838 ( rho(k+1,j,i) - rho(k-1,j,i) ) * dd2zu(k) 839 ENDDO 840 841 ELSE 842 843 DO k = nzb_diff_s_inner(j,i), nzt_diff 844 tend(k,j,i) = tend(k,j,i) - kh(k,j,i) * g / pt(k,j,i) * & 845 ( pt(k+1,j,i) - pt(k-1,j,i) ) * dd2zu(k) 846 ENDDO 847 848 IF ( use_surface_fluxes ) THEN 849 k = nzb_diff_s_inner(j,i)-1 850 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * shf(j,i) 851 ENDIF 852 853 IF ( use_top_fluxes ) THEN 854 k = nzt 855 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * tswst(j,i) 856 ENDIF 857 779 858 ENDIF 780 859 781 ELSE 782 783 DO k = nzb_diff_s_inner(j,i), nzt_diff 784 tend(k,j,i) = tend(k,j,i) - kh(k,j,i) * g / pt(k,j,i) * & 785 ( pt(k+1,j,i) - pt(k-1,j,i) ) * dd2zu(k) 786 ENDDO 787 788 IF ( use_surface_fluxes ) THEN 789 k = nzb_diff_s_inner(j,i)-1 790 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * shf(j,i) 791 ENDIF 792 793 IF ( use_top_fluxes ) THEN 794 k = nzt 795 tend(k,j,i) = tend(k,j,i) + g / pt(k,j,i) * tswst(j,i) 796 ENDIF 797 798 ENDIF 860 ENDIF 799 861 800 862 ELSE -
palm/trunk/SOURCE/prognostic_equations.f90
r96 r97 5 5 ! ----------------- 6 6 ! prognostic equation for salinity, density is calculated from equation of 7 ! state for seawater, +eqn_state_seawater_mod 8 ! new argument zw in calls of diffusion_e, calc_mean_pt_profile renamed 7 ! state for seawater and is used for calculation of buoyancy, 8 ! +eqn_state_seawater_mod 9 ! diffusion_e is called with argument rho in case of ocean runs, 10 ! new argument zw in calls of diffusion_e, new argument pt_/prho_reference 11 ! in calls of buoyancy and diffusion_e, calc_mean_pt_profile renamed 9 12 ! calc_mean_profile 10 13 ! … … 148 151 ENDIF 149 152 CALL coriolis( i, j, 1 ) 150 IF ( sloping_surface ) CALL buoyancy( i, j, pt, 1, 4 )153 IF ( sloping_surface ) CALL buoyancy( i, j, pt, pt_reference, 1, 4 ) 151 154 CALL user_actions( i, j, 'u-tendency' ) 152 155 … … 283 286 ENDIF 284 287 CALL coriolis( i, j, 3 ) 285 IF ( .NOT. humidity) THEN286 CALL buoyancy( i, j, pt, 3,4 )288 IF ( ocean ) THEN 289 CALL buoyancy( i, j, rho, prho_reference, 3, 64 ) 287 290 ELSE 288 CALL buoyancy( i, j, vpt, 3, 44 ) 291 IF ( .NOT. humidity ) THEN 292 CALL buoyancy( i, j, pt, pt_reference, 3, 4 ) 293 ELSE 294 CALL buoyancy( i, j, vpt, pt_reference, 3, 44 ) 295 ENDIF 289 296 ENDIF 290 297 CALL user_actions( i, j, 'w-tendency' ) … … 654 661 .NOT. use_upstream_for_tke ) THEN 655 662 IF ( .NOT. humidity ) THEN 656 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 657 l_grid, pt, rif, tend, zu, zw ) 663 IF ( ocean ) THEN 664 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 665 l_grid, rho, prho_reference, rif, tend, & 666 zu, zw ) 667 ELSE 668 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 669 l_grid, pt, pt_reference, rif, tend, & 670 zu, zw ) 671 ENDIF 658 672 ELSE 659 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 660 l_grid, vpt, rif, tend, zu, zw ) 673 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 674 l_grid, vpt, pt_reference, rif, tend, zu, & 675 zw ) 661 676 ENDIF 662 677 ELSE … … 680 695 IF ( .NOT. humidity ) THEN 681 696 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, & 682 km_m, l_grid, pt_m, rif_m, tend, zu,&683 zw )697 km_m, l_grid, pt_m, pt_reference, & 698 rif_m, tend, zu, zw ) 684 699 ELSE 685 700 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, & 686 km_m, l_grid, vpt_m, rif_m, tend, zu,&687 zw )701 km_m, l_grid, vpt_m, pt_reference, & 702 rif_m, tend, zu, zw ) 688 703 ENDIF 689 704 ELSE 690 705 IF ( .NOT. humidity ) THEN 691 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 692 l_grid, pt, rif, tend, zu, zw ) 706 IF ( ocean ) THEN 707 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, & 708 km, l_grid, rho, prho_reference, & 709 rif, tend, zu, zw ) 710 ELSE 711 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, & 712 km, l_grid, pt, pt_reference, rif, & 713 tend, zu, zw ) 714 ENDIF 693 715 ELSE 694 716 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 695 l_grid, vpt, rif, tend, zu, zw ) 717 l_grid, vpt, pt_reference, rif, tend, & 718 zu, zw ) 696 719 ENDIF 697 720 ENDIF … … 797 820 ENDIF 798 821 CALL coriolis( i, j, 1 ) 799 IF ( sloping_surface ) CALL buoyancy( i, j, pt, 1, 4 ) 822 IF ( sloping_surface ) CALL buoyancy( i, j, pt, pt_reference, 1, & 823 4 ) 800 824 CALL user_actions( i, j, 'u-tendency' ) 801 825 … … 896 920 ENDIF 897 921 CALL coriolis( i, j, 3 ) 898 IF ( .NOT. humidity) THEN899 CALL buoyancy( i, j, pt, 3,4 )922 IF ( ocean ) THEN 923 CALL buoyancy( i, j, rho, prho_reference, 3, 64 ) 900 924 ELSE 901 CALL buoyancy( i, j, vpt, 3, 44 ) 925 IF ( .NOT. humidity ) THEN 926 CALL buoyancy( i, j, pt, pt_reference, 3, 4 ) 927 ELSE 928 CALL buoyancy( i, j, vpt, pt_reference, 3, 44 ) 929 ENDIF 902 930 ENDIF 903 931 CALL user_actions( i, j, 'w-tendency' ) … … 1114 1142 IF ( .NOT. humidity ) THEN 1115 1143 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, & 1116 km_m, l_grid, pt_m, rif_m, tend, zu,&1117 zw )1144 km_m, l_grid, pt_m, pt_reference, & 1145 rif_m, tend, zu, zw ) 1118 1146 ELSE 1119 1147 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, & 1120 km_m, l_grid, vpt_m, rif_m, tend, zu,&1121 zw )1148 km_m, l_grid, vpt_m, pt_reference, & 1149 rif_m, tend, zu, zw ) 1122 1150 ENDIF 1123 1151 ELSE 1124 1152 IF ( .NOT. humidity ) THEN 1125 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 1126 l_grid, pt, rif, tend, zu, zw ) 1153 IF ( ocean ) THEN 1154 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, & 1155 km, l_grid, rho, prho_reference, & 1156 rif, tend, zu, zw ) 1157 ELSE 1158 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, & 1159 km, l_grid, pt, pt_reference, rif, & 1160 tend, zu, zw ) 1161 ENDIF 1127 1162 ELSE 1128 1163 CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, & 1129 l_grid, vpt, rif, tend, zu, zw ) 1164 l_grid, vpt, pt_reference, rif, tend, & 1165 zu, zw ) 1130 1166 ENDIF 1131 1167 ENDIF … … 1224 1260 ENDIF 1225 1261 CALL coriolis( 1 ) 1226 IF ( sloping_surface ) CALL buoyancy( pt, 1, 4 )1262 IF ( sloping_surface ) CALL buoyancy( pt, pt_reference, 1, 4 ) 1227 1263 CALL user_actions( 'u-tendency' ) 1228 1264 … … 1367 1403 ENDIF 1368 1404 CALL coriolis( 3 ) 1369 IF ( .NOT. humidity) THEN1370 CALL buoyancy( pt, 3,4 )1405 IF ( ocean ) THEN 1406 CALL buoyancy( rho, prho_reference, 3, 64 ) 1371 1407 ELSE 1372 CALL buoyancy( vpt, 3, 44 ) 1408 IF ( .NOT. humidity ) THEN 1409 CALL buoyancy( pt, pt_reference, 3, 4 ) 1410 ELSE 1411 CALL buoyancy( vpt, pt_reference, 3, 44 ) 1412 ENDIF 1373 1413 ENDIF 1374 1414 CALL user_actions( 'w-tendency' ) … … 1753 1793 THEN 1754 1794 IF ( .NOT. humidity ) THEN 1755 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, pt, & 1756 rif, tend, zu, zw ) 1795 IF ( ocean ) THEN 1796 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, rho, & 1797 prho_reference, rif, tend, zu, zw ) 1798 ELSE 1799 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, pt, & 1800 pt_reference, rif, tend, zu, zw ) 1801 ENDIF 1757 1802 ELSE 1758 1803 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, vpt, & 1759 rif, tend, zu, zw )1804 pt_reference, rif, tend, zu, zw ) 1760 1805 ENDIF 1761 1806 ELSE … … 1777 1822 IF ( .NOT. humidity ) THEN 1778 1823 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e_m, km_m, l_grid, & 1779 pt_m, rif_m, tend, zu, zw )1824 pt_m, pt_reference, rif_m, tend, zu, zw ) 1780 1825 ELSE 1781 1826 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e_m, km_m, l_grid, & 1782 vpt_m, rif_m, tend, zu, zw )1827 vpt_m, pt_reference, rif_m, tend, zu, zw ) 1783 1828 ENDIF 1784 1829 ELSE 1785 1830 IF ( .NOT. humidity ) THEN 1786 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, pt, & 1787 rif, tend, zu, zw ) 1831 IF ( ocean ) THEN 1832 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, & 1833 rho, prho_reference, rif, tend, zu, zw ) 1834 ELSE 1835 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, & 1836 pt, pt_reference, rif, tend, zu, zw ) 1837 ENDIF 1788 1838 ELSE 1789 1839 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, vpt, & 1790 rif, tend, zu, zw )1840 pt_reference, rif, tend, zu, zw ) 1791 1841 ENDIF 1792 1842 ENDIF -
palm/trunk/SOURCE/run_control.f90
r90 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! 6 ! Timestep format changed 7 7 ! 8 8 ! Former revisions: … … 113 113 100 FORMAT (///'Run-control output:'/ & 114 114 &'------------------'// & 115 &'RUN ITER. HH:MM:SS.SS DT(E) UMAX VMAX WMAX U*', &116 &' W* THETA* Z_I ENERG. DISTENERG DIVOLD DIVNE', &117 &'W UMAX(KJI) VMAX(KJI) WMAX(KJI) ADVECX ADVECY MG', &118 &'CYC'/&119 120 121 122 &'--')123 101 FORMAT (I3,1X,I6,1X,A8,F3.2,1X,F 7.4,A1,A1,F8.4,A1,F8.4,A1,F8.4,2X,F5.3,2X, &115 &'RUN ITER. HH:MM:SS.SS DT(E) UMAX VMAX WMAX U', & 116 &'* W* THETA* Z_I ENERG. DISTENERG DIVOLD DIV', & 117 &'NEW UMAX(KJI) VMAX(KJI) WMAX(KJI) ADVECX ADVECY ', & 118 &'MGCYC'/ & 119 &'----------------------------------------------------------------', & 120 &'----------------------------------------------------------------', & 121 &'----------------------------------------------------------------', & 122 &'-----') 123 101 FORMAT (I3,1X,I6,1X,A8,F3.2,1X,F8.4,A1,A1,F8.4,A1,F8.4,A1,F8.4,2X,F5.3,2X, & 124 124 F4.2, & 125 2X,F6.3,2X,F 5.0,1X,4(E10.3,1X),3(3(I4),1X),F8.3,1X,F8.3,5X,I3)125 2X,F6.3,2X,F6.0,1X,4(E10.3,1X),3(3(I4),1X),F8.3,1X,F8.3,5X,I3) 126 126 127 127 END SUBROUTINE run_control -
palm/trunk/SOURCE/time_integration.f90
r95 r97 4 4 ! Actual revisions: 5 5 ! ----------------- 6 ! Ghostpoint exchange for salinity and density 6 ! diffusivities is called with argument rho in case of ocean runs, 7 ! new argument pt_/prho_reference in calls of diffusivities, 8 ! ghostpoint exchange for salinity and density 7 9 ! 8 10 ! Former revisions: … … 249 251 CALL cpu_log( log_point(17), 'diffusivities', 'start' ) 250 252 IF ( .NOT. humidity ) THEN 251 CALL diffusivities( pt ) 253 IF ( ocean ) THEN 254 CALL diffusivities( rho, prho_reference ) 255 ELSE 256 CALL diffusivities( pt, pt_reference ) 257 ENDIF 252 258 ELSE 253 CALL diffusivities( vpt )259 CALL diffusivities( vpt, pt_reference ) 254 260 ENDIF 255 261 CALL cpu_log( log_point(17), 'diffusivities', 'stop' )
Note: See TracChangeset
for help on using the changeset viewer.
