# Parameter Estimation for a 2D Tidal Model with POD 4D VAR Data Assimilation

Authors: Shouguo Qian; Xianqing Lv; Yanhua Cao; Fenjing Shao
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

Combining the proper orthogonal decomposition (POD) reduced order method and 4D VAR (four-dimensional Variational) data assimilation method with a two-dimensional (2D) tidal model, a model is constructed to simulate theM
2 tide in the Bohai, Yellow, and East China Seas (BYECS). This model consists of two submodels: the POD reduced order forward model is used to simulate the tides, while its adjoint model is used to optimize the control variables. Numerical experiment is carried out to assimilate the harmonic constants, which are derived from TOPEX/Poseidon (T
/
P) altimeter data, into the 2D tidal model through optimizing the initial values and the temporally and spatially varying open boundary conditions (OBCs). The absolute mean difference between the model results and observations is 3.2 cm and 2.9∘ for amplitude and phase-lag, respectively, better than the results of Lu and Zhang (2006), suggesting that the construction of the POD reduced order model and the inversion of control variables are successful.


---

## Body

## 1. Introduction

In the recent years, with the development of large ocean observing programs and remote sensing techniques, more and more oceanic data are becoming available. Data assimilation provides a promising prospect for using these observations to improve the simulation precision of tidal models. The 4D VAR data assimilation is one of the most widely used data assimilation methods and has been widely applied to sensitivity studies and parameter estimation (see, e.g., [1–9]). By assimilating the TOPEX/Poseidon (T
/
P) altimeter data into a 2D tidal model, Lu and Zhang [10] estimated the spatially varying bottom friction coefficient (BFC). In their work, the conception of independent BFC was introduced based on the spatial distribution of water depth, whereby the simulation precision reached a high level with absolute mean simulation errors of 5.0 cm and 5.4∘ in amplitude and phase-lag, respectively. By discretizing the primitive equations of motion and continuity in spherical coordinates, Zhang and Lu [11, 12] constructed a 3D numerical adjoint barotropic tidal model to estimate three kinds of model parameters including the OBCs, BFC, and the vertical eddy coefficients. Chen et al. [13] constructed a 3D numerical internal tidal model to invert the OBCs. Gao et al. [14] studied the parameterization of an internal tide dissipation term in a 2D shallow water model and carried out some related numerical experiments in the South China Sea where the BFC and internal tide dissipation coefficient were optimized.One problem of the 4D VAR data assimilation is characterized by the great CPU time and memory load required for the computation of the cost function and its gradient. To solve this problem, some reduced order methods of high-dimensional or infinite-dimensional hydrodynamic model have been proposed. Among them, the POD method was widely used in flow fields [15–19]. Moreover, POD method had also been successfully combined with 4D VAR data assimilation to solve strong constraint problems for oceanic problems and atmospheric models (Cao et al. [20, 21], Fang et al. [22], and Daescu and Navon [23]).In this paper, based on the concept of independent BFC [10], a 2D tidal dynamical model is modified with the POD reduced order method and then applied to the Bohai, Yellow, and East China Seas (BYECS) where the OBCs are assumed to be temporally and spatially varying. The initial values and OBCs are optimized by assimilating the harmonic constants derived from T
/
P altimeter data. We get the POD reduced order model and its adjoint model. According to the division of the BYECS, we classify the OBCs as 7 types. In our experiment, the absolute mean difference between the simulation results and observations is 3.2 cm and 2.9∘ in amplitude and phase-lag, respectively, and the simulation precision has been improved by 32.7
% and 46.3
%, respectively, compared with results of Lu and Zhang (2006).This paper is organized as follows: a 2D tidal model is described in Section2; Section 3 shows the POD reduced order method; in Section 4, we apply the POD 4D VAR data assimilation to the 2D tidal model; the real experiment is carried out to optimize the initial values and the OBCs in Section 5; finally, we make a summary and draw some conclusions in Section 6.


## 2.
2D Tidal Equations Description

We consider a 2D tidal model as follows:(1)
∂
ζ
∂
t
+
∂
h
+
ζ
u
∂
x
+
∂
h
+
ζ
v
∂
y
=
0
,
∂
u
∂
t
+
g
∂
ζ
∂
x
-
f
v
+
k
u
u
2
+
v
2
h
+
ζ
=
0
,
∂
v
∂
t
+
g
∂
ζ
∂
y
+
f
u
+
k
v
u
2
+
v
2
h
+
ζ
=
0
,where t is time; x and y are Cartesian coordinates, representing the positive eastward directions and northward directions, respectively; h denotes the still water depth; ζ represents the water elevation; u
,
v stand for the eastward and northward components of the horizontal velocity, respectively; f
=
2
Ω
sin
⁡
φ, where Ω is the angular speed of earth rotation and φ is latitude; g is the acceleration due to gravity; k is the bottom friction coefficient.For the sake of convenience, we give the definitions ofC
ζ
(
m
,
n
)
,
C
u
(
m
,
n
)
,
C
v
(
m
,
n
), according to the following ways, where (
m
,
n
) is the grid point:(2)
h
m
,
n
≠
0
⟹
C
ζ
m
,
n
=
1
;
h
m
,
n
=
0
⟹
C
ζ
m
,
n
=
0
;
C
u
m
,
n
=
C
ζ
m
,
n
·
C
ζ
m
+
1
,
n
;
C
v
m
,
n
=
C
ζ
m
,
n
·
C
ζ
m
,
n
+
1
.We set(3)
Γ
=
m
,
n
∣
m
,
n
b
e
l
o
n
g
s
t
o
t
h
e
o
p
e
n
b
o
u
n
d
a
r
y
,
S
ζ
=
m
,
n
∣
C
ζ
m
,
n
=
1
,
m
,
n
∉
Γ
,
S
u
=
m
,
n
∣
C
u
m
,
n
=
1
,
S
v
=
m
,
n
∣
C
v
m
,
n
=
1
.


## 3. Proper Orthogonal Decomposition

