Transcript Document

1
-SELFE Users Trainng Course(1)
SELFE physical formulation
Joseph Zhang, CMOP, OHSU
Governing equations: Reynolds-averged Navier-Stokes




w
u 
 0, (u  (u, v))
z
Continuity equation


    udz  0
h
t
Momentum equations
Du
  u 
1
g 
ˆ
 f  g + 
;
f


f
k

u


g




p

 d   ( u)

Dt
z  z 
0 A 0 z
u

Vertical b.c.:

 τ w at z =
z
u

 CD | ub | ub , at z  h
z
z
Equation of state
   ( p, S , T )
Transport of salt and temperature
Dc   c 
     Q   ( hc),
Dt z  z 

2
MSL
x
East
c  (S ,T )
Turbulence closure: Umlauf and Burchard 2003
shear
stratification
Dk    k 
2
2
  k
  K mv M  K hv N  
Dt z  z 
D     
2
2
 
  c 1K mv M  c 3 K hv N  c 2 Fwall 
Dt z 
z  k


  c0  k m  n ,
p
Continuity equation (1)
3


DyDz   ( x  Dx)u ( x  Dx)   ( x)u ( x)  

D( V )  Dt  


DxDz   ( y  Dy )v( y  Dy )   ( y )v( y )   DxDy   ( z  Dz ) w( z  Dz )   ( z ) w( z ) 

D( V )  DDxDyDz
Eulerian form:
Lagrangian form:
w
  (  u )  (  v)  (  w)



t
x
y
z
 3  (  u)

(Dt , DV  0)

D
 3  (  u) 
 3  u  0
t
Dt
v
Dz
u
Dx
z
y
Law of mass conservation + incompressibility of water

Vertical integration
Boussinesq:

  0  3  u  0
(Volume conservation)
x


h

  udz  w  h  0
Continuity equation (2)
Vertical boundary conditions: kinematic
z   ( x, y , t )
surface
w
dz dx  dy  



 u x  v y  t
dt dt x dt y t
z   h ( x, y )
bottom

dz
w   uhx  vhy

  udtdz  w  h  0

Integrated continuity equation
h


    udz  0
h
t


h
  udz  t  u x  v y  uhx  vhy  0


h
h
   udz     udz  u (  h)
4
z
Momentum equation (1)
Newton’s 2nd law (Lagrangian)
5
F  ma
Du
1
  P + 3  ( 3u )  fk  u
Dt

Dw
1 P

 g +  (w)
Dt
 z
f
Hydrostatic assumption

1 P

 g  0  P  g   d   PA
z
 z
Separation of scales

Du
1
  u 
  (   d  PA ) +  ( u)  
  fk  u
z
Dt

z  z 
Du
  u 
1
g
  g +
p A 

  fk  u   gˆ 
Dt
z  z 
0
0
Boussinesq assumption


z
z
  d    d  


z
 d   ( u)
6
Coriolis

Coriolis: “virtual force”




Newton’s 2nd law is only valid in an inertial frame
For large scale GFD flows, the earth is not an inertial frame
Small as it is, it is reponsible for a number of unique phenomena in
GFD
Accelerations in a rotating frame (W const.):
U  u  W r
~
~
~
~
dR
or :
d

   W r
dt  dt ~  ~
~
j
Coriolis
d

A    WU  a  2 W u  W W r
~
~
~
~
~
~
 dt ~  ~ ~
 
J
u
Wxr
Centrifugal
acceleration
i
I
W
7
Coriolis

Plane coordinates: neglect curvature of the earth

x: east; y: north; z: vertical upward

Projection of earth rotation in the local frame:
W  W cosf j  W sin f k
y
e
~

~
~
Components of Coriolis acceleration:
x:
x
f* w  fv
)f
y : fu
z :  f*u
f*  2W cosf , f  2W sin f

Inertial oscillation: no external forces, advection and vertical vel.

e.g., sudden shut-down of wind in ocean.
W
du

 fv  0
 u  V sin( ft   )
dt

dv
v  V cos( ft   )
 fu  0 

dt
u  0
~
z
Momentum equation: vertical boundary condition (b.c.)
Surface


Logarithmic law:
u
 τ w at z =
z
u
ln  ( z  h) / z0 
ln( b / z0 )
Bottom

8
u
 CD | ub | ub , at z  h
z
u b , ( z0  h  z   b  h )
z0 : bottom roughness


Reynolds stress:
Turbulence closure:

u


ub
z ( z  h)ln( b / z0 )
ub
b
  2sm K 1/ 2l ,
