Optimization of Open Boundary Conditions in a 3 D Internal Tidal Model with the Adjoint Method around Hawaii

and Applied Analysis 3 In addition, the initial conditions are that the initial surface elevation, interfacial elevation, xand y-velocity at all layers are zeros. 2.2. Adjoint Model. The adjoint model is similar to that in Chen et al. [27], which can be used to optimize OBCs and/or other spatially varying parameters and displayed as follows. The cost function is defined by J (q, u, V; p) = 1

The work of this paper is an extension of Chen et al. [27], in which Fourier coefficients of  2 internal tide on four open boundaries are inverted simultaneously.During the optimization, steepest descent method with special step factor is used to minimize cost function.This paper is organized as follows.The 3D internal tidal model with adjoint method, minimization algorithm, independent point (IP) strategy and model settings are introduced in Section 2. Identical twin experiments (TEs) are carried out to test the model in Section 3. Practical experiment (PE) and its analysis are given in Section 4. Finally, the summary and conclusions are provided in Section 5.

Numerical Model and Settings
2.1.Forward Model.The governing equations of internal tide in spherical coordinate are described as follows: in which, where  = in which  is the nondimensional bottom friction coefficient), and   and   ( = 1, . . ., −1) are interfacial friction between the kth layer and the ( + 1)th layer, and they are defined as where  V is the vertical viscosity coefficient; ℎ  is the undisturbed thickness of the kth layer.
Equations ( 1) and its vertical integration, which are called the internal and external mode, respectively, compose the forward model and are used to simulate the internal tide.Miao et al. [31] used this model to investigate the characteristic of generation and propagation of  2 and  1 internal tides in the South China Sea.
The Flather condition [38] can yield good results [27,31] when applied to the open boundaries in this model.In this paper, an adaptation of the Flather condition is installed in the external mode of the forward model, which is given as follows: where  and  are Fourier coefficients of OBCs,  is the frequency of  2 tide, and Δ is time step.So the inverse problem of OBCs is turned into the inverse problem of Fourier coefficients  and .
For the internal mode, the relaxation conditions are applied.A relaxation region of several grids close to the boundary is defined.In this region, tidal values including the interface elevation and fluid velocity of both barotropic and baroclinic tides are first computed and then relaxed to final values displayed as follow  (, ) =   (, ) + [1 −  (, )] ⋅   (, ) , (7) where  represents the tidal state   ,   , or V  ,   and   represent the barotropic and baroclinic tides, respectively, and (, ) ∈ [0, 1] is the relaxation factor in grid (, ).In our model,  is increased from 0 to 1 as the grid is getting close to the open boundary.
In addition, the initial conditions are that the initial surface elevation, interfacial elevation, x-and y-velocity at all layers are zeros.

Adjoint Model.
The adjoint model is similar to that in Chen et al. [27], which can be used to optimize OBCs and/or other spatially varying parameters and displayed as follows.
The cost function is defined by where Σ is the model domain of both time and space, ζ is the observations, and p represents the control variables and is defined by the OBCs, which are denoted by a series of Fourier coefficients in this work.
The Lagrangian function is defined by where   ,   and V  are adjoint variables of   ,   , and V  , respectively.
According to the typical theory of Lagrangian multiplier method, we have the following first-order derivates of Lagrangian function with respect to all the variables and parameters: Equations ( 10) are ( 1) actually.The adjoint equations can be derived from (11).From (12), the gradients of the cost function with respect to control variables can be obtained and displayed as follows: Details on computing   , are given in [27].

Minimization Algorithm.
Considering that the aim is to make the simulated results closer to the observations, there is the optimization problem as follows: where  is the cost function;  represents Fourier coefficients  and .Suppose that the dimension of  is  and the gradient of  is ∇(   ) with respect to the ith parameter (1 ≤  ≤ ) and the kth iteration step.According to the typical steepest descent method, where    is the negative of ∇(   );   is the step factor at the kth iteration step.The initial parameter  0  should be prescribed (in the following experiments,  0  = 0).In order to minimize , the step factor   used in this paper is performed as follows: where  0 is a constant which can be obtained by prior experience; ‖∇( 1  )‖ 2 is the  2 norm of gradient at the first iteration step.Then formula ( 6) is changed into  where ∇(   )/‖∇(   )‖ 2 is the normalized gradient at the kth iteration step and ‖∇(   )‖ 2 /‖∇( 1  )‖ 2 is the ratio between  2 norm of gradient at the kth iteration step and that at the first iteration step, which could be a decreasing value as iteration step is increased.