LetV represent the model variable (e.g., ζ
,
u
,
v). The ensemble of snapshots sampled at designated time steps V
l
l
=
1
N
=
V
i
l
l
=
1
N, (
1
≤
i
≤
M
), (
N
≤
R
) can be expressed as the following M
×
N matrix A
V: (4)
A
V
=
V
1
1
V
1
2
⋯
V
1
N
V
2
1
V
2
2
⋯
V
2
N
⋮
⋮
⋮
⋮
V
M
1
V
M
2
⋯
V
M
N
,where M is the number of spatial grid points, R represents the number of time steps, and N denotes the number of selected snapshots.The average of the set of snapshotsV
¯
i
i
=
1
M is defined as(5)
V
¯
i
=
1
N
∑
l
=
1
N
V
i
l
,
1
≤
i
≤
M
.Taking the deviation from the mean of the variables yields(6)
V
^
i
l
=
V
i
l
-
V
¯
i
,
1
≤
i
≤
M
,
1
≤
l
≤
N
,which can be expressed by an M
×
N matrix A
=
V
^
i
l.The POD reduced order method can be interpreted as an approach to determine the set of orthogonal basis vectorsϕ
j
∣
j
=
1,2
,
…
,
d, which are optimal in some sense for representing each given snapshot, and the state variable can be represented by the linear combination of the POD basis vectors. The optimization problem becomes the following maximum problem:(7)
m
a
x
1
N
∑
l
=
1
N
V
l
,
∑
j
=
1
d
V
l
,
ϕ
j
ϕ
j
2
,where ·
,
· denotes the inner product between two vectors. Based on this optimization problem and the method of snapshots (see [20] for details), POD basis vectors can be calculated through the solution of the following eigenvalue problem:(8)
A
A
T
x
=
λ
x
.Since the number of spatial grid points may be much larger than the selected snapshots, that is, M
≫
N, we may first solve the eigenvectors φ
j
(
j
=
1
,
…
,
N
) corresponding to A
T
A; that is, (9)
A
T
A
φ
j
=
λ
j
φ
j
,
j
=
1
,
…
,
N
.Next, we obtain eigenvectors ϕ
j
(
j
=
1
,
…
,
N
) corresponding to the nonnull eigenvalues for A
A
T as follows:(10)
ϕ
j
=
1
λ
j
φ
j
,
j
=
1
,
…
,
N
.Finally, we define the relative information content [20], also usually referred to as “energy” I
(
d
), as (11)
I
d
=
∑
i
=
1
d
λ
i
∑
i
=
1
N
λ
iand choose d, which is small enough that I
(
d
) is near to 1, as follows [20]: (12)
d
=
a
r
g
m
i
n
I
d
∣
I
d
≥
σ
100
,where σ is a percentage.Hence, the state variable can be represented by the linear combination of the retained POD basis vectors as follows:(13)
V
i
l
≈
V
~
i
l
=
V
-
i
+
∑
j
=
1
d
a
j
l
ϕ
j
,
i
,where V
i
l is the state variable at spatial node i and temporal node l, ϕ
j
,
i is the value of the jth basis vector ϕ
j at spatial node i, and a
j
l is the POD coefficient at temporal node l corresponding to ϕ
j.


## 4. The Reduced Order Model of 2D Tidal Equations

### 4.1. POD Reduced Order Forward Equations

We assume thatζ
m
,
n
j, u
m
,
n
j, and v
m
,
n
j are the states of the 2D tidal equations (1) (see Appendix A.1) and their corresponding POD reduced order states ζ
~
m
,
n
j, u
~
m
,
n
j, v
~
m
,
n
j are set, respectively, as(14)
ζ
~
m
,
n
j
=
ζ
¯
m
,
n
+
∑
l
=
1
L
ζ
α
l
j
·
ϕ
l
,
m
,
n
ζ
,
u
~
m
,
n
j
=
u
¯
m
,
n
+
∑
l
=
1
L
u
β
l
j
·
ϕ
l
,
m
,
n
u
,
v
~
m
,
n
j
=
v
¯
m
,
n
+
∑
l
=
1
L
v
γ
l
j
·
ϕ
l
,
m
,
n
v
,where ζ
¯
m
,
n, u
¯
m
,
n, and v
¯
m
,
n are the average of the ensemble of snapshots, ϕ
l
,
m
,
n
ζ
l
=
1
L
ζ, ϕ
l
,
m
,
n
u
l
=
1
L
u, and ϕ
l
,
m
,
n
v
l
=
1
L
v are the orthogonal basis vectors of the state of ζ, u, and v, respectively, L
ζ, L
u, L
v are the number of selected POD basis vectors, and α
l
j, β
l
j, γ
l
j are the POD coefficients corresponding to the every POD basis vectors. The initial values α
l
0, β
l
0, γ
l
0 are, respectively, calculated according to(15)
α
l
0
=
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
ζ
m
,
n
0
-
ζ
¯
m
,
n
,
l
=
1
,
…
,
L
ζ
,
β
l
0
=
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
u
m
,
n
0
-
u
¯
m
,
n
,
l
=
1
,
…
,
L
u
,
γ
l
0
=
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
v
m
,
n
0
-
v
¯
m
,
n
,
l
=
1
,
…
,
L
v
.We substitute the POD reduced order states ζ
~
m
,
n
j, u
~
m
,
n
j, and v
~
m
,
n
j into the 2D tidal equations and obtain the POD reduced order forward equations (see Appendix A.2).


### 4.2. POD Reduced Order Adjoint Equations