sm  g 2 ,
1 2/3
B1 CD | ub |2
2
l   0 ( z  h)
K

Reynolds stress (const.)


Drag coefficient:
0
u

CD1/ 2 | ub | ub , ( z0  h  z   b  h)
z ln( b / z0 )
 1  
CD   ln b 
  0 z0 
2
z=h
Momentum equation (3)
9
meters
Implications
• Use of z0 seems most natural option, but over-estimates CD in shallow
area
• Strictly speaking CD should vary with bottom layer thickness (i.e. vertical
grid)
• Different from 2D, 3D results of elevation depend on vertical grid (with
constant CD)
Transport equation (1)
10
Mass conservation: advection-diffusion equation
Dc
  c 


  Q   ( hc ),
Dt
z  z 
c  (S ,T )
Fick’s law for diffusive fluxes
   ( p, S , T )
qx Dy Dz
c
Dx Dy Dz
t
qx


 qx  x Dx  Dy Dz


>
Advection
B.c.
c
  cˆ, z  
z
c
  0, z  h
n
Diffusion
Transport equation (2)
11
Precipitation and evaporation model

S
S ( E  P)

, z 
z
0
E: evaporation rate (kg/m2/s) (either measured or calculated)
P: Precipitation rate (kg/m2/s) (usually measured)
Heat exchange model
T
H

, z 
z 0C p
1 SW
Q
 0C p z

H  ( IR   IR )  S  E
IR’s are down/upwelling infrared (LW) radiation at surface
S is the turbulent flux of sensible heat (upwelling)
E is the turbulent flux of latent heat (upwelling)
SW is net downward solar radiation at surface
The solar radiation is penetrative (body force) with
attenuation (which depends upon turbidity) acts as a heat
source within the water.
Air-water exchange
1. free surface height: variations in atmospheric pressure over the
domain have a direct impact upon free surface height; set up due
to wind stresses
2. momentum: near-surface winds apply wind-stress on surface
(influences advection, location of density fronts)
3. heat: various components of heat fuxes (dependent upon many
variables) determine surface heat budget
• shortwave radiation (solar) - penetrative
• longwave radiation (infrared)
• sensible heat flux (direct transfer of heat)
• latent heat flux (heating/cooling associated with
condensation/evaporation)
4. water: evaporation, condensation, precipitation act as
sources/sinks of fresh water
12
Solar radiation
13
Downwelling SW at the surface is forecast in NWP models - a function of time of year,
time of day, weather conditions, latitude, etc
Upwelling SW is a simple function of downwelling SW
SW   SW 
 is the albedo - typically depends on solar zenith angle and sea state
Attenuation of SW radiation in the water column is a function of turbidity and depth d
(Jerlov, 1968, 1976; Paulson and Clayson, 1977):
SW ( z )  (1   ) SW  Re D / d1  (1  R)e  D / d2 
R, d1 and d2 depend on water type; D is the distance from F.S.
R
d1 (m)
d2 (m)
Jerlov I
0.58
0.35
23
Jerlov IA
0.62
0.60
20
Jerlov IB
0.67
1.00
17
Jerlov II
0.77
1.50
14
Jerlov III
0.78
1.40
7.9
Paulson and
Simpson
0.62
1.50
20
Estuary
0.80
0.90
2.1
III
Depth (m)
Type
II
I
Infrared radiation
• Downwelling IR at the surface is forecast in NWP models - a function of
air temperature, cloud cover, humidity, etc
• Upwelling IR is can be approximated as either a broadband
measurement solely within the IR wavelengths (i.e., 4-50 μm), or more
commonly the blackbody radiative flux
4
IR  s Tsfc
 is the emissivity, ~1
s is the Stefan-Boltzmann constant
Tsfc is the surface temperature
14
Turbulent Fluxes of Sensible and Latent Heat
In general, turbulent fluxes are a function of:
• Tsfc,, Tair
• near-surface wind speed
• surface atmospheric pressure
• near-surface humidity
Scales of motion responsible for these heat fluxes are much smaller than can
be resolved by any operational model - they must be parameterized (i.e., bulk
aerodynamic formulation)
S    a C pa w 'T '    a C pa u*T*
E    a Le w ' q '    a Leu*q*
where:
u* is the friction scaling velocity
T* is the temperature scaling parameter
q* is the specific humidity scaling parameter
a is the surface air density
Cpa is the specific heat of air
Le is the latent heat of vaporization
The scaling parameters are defined using Monin-Obukhov similarity theory,
and must be solved for iteratively (i.e., Zeng et al., 1998).
15
Wind Shear Stress: Turbulent Flux of Momentum
Calculation of shear stress follows naturally from calculation of
turbulent heat fluxes
Total shear stress of atmosphere upon surface:
w 
a 2
u
0 *
Alternatively, Pond and Picard’s formulation can be used as a simpler
option
τw 
a
Cds | u w | u w
0
'
0.61  0.063u w
Cds 
1000
'
uw
 max(6, min(50, u w ))