Independent Point Strategy.
Independent point strategy is considered to be a useful tool in the estimation of OBCs and spatially varying parameters.On one hand, adopting IP strategy could decrease the the number of variables so as to increase the well posedness of inverse problem.It is generally accepted that number of observations should be larger than that of variables and if number of observations is almost equal to or smaller than that of variables, the inverse problem will become ill-posed and the inverse results may be bad.On the other hand, adopting IP strategy could ensure the continuity of variables, which accords with the physical properties of variables.Lu and Zhang [28] studied the inversion problem of spatially varying bottom frication coefficients by using adjoint method and IP stragety.Results indicated that adopting IP strategy in which IPs are distributed according to depth could lead to a better simulation.Based on the theory of inverse problem, Zhang and Lu [24] simulated 3D tidal current in the Bohai and Yellow Sea by assimilating the satellite altimetry.They found that too large number of variables would lead to inaccurate inversion results, while too small number of IPs could not figure out the real condition of ocean.Guo et al. [25] and Cao et al. [26]   eddy viscosity  ℎ , and vertical eddy viscosity  V are 0.003, 1000 m 2 /s, and 0.005 m 2 /s, respectively.

Twin Experiments
TEs are designed to test the reasonability and feasibility of the model.The results and conclusions of TEs are useful for PE.The process of TE is designed as follow.(a) Run the forward model with prescribed OBCs.The simulated sea surface elevation recorded at the grid points on T/P satellite tracks is taken as the "observations." Here only the positions of T/P tracks are used and the "observations" are not the real observations but the simulated results of forward model based on the prescribed OBCs.(b) Initial values of Fourier coefficients (taken as 0 here) are assigned to the OBCs to run the forward model.The values of state variables are obtained.(c) The difference of sea surface elevation between simulated results and "observations" serves as the external force of the adjoint model.The adjoint variables are obtained through the backward integration of the adjoint equations.(d) After obtaining the values of state variables and adjoint variables, the gradient of Fourier coefficients can be computed and the Fourier coefficients can be updated.Repeat (b), (c), and (d); the Fourier coefficients will be optimized and the difference of sea surface elevation between simulated results and "observations" will decrease.Once some convergence criterion is met, the process of inversion terminates.
The iteration will be terminated once a convergence 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 [32].In our work, the criterion of TEs is that the number of iteration steps is equal to 30 exactly.

Inversion of Fourier Coefficients on Four Open Boundaries.
In this section, we study the inverse problem of Fourier coefficients on four open boundaries simultaneously, which means more variables are estimated compared with the work of Chen et al. [27] where only Fourier coefficient a on west open boundary is inverted.In addition, the feasibility and efficiency of IP strategy are tested in this section.Figure 3 is the snapshot of prescribed OBCs and the inverse values.In TE1a, IP strategy in which IPs are distributed evenly every four grid points on all four open boundaries is employed, while in TE1b no IP strategy is used.
As a whole, the shapes of inverted OBCs are similar to those of prescribed ones in both TE1a and TE1b.But it cannot be denied that inverted OBCs are disordered with fluctuations in TE1b where no IP strategy is used, implying the ill-posedness of inverse problem in TE1b and revealing feasibility and efficiency of IP strategy.Figure 4 shows the cost function and gradient norm which are normalized by their initial values at the first iteration step.As iteration step increased, the values of cost function and gradient norm are decreased obviously and at the end of iteration could reach 1% and 10% of their initial values, respectively, implying the optimizing ability of the model.In addition, as we use the step factor mentioned in Section 2, the cost function decreases efficiently without fluctuations.Tables 1 and 2 report the mean absolute errors (MAEs) and relative coefficients (RCs) between inverted OBCs and prescribed ones in TE1a and TE1b, respectively.All MAEs and RCs are less than 3 cm and greater than 0.7, respectively, for both Fourier coefficients a and b and on every  open boundary, suggesting the reasonability and feasibility of the model.In addition, the inverse results in TE1a are better than those in TE1b, revealing the feasibility and efficiency of IP strategy.Combining all these results, it can be concluded that the 3D internal tidal model involving adjoint method is suitable to invert OBCs and the IP strategy is useful in the estimation of OBCs.

Inversion by Assimilating "Inaccurate Observations".
Considering that practice in situ observations may contain noises, artificial random errors are added to "observations." In this section, based on the "observations" in TE1a, five TEs are carried out and the maximum percentages of error are 10%, 20%, 30%, 40%, and 50% in TE2a, TE2b, Te2c, TE2d and TE2e, respectively.Then these "observations" are assimilated to invert OBCs.Results are displayed in Figures 5, 6, and 7.
Figure 5 illustrates the inverse results.All the shapes of inverse OBCs are consistent with those of the prescribed ones.It is found that with increased artificial errors, the inverted OBCs become further away from the prescribed values, implying that larger errors will result in worse inversion of OBCs.The cost function and its gradient norm are plotted in Figure 6.After assimilation, the cost function is about 1% of its initial value in TE1a and about 3% in TE2e.Combined with MAEs and RCs between inverted OBCs and prescribed ones in Figure 7, where all MAEs are less than 3.5 cm and RCs are larger than 0.6, it is reasonable to believe that the prescribed OBCs have been inverted successfully though assimilating the "observations" is not accurate.But it cannot be denied that as artificial errors increased, the inversion results become worse with larger MAEs and smaller RCs.In addition, in these TEs, the inversion results of a are better than those of b.