For the sake of convenience, we assume that there are2
L
t time steps and construct the following cost function:(16)
J
=
1
2
∑
j
=
1
2
L
t
∑
m
,
n
∈
S
ζ
K
m
,
n
ζ
m
,
n
∗
,
j
-
ζ
~
m
,
n
j
2
,where ζ
m
,
n
∗
,
j is the state of the observation at j time step and m
,
n. If there is the observed value at m
,
n, K
m
,
n
=
1; otherwise K
m
,
n
=
0
.Using Lagrangian multiplier method, we construct the following Lagrangian function:(17)
L
=
J
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
ζ
λ
l
2
j
+
1
α
l
2
j
+
1
-
α
l
2
j
Δ
t
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
+
1
/
2
,
n
u
~
m
,
n
2
j
-
h
m
-
1
/
2
,
n
u
~
m
-
1
,
n
2
j
Δ
x
n
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
,
n
+
1
/
2
v
~
m
,
n
2
j
-
h
m
,
n
-
1
/
2
v
~
m
,
n
-
1
2
j
Δ
y
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
v
μ
l
2
j
+
1
γ
l
2
j
+
1
-
γ
l
2
j
Δ
t
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
k
m
,
n
+
1
/
2
s
~
m
,
n
2
j
ε
v
~
m
,
n
2
j
+
1
+
1
-
ε
v
~
m
,
n
,
2
j
h
m
,
n
+
1
/
2
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
f
n
u
~
m
,
n
2
j
¯
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
+
1
2
j
+
1
Δ
y
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
+
1
2
j
+
1
Δ
y
∑
m
,
n
∈
S
v
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
2
j
+
1
Δ
y
-
∑
m
,
n
∈
S
v
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
2
j
+
1
Δ
y
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
u
η
l
2
j
+
1
β
l
2
j
+
1
-
β
l
2
j
Δ
t
+
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
k
m
+
1
/
2
,
n
r
~
m
,
n
2
j
ε
u
~
m
,
n
2
j
+
1
+
1
-
ε
u
~
m
,
n
2
j
h
m
+
1
/
2
,
n
-
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
f
n
v
~
m
,
n
2
j
+
1
¯
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
+
1
,
n
2
j
+
1
Δ
x
n
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
+
1
,
n
2
j
+
1
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
,
n
2
j
+
1
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
,
n
2
j
+
1
Δ
x
n
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
ζ
λ
l
2
j
+
2
α
l
2
j
+
2
-
α
l
2
j
+
1
Δ
t
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
+
1
/
2
,
n
u
~
m
,
n
2
j
+
1
-
h
m
-
1
/
2
,
n
u
~
m
-
1
,
n
2
j
+
1
Δ
x
n
-
1
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
,
n
+
1
/
2
v
~
m
,
n
2
j
+
1
-
h
m
,
n
-
1
/
2
v
~
m
,
n
-
1
2
j
+
1
Δ
y
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
u
η
l
2
j
+
2
β
l
2
j
+
2
-
β
l
2
j
+
1
Δ
t
+
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
k
m
+
1
/
2
,
n
r
~
m
,
n
2
j
+
1
ε
u
~
m
,
n
2
j
+
2
+
1
-
ε
u
~
m
,
n
2
j
+
1
h
m
+
1
/
2
,
n
-
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
f
n
v
~
m
,
n
2
j
+
1
¯
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
+
1
,
n
2
j
+
2
Δ
x
n
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
+
1
,
n
2
j
+
2
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
,
n
2
j
+
2
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
,
n
2
j
+
2
Δ
x
n
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
v
μ
l
2
j
+
2
γ
l
2
j
+
2
-
γ
l
2
j
+
1
Δ
t
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
k
m
,
n
+
1
/
2
s
~
m
,
n
2
j
+
1
ε
v
~
m
,
n
2
j
+
2
+
1
-
ε
v
~
m
,
n
2
j
+
1
h
m
,
n
+
1
/
2
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
f
n
u
~
m
,
n
2
j
+
2
¯
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
+
1
2
j
+
2
Δ
y
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
+
1
2
j
+
2
Δ
y
-
∑
m
,
n
∈
S
v
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
2
j
+
2
Δ
y
-
∑
m
,
n
∈
S
v
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
2
j
+
2
Δ
y
,
where λ
l
j, η
l
j, and μ
l
j are the adjoint variables of α
l
j, β
l
j, and γ
l
j, respectively.In order to minimize the cost function, the Lagrangian function should satisfy the following equations:(18)
∂
L
∂
λ
l
j
=
0
,
l
=
1
,
…
,
L
ζ
,
j
=
1
,
…
,
2
L
t
;
(19)
∂
L
∂
η
l
j
=
0
,
l
=
1
,
…
,
L
u
,
j
=
1
,
…
,
2
L
t
;
(20)
∂
L
∂
μ
l
j
=
0
,
l
=
1
,
…
,
L
v
,
j
=
1
,
…
,
2
L
t
;
(21)
∂
L
∂
α
l
j
=
0
,
l
=
1
,
…
,
L
ζ
,
j
=
1
,
…
,
2
L
t
;
(22)
∂
L
∂
β
l
j
=
0
,
l
=
1
,
…
,
L
u
,
j
=
1
,
…
,
2
L
t
;
(23)
∂
L
∂
γ
l
j
=
0
,
l
=
1
,
…
,
L
v
,
j
=
1
,
…
,
2
L
t
;
(24)
∂
L
∂
α
l
0
=
0
,
l
=
1
,
…
,
L
ζ
;
(25)
∂
L
∂
β
l
0
=
0
,
l
=
1
,
…
,
L
u
;
(26)
∂
L
∂
γ
l
0
=
0
,
l
=
1
,
…
,
L
v
;
(27)
∂
L
∂
ζ
m
,
n
j
=
0
,
m
,
n
∈
Γ
,
j
=
1
,
…
,
2
L
t
-
1
.Equations (18)–(20) return to the POD reduced order forward equations (see Appendix A.2). By solving (21)–(23), the POD reduced order adjoint equations are obtained (see Appendix A.3). By solving (24)–(26), the gradients of the cost function with respect to α
l
0, β
l
0, and γ
l
0 can be calculated as follows:(28)
∂
J
∂
α
i
0
=
-
λ
i
1
,
i
=
1
,
…
,
L
ζ
;
(29)
∂
J
∂
β
i
0
=
-
η
i
1
+
∑
m
,
n
∈
S
ζ
Δ
t
h
m
+
1
/
2
,
n
ϕ
i
,
m
,
n
u
-
h
m
-
1
/
2
,
n
ϕ
i
,
m
-
1
,
n
u
∑
l
=
1
L
ζ
λ
l
1
ϕ
l
,
m
,
n
ζ
Δ
x
n
-
1
+
∑
m
,
n
∈
S
v
ϕ
i
,
m
,
n
u
¯
Δ
t
k
m
,
n
+
1
/
2
u
~
m
,
n
0
¯
ε
v
~
m
,
n
1
+
1
-
ε
v
~
m
,
n
0
+
f
n
h
m
,
n
+
1
/
2
s
~
m
,
n
0
∑
l
=
1
L
v
μ
l
1
ϕ
l
,
m
,
n
v
h
m
,
n
+
1
/
2
s
~
m
,
n
0
+
∑
m
,
n
∈
S
u
ϕ
i
,
m
,
n
u
Δ
t
k
m
+
1
/
2
,
n
u
~
m
,
n
0
ε
u
¯
m
,
n
+
1
-
ε
u
~
m
,
n
0
+
1
-
ε
r
~
m
,
n
0
2
∑
l
=
1
L
u
η
l
1
ϕ
l
,
m
,
n
u
h
m
+
1
/
2
,
n
r
~
m
,
n
0
,
(30)
∂
J
∂
γ
i
0
=
-
μ
i
1
+
∑
m
,
n
∈
S
ζ
Δ
t
h
m
,
n
+
1
/
2
ϕ
i
,
m
,
n
v
-
h
m
,
n
-
1
/
2
ϕ
i
,
m
,
n
-
1
v
∑
l
=
1
L
ζ
λ
l
1
ϕ
l
,
m
,
n
ζ
Δ
y
+
∑
m
,
n
∈
S
v
ϕ
i
,
m
,
n
v
Δ
t
k
m
,
n
+
1
/
2
v
~
m
,
n
0
ε
v
~
m
,
n
1
+
1
-
ε
v
~
m
,
n
0
+
1
-
ε
s
~
m
,
n
0
2
∑
l
=
1
L
v
μ
l
1
ϕ
l
,
m
,
n
v
h
m
,
n
+
1
/
2
s
~
m
,
n
0
+
∑
m
,
n
∈
S
u
ϕ
i
,
m
,
n
v
¯
Δ
t
k
m
+
1
/
2
,
n
v
~
m
,
n
0
¯
ε
u
¯
m
,
n
+
1
-
ε
u
~
m
,
n
0
∑
l
=
1
L
u
η
l
1
ϕ
l
,
m
,
n
u
h
m
+
1
/
2
,
n
r
~
m
,
n
0
.By solving (27), the gradients of the cost function with respect to OBCs are obtained. For m
,
n belonging to the open boundary,(31)
i
f
m
,
n
∈
S
u
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
-
∑
l
=
1
L
u
η
l
j
ϕ
l
,
m
,
n
u
Δ
t
g
Δ
x
n
;
(32)
i
f
m
-
1
,
n
∈
S
u
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
∑
l
=
1
L
u
η
l
j
ϕ
l
,
m
-
1
,
n
u
Δ
t
g
Δ
x
n
-
1
;
(33)
i
f
m
,
n
∈
S
v
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
-
∑
l
=
1
L
v
μ
l
j
ϕ
l
,
m
,
n
v
Δ
t
g
Δ
y
;
(34)
i
f
m
,
n
-
1
∈
S
v
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
∑
l
=
1
L
v
μ
l
j
ϕ
l
,
m
,
n
-
1
v
Δ
t
g
Δ
y
;if m
,
n belongs to the above variety of situations, the gradient of the cost function with respect to OBC is equal to the sum of the corresponding gradients.Using (28)–(34), the initial values and OBCs are adjusted, respectively, as follows:(35)
α
^
i
0
=
α
i
0
-
ρ
α
i
0
∂
J
∂
α
i
0
,
i
=
1
,
…
,
L
ζ
;
(36)
β
^
i
0
=
β
i
0
-
ρ
β
i
0
∂
J
∂
β
i
0
,
i
=
1
,
…
,
L
u
;
γ
^
i
0
=
γ
i
0
-
ρ
γ
i
0
∂
J
∂
γ
i
0
,
i
=
1
,
…
,
L
v
;
i
f
m
,
n
∈
Γ
,
t
h
e
n
ζ
^
m
,
n
j
=
ζ
m
,
n
j
-
ρ
ζ
m
,
n
j
∂
J
∂
ζ
m
,
n
j
,
j
=
0
,
…
,
2
L
t
-
1
,where ρ is the step length and α
^
i
0, β
^
i
0, γ
^
i
0, ζ
^
m
,
n
j are the updated values of α
i
0, β
i
0, γ
i
0, ζ
m
,
n
j, respectively.


