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