16
Inputs to heat/momentum exchange model
17
• Forecasts of uw vw @ 10m, PA @ MSL,Tair and qair @2m
• used by the bulk aerodynamic model to calculate the scaling parameters
• Forecasts of downward IR and SW
Atmospheric properties
supplying
agency
time period
spatial resolution
MRF
NCEP
04/01/200102/29/2004
1° x 1°
GFS
NCEP
07/03/2003-present
OSU-ARPS
OSU
ETA/NAM
NCAR/NCEP
Reanalysis (NARR)
data source
temporal
resolution
area of coverage
data type
12-hour
snapshots
129W-120W, 35N-51N
forecast
1° x 1°
3-hour snapshots
180W-70W, 90S-90N
forecast
05/04/200102/25/2004
12 km
1-hour snapshots
128W-119W, 41N-47N
(approx)
forecast
NCEP
07/03/2003-present
12 km
3-hour snapshots
Eta Grid 218, west of
100W
forecast
NCAR
01/01/197912/31/2006
12km
6-hour snapshots
North America
reanalysis
Heat fluxes
supplying
agency
time period
spatial resolution
AVN (lo-res)
NCEP
04/01/200110/28/2002
0.7° x 0.7° (approx)
AVN (hi-res)
NCEP
10/29/200202/29/2004
GFS
NCEP
ETA/NAM
NCAR/NCEP
Reanalysis (NARR)
data source
temporal
resolution
area of coverage
data type
3-hour averages
129W-120W, 35N-51N
forecase
0.5° x 0.5° (approx)
3-hour averages
129W-120W, 35N-51N
forecast
07/03/2003-present
0.5° x 0.5° (approx)
3-hour averages
180W-70W, 90S-90N
forecast
NCEP
07/03/2003-present
12 km
3-hour snapshots
NCAR
01/01/197912/31/2006
12km
6-hour averages
Eta Grid 218, west of
100W
North America
forecast
reanalysi
s
Turbulence closure: Pacanowski and Philander (1982)
• Kelvin-Helmholtz instability: occurs at interface of two fluids or
fluids of different density (stratification)
• When stratification is present, the onset of instability is govern
by the Richardson number

Ri 
g 
 0 z
 u    v 




 z 
 z 
2
2

buoyancy
shear
• Theory predicts that the fluid is unstable when Ri<0.25
 

 max
(1  5 Ri ) 2

1  5 Ri
  min
  min
18
Turbulence closure: Mellor-Yamada-Galperin (1988)
19
• Theorize that the turbulence mixing is related to the growth and dissipation of
turbulence kinetic energy (k) and Monin-Obukhov mixing length (l)
• k and l are derived from the 2nd moment turbulence theory
Dk
 
k 

 k
  M
Dt
z 
z 
2
lE1
D( kl )
 
 ( kl ) 



 k

Dt
z 
z 
2

Dissipation rate
M
Wall proximity function
2
 N2
(2k ) 3 / 2

B1l
2
 u 
 v 

 

 z 
 z 
  1  E2
2
N2 


l
 d 

0
b 
  sm (2k )1/ 2 l
Viscosity and diffusivities
Galperin clipping
 M
2
Shear & buoyancy freq
Stability functions
 N2 
2
   2l
g 
0 z


l
 E3 

 0ds 
2
  sh (2k )1/ 2 l
 k  0.2(2k )1/ 2 l
sm 
g 2  g3Gh
g6
sh 
(1  g 4Gh )(1  g 5Gh )
1  g 4Gh
0.28  Gh
l2

2k
g 
 0.0233
 0 z
Turbulence closure: Mellor-Yamada-Galperin (1988)
• Vertical b.c.
• TKE related to stresses
• Mixing length related to the distance from the “walls” (note the
singularity)
B12 / 3 u
k 

2
z
l / 0  d
• Initial condition
k  5 106 m 2 / s
 =108
20
Turbulence closure: Umlauf and Burchard (2003)
21
• Use a generic length-scale variable to represent various closure schemes
• Mellor-Yamada-Galperin (with modification to wall-proximity function
 k-kl)