## 4.1. POD Reduced Order Forward Equations

We assume thatζ
m
,
n
j, u
m
,
n
j, and v
m
,
n
j are the states of the 2D tidal equations (1) (see Appendix A.1) and their corresponding POD reduced order states ζ
~
m
,
n
j, u
~
m
,
n
j, v
~
m
,
n
j are set, respectively, as(14)
ζ
~
m
,
n
j
=
ζ
¯
m
,
n
+
∑
l
=
1
L
ζ
α
l
j
·
ϕ
l
,
m
,
n
ζ
,
u
~
m
,
n
j
=
u
¯
m
,
n
+
∑
l
=
1
L
u
β
l
j
·
ϕ
l
,
m
,
n
u
,
v
~
m
,
n
j
=
v
¯
m
,
n
+
∑
l
=
1
L
v
γ
l
j
·
ϕ
l
,
m
,
n
v
,where ζ
¯
m
,
n, u
¯
m
,
n, and v
¯
m
,
n are the average of the ensemble of snapshots, ϕ
l
,
m
,
n
ζ
l
=
1
L
ζ, ϕ
l
,
m
,
n
u
l
=
1
L
u, and ϕ
l
,
m
,
n
v
l
=
1
L
v are the orthogonal basis vectors of the state of ζ, u, and v, respectively, L
ζ, L
u, L
v are the number of selected POD basis vectors, and α
l
j, β
l
j, γ
l
j are the POD coefficients corresponding to the every POD basis vectors. The initial values α
l
0, β
l
0, γ
l
0 are, respectively, calculated according to(15)
α
l
0
=
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
ζ
m
,
n
0
-
ζ
¯
m
,
n
,
l
=
1
,
…
,
L
ζ
,
β
l
0
=
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
u
m
,
n
0
-
u
¯
m
,
n
,
l
=
1
,
…
,
L
u
,
γ
l
0
=
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
v
m
,
n
0
-
v
¯
m
,
n
,
l
=
1
,
…
,
L
v
.We substitute the POD reduced order states ζ
~
m
,
n
j, u
~
m
,
n
j, and v
~
m
,
n
j into the 2D tidal equations and obtain the POD reduced order forward equations (see Appendix A.2).


## 4.2. POD Reduced Order Adjoint Equations

