Inverse Estimation of Open Boundary Conditions in the Bohai Sea

This paper presents an algorithm for the estimation of open boundary conditions OBCs which force tides in the interior region by an adjoint data assimilation approach. Assuming that OBCs are position dependent, OBCs can be approximated by linear interpolation among values at certain independent points IPs . Twin experiments are performed to examine the sensitivity of the model to the IP distribution and interpolation radius. It is proved that the prescribed OBCs can be well recovered with appropriate number of IP and interpolation radius. In the Bohai Sea model domain with horizontal resolution of 10′×10′, the appropriate number of IP is 3 and the interpolation radius is 60′. In the practical experiment, the M2 constituent in the Bohai Sea is simulated by assimilating the T/P data and tidal gauge data. The mean absolute errors in amplitude and phase are 5.0 cm and 5.7◦, respectively, and the cochart obtained shows the character of M2 constituent in the Bohai Sea.


Introduction
Open boundary conditions OBCs are crucial for the representation of tidal processes in the regional ocean model.Traditionally, they can be obtained from global tidal estimates, larger scale models, or extrapolation from available observations.However, experience is required when the methods mentioned above are used.In addition, OBCs obtained need to be adjusted by comparison between model results and observations, which is time consuming.
The adjoint method is a powerful tool for parameter estimation.Inverse estimation of OBCs in meteorology and oceanography dates back to Sasaki et

Optimization of Independent Points in the OBCs and Interpolation Scheme
We use the same adjoint tidal model as in Lu and Zhang 13 .In this study, the amplitude and phase at IPs are transformed into Fourier coefficients.The correction of the independent Fourier coefficients is as follows: where w j,jj is the weight coefficient in the Cressman form 14 where R is the interpolation radius and r j,jj is the distance from grid point j to jj on the boundary.

Twin Experiments
In this study, twin experiments are implemented in the Bohai Sea Figure 1 from 117.5 • E to 122.5 • E and from 37 • N to 41 • N. The horizontal resolution of this model is 10 × 10 .The open boundary is installed on the right edge.The bottom friction coefficients are taken as constant 1.2 × 10 −3 .The iteration process of twin experiment TE is designed as follow: 1 Run the original dynamic forward model with prescribed OBCs.Take the simulated results at the grid points on T/P satellite tracks and tidal gauges as the "observations." 2 A value is assigned to each IP the initial guess is taken as zero here and values of other grid points on the open boundary are obtained by interpolation.The forward simulation is performed under this condition, which gives the values of state variables.3 The differences between simulated results and observations serve as the external force of the adjoint model.And the adjoint variables in a period of M 2 constituent are calculated through backward integration of the adjoint equations.4 Update OBCs by correction equation and interpolation.Repeat steps 2 -4 until some stopping criterion is satisfied.As the OBCs are increasingly optimized, the differences between simulated results and observations will decrease.Meanwhile, differences between the prescribed OBCs and inverted ones are usually decreased as well.The iteration will be terminated once a stopping criterion is met.The criterion could be that the last two values of the cost function are sufficiently close, the magnitude of the gradient is sufficiently small, the discrepancy between the updated and old parameters is sufficiently small, or a combination of these.In this work, the criterion is that the number of iteration steps is equal to 100 exactly in both twin experiments and practical experiment.

Influence of the Interpolation Radius on Inversion Results
The model is sensitive to interpolation radius when the algorithm presented in this paper is implemented.In this part, TE1∼TE4 are conducted with prescribed OBCs Figure 2 , in which 3 IPs are uniformly distributed Figure 4 a and the interpolation radius is varying, that is, R 60 , R 80 , R 100 , and R 120 .The results are shown as follow.
The decreasing process of cost function value is an important criterion to access the inversion.The ratio of the cost function and its initial value is used to describe this process.The inverted OBCs, cost function, and gradient are shown in Figures 2 and 3, respectively.From Figure 3, we can find that after 100 iteration steps, the cost function and gradient can reach 1% or even 0.1% and 10 −3 ∼10 −5 of their initial value.From Figure 2, it can be concluded that by the use of IP scenarios, the prescribed OBCs can be inverted successfully by assimilating the "observations."In addition, we can find that when R is larger than the maximum distance between adjacent IPs the cases of R 80 , R 100 and R 120 , the inversion results may have strange values between adjacent IPs which is related to the interpolation in the Cressman form.
In this paper, mean absolute error MAE is used to test the simulated results.MAE is defined as follow: where A is amplitude or phase, ΔA is the MAE of A, n is the total number of observations, and A obs and A simu represent the observation and simulated results at observation locations, respectively.Considering the MAEs shown in Table 1, we can conclude that better inversion results can be obtained when the interpolation radius is equal to the maximum distance between adjacent IPs.Furthermore, the difference increases with the growth of interpolation radius.