• k- (Rodi)
• k-w (Wilcox)
• KPP (Large et al.; Durski et al)
Dk
   k 
2
2

 k
  M   N  
Dt
z 
z 
D

  
2
2
 
   c 1 M  c 3 N  c 2 Fwall  
Dt z 
z  k
2
Wall function
Fwall
 l 
 l 
 1  E2 
 E4 


  0 db 
 0ds 
1/ 2
Diffusivity   c k
  c' k1/ 2
 k 

s k
   c  k m
0
p
   c  k 3/ 2
0 3
n
1
2
 

s
Stability: Kantha and Clayson (smoothes Galperin’s stability function as it approaches max)
c  2sm , c'  2sh
sh 
0.4939
1  30.19Gh
sm 
0.392  17.07 shGh
1  6.127Gh
0.28  Gh 
Gh _ c  0.02
Gh _ u  (Gh _ u  Gh _ c ) 2
Gh _ u  Gh 0  2Gh _ c
 0.023
Turbulence closure: Umlauf and Burchard (2003)
22
Constants
Wall proximity function
c 3
 1, N 2  0
 
2
c 3 , N  0
c 3
depends on choice of scheme and
stability function
Geophysical fluid dynamics in ocean

Coriolis plays an important role in GFD
u
u
u
u
1 p   u 
u
 v  w  fv  
 

t
x
y
z
 0 x z  z 
U
T
U2
L
U2
L
UW
H
WU
P
ρ0 L
1
WT
U
WL
U
WL
WL
UH
1
P
ρ0 WLU
RoT  1 Ro  1


U
H2

WH 2
Ek  1
Ekman and Rossby numbers
Although Ek<<1, the friction term is important in boundary layer, and is the
highest-order term in the eq.
Reynolds number:
R
Re 
 o

Ek
UL

2
L
   1
H
Eddy viscosity >> molecular viscosity
23
Geostrophy


24
Rapidly rotating (Coriolis dominant), homogeneous inviscid fluids:
N-S eq:
W
fv 
RoT , Ro , Ek  1,   0
1 p
,
 0 x
fu  
0
1 p
,
 0 y
1 p
,
 0 z
  u  0,


(a) 

(b) 



(c ) 

(d ) 

low
>
>
u
~
>
p
Taylor-Proudman theorem (vertical rigidity):
high
high
<
<
low
u v

0
z z

Flow is along isobars, i.e., streamline = isobar


No work done by pressure; no energy required to sustain the flow (inertial)
Direction along the isobars
25
Surface Ekman layer


EL exists whenever viscosity/shear is important
Steady, homogeneous, uniform interior:
d u
 f (v  v )  
,
z
dz
d v
f (u  u )  
,
wind
dz
2
2
2
2
du
dv
  x ,  0
  y , at z  0
dz
dz
u  u, v  v, z  
 0

EL
Solution:
(u , v )
u u 
2 z/d  x
z 
 z  
e  cos     y sin   
 0 fd
d 4
 d 4 

vv
2 z/d  x  z   y
 z  
e  sin      cos  
 0 fd
d 4
 d 4 

d
2
f
interior
d
26
Ekman transport

Properties:




Velocities can be very large for nearly inviscid fluid
Wind-generated vel. makes a 45o angle to the wind direction
Ekman spiral
z
Net horizontal transport:
0
1

0 f
U   (u  u )dz 
 y,
0
1

0 f
V   (v  v )dz  
 x,
(U ,V )  ( x , y )  0


Ekman transport is perpendicular to the direction of wind stress (to the
right in the Northern Hemisphere)
Bottom Ekman layer
) 45o
Surface current
Ekman transport
Barotropic waves

27
Linear wave theory

crest
2D monochromatic wave:
  A cos(lx  m y  wt  f )  Re Aei (lx my wt f ) 

A : amplitude
  (l,m) : wave number vector
~
w : angular frequency
f : referencephase
phaseline : lx  m y  wt  f  const.
2
wave length:  
| |
trough
y
~

Fixed (x,y), varying t:
wT=2  T=2/w (period)
L2
(x,y) fixed
L1
t fixed


T
=A
=0 =A=0 =A
x
28
Phase speed


Vary (x,y,t), but follow a phase line L: lx+my-wt=C.
The speed of a phase
w
c
| |
in thedirectionof 
~
~

In many dynamic systems, w and  are related in the dispersion relation:
w  w ( )
~


e.g., surface gravity waves: w  g tanhh
Dispersion: a group of waves with different wavelength travel with
different speed
Non-dispersive waves: tides
Group velocity