For the sake of convenience, we assume that there are2
L
t time steps and construct the following cost function:(16)
J
=
1
2
∑
j
=
1
2
L
t
∑
m
,
n
∈
S
ζ
K
m
,
n
ζ
m
,
n
∗
,
j
-
ζ
~
m
,
n
j
2
,where ζ
m
,
n
∗
,
j is the state of the observation at j time step and m
,
n. If there is the observed value at m
,
n, K
m
,
n
=
1; otherwise K
m
,
n
=
0
.Using Lagrangian multiplier method, we construct the following Lagrangian function:(17)
L
=
J
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
ζ
λ
l
2
j
+
1
α
l
2
j
+
1
-
α
l
2
j
Δ
t
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
+
1
/
2
,
n
u
~
m
,
n
2
j
-
h
m
-
1
/
2
,
n
u
~
m
-
1
,
n
2
j
Δ
x
n
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
,
n
+
1
/
2
v
~
m
,
n
2
j
-
h
m
,
n
-
1
/
2
v
~
m
,
n
-
1
2
j
Δ
y
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
v
μ
l
2
j
+
1
γ
l
2
j
+
1
-
γ
l
2
j
Δ
t
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
k
m
,
n
+
1
/
2
s
~
m
,
n
2
j
ε
v
~
m
,
n
2
j
+
1
+
1
-
ε
v
~
m
,
n
,
2
j
h
m
,
n
+
1
/
2
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
f
n
u
~
m
,
n
2
j
¯
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
+
1
2
j
+
1
Δ
y
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
+
1
2
j
+
1
Δ
y
∑
m
,
n
∈
S
v
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
2
j
+
1
Δ
y
-
∑
m
,
n
∈
S
v
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
2
j
+
1
Δ
y
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
u
η
l
2
j
+
1
β
l
2
j
+
1
-
β
l
2
j
Δ
t
+
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
k
m
+
1
/
2
,
n
r
~
m
,
n
2
j
ε
u
~
m
,
n
2
j
+
1
+
1
-
ε
u
~
m
,
n
2
j
h
m
+
1
/
2
,
n
-
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
f
n
v
~
m
,
n
2
j
+
1
¯
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
+
1
,
n
2
j
+
1
Δ
x
n
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
+
1
,
n
2
j
+
1
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
,
n
2
j
+
1
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
,
n
2
j
+
1
Δ
x
n
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
ζ
λ
l
2
j
+
2
α
l
2
j
+
2
-
α
l
2
j
+
1
Δ
t
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
+
1
/
2
,
n
u
~
m
,
n
2
j
+
1
-
h
m
-
1
/
2
,
n
u
~
m
-
1
,
n
2
j
+
1
Δ
x
n
-
1
+
∑
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
ζ
h
m
,
n
+
1
/
2
v
~
m
,
n
2
j
+
1
-
h
m
,
n
-
1
/
2
v
~
m
,
n
-
1
2
j
+
1
Δ
y
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
u
η
l
2
j
+
2
β
l
2
j
+
2
-
β
l
2
j
+
1
Δ
t
+
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
k
m
+
1
/
2
,
n
r
~
m
,
n
2
j
+
1
ε
u
~
m
,
n
2
j
+
2
+
1
-
ε
u
~
m
,
n
2
j
+
1
h
m
+
1
/
2
,
n
-
∑
m
,
n
∈
S
u
ϕ
l
,
m
,
n
u
f
n
v
~
m
,
n
2
j
+
1
¯
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
+
1
,
n
2
j
+
2
Δ
x
n
+
∑
m
,
n
∈
S
u
,
m
+
1
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
+
1
,
n
2
j
+
2
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
~
m
,
n
2
j
+
2
Δ
x
n
-
∑
m
,
n
∈
S
u
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
u
g
ζ
m
,
n
2
j
+
2
Δ
x
n
+
∑
j
=
0
L
t
-
1
∑
l
=
1
L
v
μ
l
2
j
+
2
γ
l
2
j
+
2
-
γ
l
2
j
+
1
Δ
t
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
k
m
,
n
+
1
/
2
s
~
m
,
n
2
j
+
1
ε
v
~
m
,
n
2
j
+
2
+
1
-
ε
v
~
m
,
n
2
j
+
1
h
m
,
n
+
1
/
2
+
∑
m
,
n
∈
S
v
ϕ
l
,
m
,
n
v
f
n
u
~
m
,
n
2
j
+
2
¯
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
+
1
2
j
+
2
Δ
y
+
∑
m
,
n
∈
S
v
,
m
,
n
+
1
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
+
1
2
j
+
2
Δ
y
-
∑
m
,
n
∈
S
v
,
m
,
n
∈
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
~
m
,
n
2
j
+
2
Δ
y
-
∑
m
,
n
∈
S
v
,
m
,
n
∉
S
ζ
ϕ
l
,
m
,
n
v
g
ζ
m
,
n
2
j
+
2
Δ
y
,
where λ
l
j, η
l
j, and μ
l
j are the adjoint variables of α
l
j, β
l
j, and γ
l
j, respectively.In order to minimize the cost function, the Lagrangian function should satisfy the following equations:(18)
∂
L
∂
λ
l
j
=
0
,
l
=
1
,
…
,
L
ζ
,
j
=
1
,
…
,
2
L
t
;
(19)
∂
L
∂
η
l
j
=
0
,
l
=
1
,
…
,
L
u
,
j
=
1
,
…
,
2
L
t
;
(20)
∂
L
∂
μ
l
j
=
0
,
l
=
1
,
…
,
L
v
,
j
=
1
,
…
,
2
L
t
;
(21)
∂
L
∂
α
l
j
=
0
,
l
=
1
,
…
,
L
ζ
,
j
=
1
,
…
,
2
L
t
;
(22)
∂
L
∂
β
l
j
=
0
,
l
=
1
,
…
,
L
u
,
j
=
1
,
…
,
2
L
t
;
(23)
∂
L
∂
γ
l
j
=
0
,
l
=
1
,
…
,
L
v
,
j
=
1
,
…
,
2
L
t
;
(24)
∂
L
∂
α
l
0
=
0
,
l
=
1
,
…
,
L
ζ
;
(25)
∂
L
∂
β
l
0
=
0
,
l
=
1
,
…
,
L
u
;
(26)
∂
L
∂
γ
l
0
=
0
,
l
=
1
,
…
,
L
v
;
(27)
∂
L
∂
ζ
m
,
n
j
=
0
,
m
,
n
∈
Γ
,
j
=
1
,
…
,
2
L
t
-
1
.Equations (18)–(20) return to the POD reduced order forward equations (see Appendix A.2). By solving (21)–(23), the POD reduced order adjoint equations are obtained (see Appendix A.3). By solving (24)–(26), the gradients of the cost function with respect to α
l
0, β
l
0, and γ
l
0 can be calculated as follows:(28)
∂
J
∂
α
i
0
=
-
λ
i
1
,
i
=
1
,
…
,
L
ζ
;
(29)
∂
J
∂
β
i
0
=
-
η
i
1
+
∑
m
,
n
∈
S
ζ
Δ
t
h
m
+
1
/
2
,
n
ϕ
i
,
m
,
n
u
-
h
m
-
1
/
2
,
n
ϕ
i
,
m
-
1
,
n
u
∑
l
=
1
L
ζ
λ
l
1
ϕ
l
,
m
,
n
ζ
Δ
x
n
-
1
+
∑
m
,
n
∈
S
v
ϕ
i
,
m
,
n
u
¯
Δ
t
k
m
,
n
+
1
/
2
u
~
m
,
n
0
¯
ε
v
~
m
,
n
1
+
1
-
ε
v
~
m
,
n
0
+
f
n
h
m
,
n
+
1
/
2
s
~
m
,
n
0
∑
l
=
1
L
v
μ
l
1
ϕ
l
,
m
,
n
v
h
m
,
n
+
1
/
2
s
~
m
,
n
0
+
∑
m
,
n
∈
S
u
ϕ
i
,
m
,
n
u
Δ
t
k
m
+
1
/
2
,
n
u
~
m
,
n
0
ε
u
¯
m
,
n
+
1
-
ε
u
~
m
,
n
0
+
1
-
ε
r
~
m
,
n
0
2
∑
l
=
1
L
u
η
l
1
ϕ
l
,
m
,
n
u
h
m
+
1
/
2
,
n
r
~
m
,
n
0
,
(30)
∂
J
∂
γ
i
0
=
-
μ
i
1
+
∑
m
,
n
∈
S
ζ
Δ
t
h
m
,
n
+
1
/
2
ϕ
i
,
m
,
n
v
-
h
m
,
n
-
1
/
2
ϕ
i
,
m
,
n
-
1
v
∑
l
=
1
L
ζ
λ
l
1
ϕ
l
,
m
,
n
ζ
Δ
y
+
∑
m
,
n
∈
S
v
ϕ
i
,
m
,
n
v
Δ
t
k
m
,
n
+
1
/
2
v
~
m
,
n
0
ε
v
~
m
,
n
1
+
1
-
ε
v
~
m
,
n
0
+
1
-
ε
s
~
m
,
n
0
2
∑
l
=
1
L
v
μ
l
1
ϕ
l
,
m
,
n
v
h
m
,
n
+
1
/
2
s
~
m
,
n
0
+
∑
m
,
n
∈
S
u
ϕ
i
,
m
,
n
v
¯
Δ
t
k
m
+
1
/
2
,
n
v
~
m
,
n
0
¯
ε
u
¯
m
,
n
+
1
-
ε
u
~
m
,
n
0
∑
l
=
1
L
u
η
l
1
ϕ
l
,
m
,
n
u
h
m
+
1
/
2
,
n
r
~
m
,
n
0
.By solving (27), the gradients of the cost function with respect to OBCs are obtained. For m
,
n belonging to the open boundary,(31)
i
f
m
,
n
∈
S
u
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
-
∑
l
=
1
L
u
η
l
j
ϕ
l
,
m
,
n
u
Δ
t
g
Δ
x
n
;
(32)
i
f
m
-
1
,
n
∈
S
u
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
∑
l
=
1
L
u
η
l
j
ϕ
l
,
m
-
1
,
n
u
Δ
t
g
Δ
x
n
-
1
;
(33)
i
f
m
,
n
∈
S
v
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
-
∑
l
=
1
L
v
μ
l
j
ϕ
l
,
m
,
n
v
Δ
t
g
Δ
y
;
(34)
i
f
m
,
n
-
1
∈
S
v
,
t
h
e
n
∂
J
∂
ζ
m
,
n
j
=
∑
l
=
1
L
v
μ
l
j
ϕ
l
,
m
,
n
-
1
v
Δ
t
g
Δ
y
;if m
,
n belongs to the above variety of situations, the gradient of the cost function with respect to OBC is equal to the sum of the corresponding gradients.Using (28)–(34), the initial values and OBCs are adjusted, respectively, as follows:(35)
α
^
i
0
=
α
i
0
-
ρ
α
i
0
∂
J
∂
α
i
0
,
i
=
1
,
…
,
L
ζ
;
(36)
β
^
i
0
=
β
i
0
-
ρ
β
i
0
∂
J
∂
β
i
0
,
i
=
1
,
…
,
L
u
;
γ
^
i
0
=
γ
i
0
-
ρ
γ
i
0
∂
J
∂
γ
i
0
,
i
=
1
,
…
,
L
v
;
i
f
m
,
n
∈
Γ
,
t
h
e
n
ζ
^
m
,
n
j
=
ζ
m
,
n
j
-
ρ
ζ
m
,
n
j
∂
J
∂
ζ
m
,
n
j
,
j
=
0
,
…
,
2
L
t
-
1
,where ρ is the step length and α
^
i
0, β
^
i
0, γ
^
i
0, ζ
^
m
,
n
j are the updated values of α
i
0, β
i
0, γ
i
0, ζ
m
,
n
j, respectively.