Practical Experiments
The 3D internal tidal model is used to simulate the  2 internal tide around Hawaii by assimilating T/P data with OBCs optimized.IP strategy is used, which is the same as that used in TE1a.Both locations and values of observations are used in PE, which is different from TEs.In PE, the iterative assimilation terminates once the number of iteration steps reaches 50 exactly.
After assimilation, the cost function is about 3% of its initial value.The simulated surface cochart is illustrated in Figure 8. Considering the co-phase lines are disordered due to the influence of internal tide, only coamplitude lines are contoured.On the whole, the simulated cochart shows good agreement with that obtained from Oregon State University (OSU) tidal model in which only barotropic tide is simulated.The sea surface amplitude is about 25 cm at the southeast of Hawaii and smaller than 15 cm at the northwest.Careful inspection shows that the major difference between simulated results and OSU model results is in the short-scale fluctuations superposed on co-amplitude lines.These fluctuations can be considered as the surface manifestation modulated by the internal tide, from which the propagation of internal tide could be found.The  2 internal tide generated along the Hawaiian Ridge mainly propagates towards to northward and southeastward.
Figure 9 illustrates the difference of surface amplitude along tracks 006, 045, 082, 147, 184, and 223.Simulated results show good agreement with observations.Table 2 shows the corresponding MAEs.All MAEs of amplitude and phase are less than 2 cm and 8 ∘ (Table 3), respectively.Relative errors of amplitude are about 0.1, smaller than those in Ray and Mitchum [3], suggesting the advantage of this 3D internal tidal model involving adjoint method in simulating internal tides.
Compared with the barotropic results from OSU tidal model, we can calculate the surface manifestation modulated by the internal tide, the maximum of which are 4.0, 4.2, 5.6, 3.6, 9.1, and 6.7 cm, corresponding to tracks 006, 045, 082, 147, 184, and 223, respectively.Except the manifestation of 9.1 cm

Conclusion and Summary
In this paper, the optimization of OBCs in a 3D internal tidal model is investigated with the adjoint method.Fourier coefficients of  2 internal tide on four open boundaries are inverted simultaneously, which means more variables are estimated in the inverse problem.Special step factor is used in the steepest descent method in order to minimize cost function efficiently without fluctuations.The reasonability and feasibility of the model and IP strategy are tested by TE1.In   TE2, artificial random errors, of which the percentage is from 10% to 50%, are added to observations.After assimilation, the cost function reaches less than 3% of its initial values MAEs and RCs between inverted OBCs and prescribed ones are less than 3.5 cm and larger than 0.6, respectively, indicating that the prescribed OBCs are inverted successfully.But it cannot be denied that as artificial errors increased, the inversion results become worse gradually.
In PE, the OBCs are optimized and the  2 internal tide around Hawaii is simulated by assimilating T/P data.The simulated co-chart shows good agreement with that from Oregon State University (OSU) tidal model and T/P observations on the whole.Careful inspection shows that the major difference between simulated results and OSU model results is in the short-scale fluctuations superposed on co-amplitude lines, which can be treated as the surface manifestation modulated by the internal tide.The maximum manifestation along tracks 006, 045, 082, 147, 184, and 223 are 4.0, 4.2, 5.6, 3.6, 9.1, and 6.7 cm, receptively.Except the manifestation of 9.1 cm along track 184, the values along other tracks are consistent with the estimation of Ray and Mitchum [3,4].

Figure 1 :
Figure 1: Bathymetry and T/P tracks around Hawaii ("o" represents grid point on the open boundary).

Figure 5 :
Figure 5: (a) Prescribed OBCs and inverted ones on the west open boundary in TE2.(b) Prescribed OBCs and inverted ones on the south open boundary in TE2.(c) Prescribed OBCs and inverted ones on the east open boundary in TE2.(d) Prescribed OBCs and inverted ones on the north open boundary in TE2.

Figure 6 :
Figure 6: Cost function and gradient norm in TE2.

Figure 7 :
Figure 7: MAEs and RCs in TE2 (aw, bw, as, bs, ae, be, an, and bn represent Fourier coefficients on west, south, east, and north open boundaries, respectively).

Table 2 :
RCs between inverted OBCs and prescribed ones in TE1.

Table 3 :
MAEs of amplitude and phase between observations and simulated results (values in brackets are mean relative errors of amplitude).