[95] | 1 | SUBROUTINE init_ocean |
---|
| 2 | |
---|
| 3 | !------------------------------------------------------------------------------! |
---|
| 4 | ! Actual revisions: |
---|
| 5 | ! ----------------- |
---|
[388] | 6 | ! Bugfix: Initial profiles of hydrostatic pressure and density are calculated |
---|
| 7 | ! iteratively. First calculation of hyp(0) changed. |
---|
[95] | 8 | ! |
---|
| 9 | ! Former revisions: |
---|
| 10 | ! ------------------ |
---|
[96] | 11 | ! $Id: init_ocean.f90 388 2009-09-23 09:40:33Z raasch $ |
---|
[95] | 12 | ! |
---|
[139] | 13 | ! 124 2007-10-19 15:47:46Z raasch |
---|
| 14 | ! Bugfix: Initial density rho is calculated |
---|
| 15 | ! |
---|
[98] | 16 | ! 97 2007-06-21 08:23:15Z raasch |
---|
| 17 | ! Initial revision |
---|
[95] | 18 | ! |
---|
| 19 | ! Description: |
---|
| 20 | ! ------------ |
---|
| 21 | ! Initialization of quantities needed for the ocean version |
---|
| 22 | !------------------------------------------------------------------------------! |
---|
| 23 | |
---|
| 24 | USE arrays_3d |
---|
| 25 | USE control_parameters |
---|
| 26 | USE eqn_state_seawater_mod |
---|
[336] | 27 | USE pegrid |
---|
[95] | 28 | USE grid_variables |
---|
| 29 | USE indices |
---|
| 30 | |
---|
| 31 | IMPLICIT NONE |
---|
| 32 | |
---|
[336] | 33 | INTEGER :: k, n |
---|
[95] | 34 | |
---|
[388] | 35 | REAL :: sa_l, pt_l |
---|
[95] | 36 | |
---|
[336] | 37 | REAL, DIMENSION(nzb:nzt+1) :: rho_init |
---|
| 38 | |
---|
[95] | 39 | ALLOCATE( hyp(nzb:nzt+1) ) |
---|
| 40 | |
---|
| 41 | ! |
---|
| 42 | !-- Set water density near the ocean surface |
---|
| 43 | rho_surface = 1027.62 |
---|
| 44 | |
---|
| 45 | ! |
---|
| 46 | !-- Calculate initial vertical profile of hydrostatic pressure (in Pa) |
---|
[96] | 47 | !-- and the reference density (used later in buoyancy term) |
---|
[388] | 48 | !-- First step: Calculate pressure using reference density |
---|
[95] | 49 | hyp(nzt+1) = surface_pressure * 100.0 |
---|
| 50 | |
---|
[97] | 51 | hyp(nzt) = hyp(nzt+1) + rho_surface * g * 0.5 * dzu(nzt+1) |
---|
[336] | 52 | rho_init(nzt) = rho_surface |
---|
[95] | 53 | |
---|
[366] | 54 | DO k = nzt-1, 1, -1 |
---|
[336] | 55 | hyp(k) = hyp(k+1) + rho_surface * g * dzu(k) |
---|
| 56 | ENDDO |
---|
[366] | 57 | hyp(0) = hyp(1) + rho_surface * g * dzu(1) |
---|
[95] | 58 | |
---|
[388] | 59 | ! |
---|
| 60 | !-- Second step: Iteratively calculate in situ density (based on presssure) |
---|
| 61 | !-- and pressure (based on in situ density) |
---|
[336] | 62 | DO n = 1, 5 |
---|
[95] | 63 | |
---|
[336] | 64 | rho_reference = rho_surface * 0.5 * dzu(nzt+1) |
---|
[95] | 65 | |
---|
[336] | 66 | DO k = nzt-1, 0, -1 |
---|
| 67 | |
---|
| 68 | sa_l = 0.5 * ( sa_init(k) + sa_init(k+1) ) |
---|
| 69 | pt_l = 0.5 * ( pt_init(k) + pt_init(k+1) ) |
---|
| 70 | |
---|
| 71 | rho_init(k) = eqn_state_seawater_func( hyp(k), pt_l, sa_l ) |
---|
| 72 | |
---|
| 73 | rho_reference = rho_reference + rho_init(k) * dzu(k+1) |
---|
| 74 | |
---|
| 75 | ENDDO |
---|
| 76 | |
---|
| 77 | rho_reference = rho_reference / ( zw(nzt) - zu(nzb) ) |
---|
| 78 | |
---|
| 79 | DO k = nzt-1, 0, -1 |
---|
| 80 | hyp(k) = hyp(k+1) + g * 0.5 * ( rho_init(k) + rho_init(k+1 ) ) * & |
---|
| 81 | dzu(k+1) |
---|
| 82 | ENDDO |
---|
| 83 | |
---|
[95] | 84 | ENDDO |
---|
| 85 | |
---|
[97] | 86 | ! |
---|
| 87 | !-- Calculate the reference potential density |
---|
| 88 | prho_reference = 0.0 |
---|
| 89 | DO k = 0, nzt |
---|
[96] | 90 | |
---|
[97] | 91 | sa_l = 0.5 * ( sa_init(k) + sa_init(k+1) ) |
---|
| 92 | pt_l = 0.5 * ( pt_init(k) + pt_init(k+1) ) |
---|
| 93 | |
---|
| 94 | prho_reference = prho_reference + dzu(k+1) * & |
---|
[336] | 95 | eqn_state_seawater_func( 0.0, pt_l, sa_l ) |
---|
[97] | 96 | |
---|
| 97 | ENDDO |
---|
| 98 | |
---|
| 99 | prho_reference = prho_reference / ( zu(nzt) - zu(nzb) ) |
---|
| 100 | |
---|
[124] | 101 | ! |
---|
[388] | 102 | !-- Calculate the 3d array of initial in situ and potential density, |
---|
| 103 | !-- based on the initial temperature and salinity profile |
---|
[124] | 104 | CALL eqn_state_seawater |
---|
[97] | 105 | |
---|
[124] | 106 | |
---|
[95] | 107 | END SUBROUTINE init_ocean |
---|