## 5. Numerical Experiment and Results

4D VAR data assimilation is a powerful tool for parameter estimation. The cost function measuring the data misfit between the model output and observations is minimized through optimizing the parameters. This means that the physical process in oceans is simulated while the observation data is assimilated into the hydrodynamic model. In this section, the performance and feasibility of the POD 4D VAR data assimilation when applied to estimating parameters such as the initial values and OBCs are verified through the real experiment.


### 5.1. The Design of Experiment

In the present study, the computing area is BYECS (117
.
5
∘
E
–
13
1
∘
E, 2
4
∘
N
–
4
1
∘
N, Figure 1). From Figure 1, we can find that the water is shallower in the Bohai Sea and the Yellow Sea and deeper in the East China Sea. The deepest water appears in and around the Okinawa Trough of the East China Sea. The open boundaries and the T
/
P satellite tracks are also shown in Figure 1. The satellite tracks contain six descending tracks (036, 062, 138, 164, 214, and 240) and six ascending tracks (051, 077, 127, 153, 203, and 229). The open boundaries are set along the Taiwan Strait, the first island chains, the Korea Strait, and the Tsushima Strait. The horizontal resolution is 10
′
×
10
′. The angular frequency of M
2 tide is 1.405189
×
1
0
-
4
s
-
4 and the time step is 62.103 s (1/720 of the period of M
2 tide). The initial conditions are ζ
=
u
=
v
=
0 at t
=
0. The closed boundary condition is that the normal velocity is zero. BFC and OBCs are derived from the results of Lu and Zhang (2006).Figure 1
The bathymetric map of BYECS and the position ofT
/
P satellite tracks and open boundaries. The solid lines are T
/
P satellite tracks. The symbol “+” stands for the open boundary.According to the division of BYECS, we classify the OBCs into 7 types, for a grid pointm
,
n belonging to the open boundary:(1)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the first type Γ
1; that is, m
,
n
∈
Γ
1.
(2)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
0, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the second type Γ
2; that is, m
,
n
∈
Γ
2.
(3)
IfC
u
m
,
n
=
0, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the third type Γ
3; that is, m
,
n
∈
Γ
3.
(4)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
0, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the fourth type Γ
4; that is, m
,
n
∈
Γ
4.
(5)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
0, C
v
m
,
n
=
0, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the fifth type Γ
5; that is, m
,
n
∈
Γ
5.
(6)
IfC
u
m
,
n
=
0, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
0, m
,
n belongs to the sixth type Γ
6; that is, m
,
n
∈
Γ
6.
(7)
IfC
u
m
,
n
=
0, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
0, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the seventh type Γ
7; that is, m
,
n
∈
Γ
7.The number of grid points belonging toΓ
i is denoted by N
i
Γ, i
=
1
,
…
,
7.The 2D tidal models were run 12 times when they reached steady states, and we choose the states every 7 times as the snapshots. The numbers of POD basis vectorsL
ζ
=
L
u
=
L
v
=
6 were chosen to maintain 99% of energy as stated in (12).


### 5.2. Optimization Steps