Influence of IP Distribution on Inversion Results
The number and distribution of IPs also have an important influence on the inverse estimation.To explore the influence, three strategies are designed, that is, 3 IPs in strategy A , 5 IPs in strategy B , and 7 IPs in strategy C , and applied to TE5∼TE7 where R is equal to the maximum distance between adjacent IPs, respectively.Figure 5 displays the prescribed OBCs and inversion results of TE5∼TE7.
Cost function and gradient are shown in Figure 6 from which we can find that after 100 iteration steps, the cost function and gradient can reach 1% or even 0.1% and 10 −3 ∼ 10 −5 of  their initial value, respectively.Considering the inverted OBCs shown in Figure 5 and MAEs between prescribed OBCs and inverted ones shown in Table 2, we can conclude that the OBCs can be inverted successfully by using IP scenarios, which demonstrates the reasonability and feasibility of the model.Though the prescribed OBCs can be inverted by using all these 3 strategies, strategy A with 3 IPs and R equal to 60 leads to better results smaller values of MAE, cost function, and gradient .

Practical Experiments
According to the conclusion of Sections 3.1 and 3.2, in the practical experiment, only strategy A is used to estimate the OBCs and simulate the M 2 constituent in the Bohai Sea.Tidal gauge data are used during the assimilation, while T/P data are used to test the simulated results.The bottom friction coefficients are similar to those in Lu and Zhang 13 .The MAEs of amplitude and phase are 5.0 cm and 5.7 • , respectively.Cochart based on the model results shows the characteristic of M 2 constituent in the Bohai Sea Figure 7 .There are two amphidromic points in the Bohai Sea, of which one is near Qinhuangdao and the other is near the Yellow River delta.

Summary
OBCs have crucial impacts on the representation of tidal processes in regional ocean models.They can be inverted by the adjoint method.In this study, OBCs are assumed to be positiondependent and therefore can be approximated by linear interpolation among certain nodal values.In this way, we reduce the number of control variables.
Twin experiments are performed to verify the feasibility of the method and examine the sensitivity of the model to the number of IPs and interpolation radius.It is proved that the prescribed OBCs can be well recovered with appropriate number of independent points and interpolation radius.In the Bohai Sea model domain with horizontal resolution of 10 × 10 , the appropriate number of IPs is 3 and the interpolation radius is 60 .
Based on the strategy of appropriate number of IPs and interpolation radius obtained in the twin experiment, the M 2 constituent in the Bohai Sea is simulated by assimilating the T/P data and tidal gauge data.The MAEs of amplitude and phase are 5.0 cm and 5.7  respectively, and the cochart obtained can describe the character of M 2 constituent in the Bohai Sea.Used to more general and open domain, the IP strategy is still feasible and reasonable; however, the results will vary with specific conditions while two conclusions stay unchanged.On one hand, under the situation that IPs are distributed evenly, better inversion results can be obtained with interpolation radius equal to the distance between adjacent IPs.On the other hand, the IP number should be chosen carefully, avoiding a too large one.So when this method is used to a more general and open domain, only several experiments are needed in order to obtain the best simulated results.
al. 1 .Shulman et al. 2 , Nguyen 3 , Das and Lardner 4, 5 , Lardner et al. 6 , Seiler 7 , and Zhang et al. 8, 9 made some studies on optimizing parameters or estimating OBCs in different models over the past few decades.Efforts also have been made with regard to inverse estimation of OBCs in the China Seas.He et al. 10 investigated the shallow water tidal constituents in the Bohai and Yellow Sea by assimilating T/P altimeter data with the adjoint method.Han et al. 11 employed a 2-dimensional nonlinear numerical POM to describe tide in the China Seas, in which the OBCs and other unknown internal model parameters were inversed in the meanwhile.With the help of a 3D numerical barotropic adjoint tidal model, Zhang and Lu 12 studied the estimation of OBCs more deeply and simulated the 3D M 2 tides and tidal currents in the Bohai and Yellow Seas by assimilating the satellite altimetry data.In the above-mentioned studies about inverse estimation of OBCs, all grid points on the open boundary are inverted independently.In order to reduce the number of control variables, amplitude and phase of each tidal constituent at a given open boundary grid were taken as a quadratic polynomial based on the results of Schwiderski's global tidal model, which decreased the control variables to six for each tidal constituent 9 .Similar to what had been done in Das and Lardner 4, 5 , we assume several independent points IPs on the open boundary.Values at the IPs can be optimized by correction equation, and those at other grid points on the open boundary are determined by linearly interpolating the values at IPs.

Figure 1 :
Figure 1: The bathymetry map of the Bohai Sea and positions of T/P altimeter tracks "•" and tidal gauges "o" .

Figure 2 :Figure 3 :
Figure 2: Prescribed OBCs and inverted ones corresponding to different interpolation radius.

Figure 4 :
Figure 4: IP distribution of the 3 strategies large and small dots denote IPs and the other grid points on the open boundary, resp .

Figure 5 :
Figure 5: Prescribed OBCs and inverted ones corresponding to different strategies.

Figure 7 :
Figure 7: Cotidal chart obtained in the practical experiment using strategy A the dashed line denotes coamplitude line meter and the solid line denotes cophase line degree .
1where jj is the index of IPs and j is the index of grids on the entire open boundary.a jj and b jj are optimized values of Fourier coefficients at IPs, while a jj and b jj are prior values.W j,jj is the weight of linear interpolation.K is the undetermined constant.T i varies with specific conditions, that is,

Table 1 :
Differences between prescribed OBCs and inverted ones.

Table 2 :
Differences between prescribed OBCs and inverted ones.