-
Notifications
You must be signed in to change notification settings - Fork 0
/
simple_ocean.f90
147 lines (112 loc) · 3.86 KB
/
simple_ocean.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
module simple_ocean
!------------------------------------------------------------
! Purpose:
!
! A collection of routines used to specify fixed
! or compute interactive SSTs, like slab-ocean model, etc.
!
! Author: Marat Khairoutdinov
! Based on dynamic ocean impelemntation from the UW version of SAM.
!------------------------------------------------------------
use grid
implicit none
public set_sst ! set SST
public sst_evolve ! evolve SST according to a model set by the ocean_type
CONTAINS
SUBROUTINE set_sst()
use vars, only: sstxy,t00
use params, only: tabs_s, delta_sst, ocean_type
! parameters of the sinusoidal SST destribution
! along the X for Walker-type simulatons( ocean-type = 1):
real(8) tmpx(nx), pii, lx, yy
integer i,j, it,jt
select case (ocean_type)
case(0) ! fixed constant SST
sstxy = tabs_s - t00
case(1) ! Sinusoidal distribution along the x-direction:
lx = float(nx_gl)*dx
do i = 1,nx
tmpx(i) = float(mod(rank,nsubdomains_x)*nx+i-1)*dx
end do
pii = atan2(0.d0,-1.d0)
do j=1,ny
do i=1,nx
sstxy(i,j) = tabs_s-delta_sst*cos(2.*pii*tmpx(i)/lx) - t00
end do
end do
case(2) ! Sinusoidal distribution along the y-direction:
call task_rank_to_index(rank,it,jt)
pii = atan2(0.d0,-1.d0)
lx = float(ny_gl)*dy
do j=1,ny
yy = dy*(j+jt-(ny_gl+YES3D-1)/2-1)
do i=1,nx
sstxy(i,j) = tabs_s+delta_sst*(2.*cos(pii*yy/lx)-1.) - t00
end do
end do
case default
if(masterproc) then
print*, 'unknown ocean type in set_sst. Exitting...'
call task_abort
end if
end select
end subroutine set_sst
SUBROUTINE sst_evolve
use vars, only: sstxy, t00, fluxbt, fluxbq, rhow,qocean_xy
use params, only: cp, lcond, tabs_s, ocean_type, dossthomo, &
depth_slab_ocean, Szero, deltaS, timesimpleocean
use rad, only: swnsxy, lwnsxy
real, parameter :: rhor = 1000. ! density of water (kg/m3)
real, parameter :: cw = 4187. ! Liquid Water heat capacity = 4187 J/kg/K
real factor_cp, factor_lc, qoceanxy
real tmpx(nx), lx
real(8) sss,ssss, tmp1(1), tmp2(1)
integer i,j
if(time.lt.timesimpleocean) return
lx = float(nx_gl)*dx
do i = 1,nx
tmpx(i) = float(mod(rank,nsubdomains_x)*nx+i-1)*dx
end do
! Define weight factors for the mixed layer heating due to
! the model's sensible and latent heat flux.
factor_cp = rhow(1)*cp
factor_lc = rhow(1)*lcond
! Use forward Euler to integrate the differential equation
! for the ocean mixed layer temperature: dT/dt = S - E.
! The source: CPT?GCSS WG4 idealized Walker-circulation
! RCE Intercomparison proposed by C. Bretherton.
do j=1,ny
do i=1,nx
qoceanxy = Szero + deltaS*abs(2.*tmpx(i)/lx - 1)
qocean_xy(i,j) = qocean_xy(i,j) + qoceanxy
sstxy(i,j) = sstxy(i,j) &
+ dtn*(swnsxy(i,j) & ! SW Radiative Heating
- lwnsxy(i,j) & ! LW Radiative Heating
- factor_cp*fluxbt(i,j) & ! Sensible Heat Flux
- factor_lc*fluxbq(i,j) & ! Latent Heat Flux
+ qoceanxy) & ! Ocean Heating
/(rhor*cw*depth_slab_ocean) ! Convert W/m^2 Heating to K/s
end do
end do
if(dossthomo) then
sss = 0.
do j=1,ny
do i=1,nx
sss = sss + sstxy(i,j)
end do
end do
sss = sss / dble(nx*ny)
if(dompi) then
tmp1(1) = sss
call task_sum_real8(tmp1,tmp2,1)
sss = tmp2(1) /float(nsubdomains)
end if ! dompi
if(ocean_type.eq.2) then
tabs_s = sss + t00
call set_sst()
else
sstxy(:,:) = sss
end if
end if
end subroutine sst_evolve
end module simple_ocean