The initial valuesα
i
0, β
i
0, γ
i
0 and OBCs can be optimized with the following steps.Step 1.
Calculateα
i
0, β
i
0, γ
i
0 by (15).Step 2.
Integrate the POD reduced order forward equations and adjoint equations.Step 3.
Calculate the cost function and its gradients with respect to the control variables.Step 4.
Optimizeα
i
0, β
i
0, γ
i
0 and OBCs by (35).Step 5.
Substitute the optimized parameters into the POD reduced order forward equations and then perform them.Step 6.
Calculate the POD reduced order state variables by (14).Step 7.
If some stopping criterion is met, stop and return the optimized parameters; otherwise, replace the control parameters with the new values and repeat Steps  2–6.The algorithm shows that the process optimizing the parameters is the process assimilating the observation data. In this way, the combination of the hydrodynamic model and the observations is realized in this model.Iterative minimization is stopped once a convergence criterion is met. The criterion could be the cost function being sufficiently small, the last two values of the cost function being sufficiently close, the gradient being sufficiently small, or a combination of these. However, for the sake of convenience, we give the optimization result after 200 iteration steps in the experiment.


### 5.3. Experimental Results and Discussion

We defines representing the simulated value, o representing the observed value, and N
o representing the number of observation sites. Using harmonic analysis, we gain the phase-lag θ
s and amplitude Z
s based on the POD reduced order tidal elevations. The deviations are used to represent the discrepancy between the simulated and observed harmonic constants, which are calculated by the following formulas.Absolute mean error of amplitude is as follows:(37)
Δ
Z
=
1
N
o
∑
n
=
1
N
o
Z
s
,
n
-
Z
o
,
n
.Absolute mean error of phase-lag is as follows: (38)
Δ
θ
=
1
N
o
∑
n
=
1
N
o
θ
s
,
n
-
θ
o
,
n
.Lu and Zhang [10] studied the spatially varying bottom friction coefficient (BFC) by using 4D VAR data assimilation. In the EXP-PE4 (see [10] for details), the absolute mean difference of amplitude and phase-lag is 7.3 cm and 6.3∘, respectively. In the EXP-PE2 (see [10] for details), the absolute mean difference of amplitude and phase-lag is 5.0 cm and 5.4∘, respectively. Based on the results of independent BFC [10], we invert the initial values and the OBCs using the POD 4D VAR data assimilation and the absolute mean difference of amplitude and phase-lag is 3.2 cm and 2.9∘, respectively. The absolute mean difference of amplitude and phase-lag is declined by 32.7% and 46.3% compared to the results of EXP-PE2 (see [10]), respectively. Compared with the results of Lu and Zhang (2006), the results of the experiment in this paper appear to coincide with the observed M
2 tide in BYECS fairly well (Figure 2). The results of the absolute mean difference and Figure 2 are satisfactory and indicate that the observations have been successfully assimilated into the model.Figure 2
The computed cotidal charts ofM
2 tide. The dotted and solid lines denote the amplitude and phase-lag (deg.), respectively.The cost function and its gradients with respect to the control variables are two important convergence criteria for data assimilation in this model. Figure3 shows the cost values J (normalized by their initial value J
0). Figures 4–6 show the gradients of the cost function J with respect to α
i
0, β
i
0, γ
i
0, (normalized by their first values ∂
J
1
/
∂
α
i
0, ∂
J
1
/
∂
β
i
0, ∂
J
1
/
∂
γ
i
0, where J
1 is the cost value J of the first iteration), respectively, versus the number of iterations. We define the gradients of the cost function J with respect to the seven types of the OBCs as follows:(39)
∂
J
∂
Γ
i
=
1
2
L
t
-
1
N
i
Γ
∑
j
=
1
2
L
t
-
1
∑
m
,
n
∈
Γ
i
∂
J
∂
ζ
m
,
n
j
2
,
i
=
1
,
…
,
7
.Figure 7 shows the gradients of the cost function J with respect to the seven types of the OBCs (normalized by their first values ∂
J
1
/
∂
Γ
i, i
=
1
,
…
,
7). From Figures 3–7, we can see that both the cost function and gradients drop dramatically in the optimization process, which demonstrate the strong inversion ability of this model.Figure 3
The valuesJ
/
J
0 versus iteration number, where J
0 represents the initial value of the cost function.Figure 4
The curveF
i represents the gradients of J with respect to α
i
0 versus iteration number, i
=
1
,
…
,
6, respectively.Figure 5
The curveP
i represents the gradients of J with respect to β
i
0 versus iteration number, i
=
1
,
…
,
6, respectively.Figure 6
The curveZ
i represents the gradients of J with respect to γ
i
0 versus iteration number, i
=
1
,
…
,
6, respectively.Figure 7
The curveG
(
i
) represents the gradients of J with respect to ith type of the OBCs versus iteration number, i
=
1
,
…
,
7, respectively.


## 5.1. The Design of Experiment

In the present study, the computing area is BYECS (117
.
5
∘
E
–
13
1
∘
E, 2
4
∘
N
–
4
1
∘
N, Figure 1). From Figure 1, we can find that the water is shallower in the Bohai Sea and the Yellow Sea and deeper in the East China Sea. The deepest water appears in and around the Okinawa Trough of the East China Sea. The open boundaries and the T
/
P satellite tracks are also shown in Figure 1. The satellite tracks contain six descending tracks (036, 062, 138, 164, 214, and 240) and six ascending tracks (051, 077, 127, 153, 203, and 229). The open boundaries are set along the Taiwan Strait, the first island chains, the Korea Strait, and the Tsushima Strait. The horizontal resolution is 10
′
×
10
′. The angular frequency of M
2 tide is 1.405189
×
1
0
-
4
s
-
4 and the time step is 62.103 s (1/720 of the period of M
2 tide). The initial conditions are ζ
=
u
=
v
=
0 at t
=
0. The closed boundary condition is that the normal velocity is zero. BFC and OBCs are derived from the results of Lu and Zhang (2006).Figure 1
The bathymetric map of BYECS and the position ofT
/
P satellite tracks and open boundaries. The solid lines are T
/
P satellite tracks. The symbol “+” stands for the open boundary.According to the division of BYECS, we classify the OBCs into 7 types, for a grid pointm
,
n belonging to the open boundary:(1)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the first type Γ
1; that is, m
,
n
∈
Γ
1.
(2)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
0, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the second type Γ
2; that is, m
,
n
∈
Γ
2.
(3)
IfC
u
m
,
n
=
0, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the third type Γ
3; that is, m
,
n
∈
Γ
3.
(4)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
0, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the fourth type Γ
4; that is, m
,
n
∈
Γ
4.
(5)
IfC
u
m
,
n
=
1, C
u
m
-
1
,
n
=
0, C
v
m
,
n
=
0, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the fifth type Γ
5; that is, m
,
n
∈
Γ
5.
(6)
IfC
u
m
,
n
=
0, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
1, and C
v
m
,
n
-
1
=
0, m
,
n belongs to the sixth type Γ
6; that is, m
,
n
∈
Γ
6.
(7)
IfC
u
m
,
n
=
0, C
u
m
-
1
,
n
=
1, C
v
m
,
n
=
0, and C
v
m
,
n
-
1
=
1, m
,
n belongs to the seventh type Γ
7; that is, m
,
n
∈
Γ
7.The number of grid points belonging toΓ
i is denoted by N
i
Γ, i
=
1
,
…
,
7.The 2D tidal models were run 12 times when they reached steady states, and we choose the states every 7 times as the snapshots. The numbers of POD basis vectorsL
ζ
=
L
u
=
L
v
=
6 were chosen to maintain 99% of energy as stated in (12).