29
Wave energy ~ A2
Two 1D waves of equal amplitude and almost equal wave numbers:
C
  A cos( 1 x  w1t )  A cos( 2 x  w 2t ),
w i  w i ( i ) (i  1,2),
 
  1 2 , D   1   2 , | D ||  |,
Cg
2
w  w2
w 1
, Dw  w1  w 2 , | Dw || w |,
2
Dw 
 D
  2 A cos
x
t  cos(x  wt )
2 
 2
’
average wave
Modulation/envelop
of average wave
Dw : beatingfreq.

The speed at which the envelop/energy travels is the group velocity
cg 
Dw dw

D d
Linearization


Linearization of advective terms: Ro << 1
Shallow water model:
u

 fv   g
(a)
t
x
v

 fu   g
(b)
t
y
h  (hu)  (hv)


0
t
x
y

30
z
h(x,y,t)
b(x,y)
(wave continuity equation)
O
For flat bottom, h=+H (H: mean depth). Linearizing for small-amplitude
waves ~A<<H leads to:
 u v 

 H     0 (c)
t
 x y 


(a-c) govern the dynamics of linear barotropic waves on a flat bottom
Neglecting Coriolis, the 1D wave equation is recovered
tt  gH xx
phase speed
Cp 
gH
Kelvin waves



31
Semi-infinite domain of flat bottom, free slip at coast
(a-c) become (3 eqs. for 2 unknowns; u=0):

x
v

 g
t
y

v
H
0
u0
t
y
c
v  gH F ( y  ct )e  x / R
R  , c  gH
f
   HF ( y  ct )e  x / R
fv  g

y
Solution:
u=0
x
where F is an arbitrary function
Properties:

Kelvin wave is trapped near the coast; Rossby radius of deformation R is a measurement
of trapping;

The wave travels alongshore w/o deformation, at the speed of surface gravity waves c,
with the coast on its right in the N-H;
W

The direction of the alongshore current (v) is arbitrary:

Upwelling wave >0  v<0, same direction as phase speed;

Downwelling wave <0  v>0, opposite direction to phase speed

As f0, Kelvin wave degenerates to plane gravity wave;

Generated by ocean tides and wind near coastal areas
Natural frequency of a stratified system

A simple 1D model for continuously stratified fluid



Horizontally homogenous: (z)
Static equilibrium initially
Wave motion due to buoyancy
d 2h
 ( z )V 2  gV  ( z  h)   ( z )
dz
2
d h g d

h0
dz2  0 dz

(z+h)
z
(z)
(z-h)
The freq. of oscillation is
N2  


g d
 0 dz
Brunt-Vaisala frequency
If N2 > 0, stable stratification  internal wave motion possible
If N2 < 0, unstable stratification  small disturbances leads to
overturning; no internal waves possible
32
Internal waves


Surface & internal gravity waves: same origin (interface); interchange of kinetic
and potential energy
Quintessential internal waves:





Over an infinite domain (even in the vertical!) with vertical stratification
No ambient rotation (gravity waves)
Neglect viscosity; small-amplitude motion
Hydrodynamic effects are important
u
1 p '

t
 0 x
Perturbation expansion:
   0   ( z )   ' ( x, t ),
|  ' ||  |  0
p
| p' || p |, 
~
p ( z )  p' ( x, t ),
~
1 p
g
 0   z
(u,v,w), generated by perturbation, are small

33
Wave solution ~ exp[i(lx+my+nzwt)] 
2
2
l

m
w2  N 2 2
l  m2  n 2
v
1 p '

t
 0 y
w
1 p '
'

g
t
 0 y
0
u v w
 
0
x y z
 '
d
 ' N 2  0
w
 0, or :

w0
t
dz
t
g
Transport eq.

Properties of internal waves


)q
For a given N, wmaxN
Dispersion relation only depends on the angle between the wave number and horizontal
w   N cos q
plane:



34
For a given frequency, all waves propagate at a a fixed angle to the horizontal
As w0, q90o
Wave motion:
l
wN
l  n2
u  U sin(lx  nz  wt )
2
, C
)q
l
w u
n
' 
0 N l 2  n 2
g
(u , w)    0
n
U cos(lx  nz  wt )
~






Cg
Transverse wave
Group vel. inside the phase plane
Purely vertical motion when n=0; Cg=0  no energy radiation
Purely horizontal motion when l=0; w=0 steady state; blocking
For continuously stratified fluid on a finite depth, the phase speed is bounded by a maximum but
there are infinitely many modes of internal waves
z
x
Inertial-internal waves:
f w  N