## 5.2. Optimization Steps

The initial valuesα
i
0, β
i
0, γ
i
0 and OBCs can be optimized with the following steps.Step 1.
Calculateα
i
0, β
i
0, γ
i
0 by (15).Step 2.
Integrate the POD reduced order forward equations and adjoint equations.Step 3.
Calculate the cost function and its gradients with respect to the control variables.Step 4.
Optimizeα
i
0, β
i
0, γ
i
0 and OBCs by (35).Step 5.
Substitute the optimized parameters into the POD reduced order forward equations and then perform them.Step 6.
Calculate the POD reduced order state variables by (14).Step 7.
If some stopping criterion is met, stop and return the optimized parameters; otherwise, replace the control parameters with the new values and repeat Steps  2–6.The algorithm shows that the process optimizing the parameters is the process assimilating the observation data. In this way, the combination of the hydrodynamic model and the observations is realized in this model.Iterative minimization is stopped once a convergence criterion is met. The criterion could be the cost function being sufficiently small, the last two values of the cost function being sufficiently close, the gradient being sufficiently small, or a combination of these. However, for the sake of convenience, we give the optimization result after 200 iteration steps in the experiment.


## 5.3. Experimental Results and Discussion

We defines representing the simulated value, o representing the observed value, and N
o representing the number of observation sites. Using harmonic analysis, we gain the phase-lag θ
s and amplitude Z
s based on the POD reduced order tidal elevations. The deviations are used to represent the discrepancy between the simulated and observed harmonic constants, which are calculated by the following formulas.Absolute mean error of amplitude is as follows:(37)
Δ
Z
=
1
N
o
∑
n
=
1
N
o
Z
s
,
n
-
Z
o
,
n
.Absolute mean error of phase-lag is as follows: (38)
Δ
θ
=
1
N
o
∑
n
=
1
N
o
θ
s
,
n
-
θ
o
,
n
.Lu and Zhang [10] studied the spatially varying bottom friction coefficient (BFC) by using 4D VAR data assimilation. In the EXP-PE4 (see [10] for details), the absolute mean difference of amplitude and phase-lag is 7.3 cm and 6.3∘, respectively. In the EXP-PE2 (see [10] for details), the absolute mean difference of amplitude and phase-lag is 5.0 cm and 5.4∘, respectively. Based on the results of independent BFC [10], we invert the initial values and the OBCs using the POD 4D VAR data assimilation and the absolute mean difference of amplitude and phase-lag is 3.2 cm and 2.9∘, respectively. The absolute mean difference of amplitude and phase-lag is declined by 32.7% and 46.3% compared to the results of EXP-PE2 (see [10]), respectively. Compared with the results of Lu and Zhang (2006), the results of the experiment in this paper appear to coincide with the observed M
2 tide in BYECS fairly well (Figure 2). The results of the absolute mean difference and Figure 2 are satisfactory and indicate that the observations have been successfully assimilated into the model.Figure 2
The computed cotidal charts ofM
2 tide. The dotted and solid lines denote the amplitude and phase-lag (deg.), respectively.The cost function and its gradients with respect to the control variables are two important convergence criteria for data assimilation in this model. Figure3 shows the cost values J (normalized by their initial value J
0). Figures 4–6 show the gradients of the cost function J with respect to α
i
0, β
i
0, γ
i
0, (normalized by their first values ∂
J
1
/
∂
α
i
0, ∂
J
1
/
∂
β
i
0, ∂
J
1
/
∂
γ
i
0, where J
1 is the cost value J of the first iteration), respectively, versus the number of iterations. We define the gradients of the cost function J with respect to the seven types of the OBCs as follows:(39)
∂
J
∂
Γ
i
=
1
2
L
t
-
1
N
i
Γ
∑
j
=
1
2
L
t
-
1
∑
m
,
n
∈
Γ
i
∂
J
∂
ζ
m
,
n
j
2
,
i
=
1
,
…
,
7
.Figure 7 shows the gradients of the cost function J with respect to the seven types of the OBCs (normalized by their first values ∂
J
1
/
∂
Γ
i, i
=
1
,
…
,
7). From Figures 3–7, we can see that both the cost function and gradients drop dramatically in the optimization process, which demonstrate the strong inversion ability of this model.Figure 3
The valuesJ
/
J
0 versus iteration number, where J
0 represents the initial value of the cost function.Figure 4
The curveF
i represents the gradients of J with respect to α
i
0 versus iteration number, i
=
1
,
…
,
6, respectively.Figure 5
The curveP
i represents the gradients of J with respect to β
i
0 versus iteration number, i
=
1
,
…
,
6, respectively.Figure 6
The curveZ
i represents the gradients of J with respect to γ
i
0 versus iteration number, i
=
1
,
…
,
6, respectively.Figure 7
The curveG
(
i
) represents the gradients of J with respect to ith type of the OBCs versus iteration number, i
=
1
,
…
,
7, respectively.


## 6. Summary and Conclusions

One key of 4D VAR data assimilation is optimization of the model parameters such as the initial values, OBCs, BFC, interface friction coefficient, and horizontal eddy viscosity coefficient [13]. It has been shown that the POD method by combining with some numerical methods for partial differential equations can provide an efficient means of generating reduced order formulations and alleviating the computational load and memory requirements [20]. As a preliminary feasibility study, we combine the POD method and 4D VAR data assimilation with a 2D tidal model to simulate the M
2 tide in the BYECS. The numerical experiment is carried out to optimize the initial values and the temporally and spatially varying open boundary conditions by assimilating T
/
P altimeter data into the 2D tidal model. This model consists of two submodels: the POD reduced order forward model is used to simulate the tides, while its adjoint model is used to optimize the control variables. According to the division of the BYECS, we classify the OBCs as 7 types.We can see that, after 200 iterations of assimilation, both the cost function and gradients with respect to the initial values and OBCs have been reduced dramatically. The absolute mean difference between the model results and observations is 3.2 cm and 2.9∘ for amplitude and phase-lag, respectively, better than the results of Lu and Zhang (2006), suggesting that the construction of the POD reduced order model and the inversion of control variables are successful. In future research we will extend the application of our model to more cases such as inverting other model parameters. The relevant work is underway.



---

*Source: 6751537-2016-01-20.xml*