Analysis of the Properties of Adjoint Equations and Accuracy Verification of Adjoint Model Based on FVM

and Applied Analysis 3


Introduction
Numerical simulation becomes an effective tool for the analysis and forecasting in hydrodynamics, hydrology, oceanic and atmospheric.However, the uncertainties in numerical model, for example, initial conditions, open boundary conditions, and physical parameters (wind drag coefficients, bottom friction coefficient, and so on), make the underlying problem of calibration very challenging.
Adjoint data assimilation method (shortened as "adjoint method" hereafter), based on optimal control theory and 4dimensional variational theory, is one of the most effective data assimilation approaches with extensive applications in the studies of meteorological and oceanic numerical simulation [1,2].It has the advantage of assimilating observations that are distributed in time and space in the numerical model, while maintaining dynamical and physical consistency with the model.The adjoint method is primarily focused on reducing model errors induced by uncertainties in the model [3][4][5][6].
The derivation of adjoint numerical model (ANM) is one of the essential issues for developing an adjoint data assimilation system.There are two different approaches on how to get ANM [7].The first approach is "adjoint of model, " that is, getting numerical solution of the forward continuous equations and then deriving ANM from this numerical model.The second approach is "adjoint of equation, " that is, deriving the adjoint equations directly from the forward continuous equations and then determining the adjoint numerical model by modeling the adjoint equations.
It is conventionally believed that the first method should be adopted, so most studies employed the "adjoint of model" to construct ANM [1,7,8].It has to be stressed out that using the approach of "adjoint of model" to get ANM from forward numerical model allows no selection of appropriate numerical methods according to specific needs, such as the problem of discontinuous solution.In fact, some studies suggest that "adjoint of equation" will have similar effectiveness and can also get appropriate descend gradient [9,10]."Adjoint of equation" avoids the limitations of forward numerical model and is capable of designing numerical scheme according to mathematical properties of adjoint equation.
Most adjoint numerical models are based on finite difference method (FDM).FDM has good computational performance and coding efficiency.However, the classic FDM based on rectangular meshes could not approximate boundary accurately in bays and estuaries, which have complex geometries [11,12].To improve the fitting ability of numerical model on the boundary in bays and estuaries, finite element method (FEM) based on adaptive meshes was introduced into adjoint numerical model [13,14].However, in FEM the established algebraic equations are nonlinear and have to be solved by iteration process [15].Moreover, computation of the adjoint method is huge itself: the forward model and adjoint model would be integrated forward and integrated backward once, respectively, in each assimilation iteration.So, FEM will certainly increase the computational cost of the assimilation greatly.
In fact, the finite volume method (FVM), which has the merits of FDM and FEM, combines the best attributes of FDM for simple discrete coding and computational efficiency and FEM for geometric flexibility.Furthermore, the integral equations can be solved numerically by flux calculation; the FVM is better suited to guarantee mass conservation [16].
To develop adjoint numerical model based on FVM, the properties of equation will be discussed first.Although some researches used "adjoint of equation" to develop adjoint model [9,10,13], few studies deal with the properties of adjoint equation, such as hyperbolicity and homogeneity of the adjoint equation.In addition, What are the reasons that FVM can be used to solve adjoint equations?If the new adjoint model based on FVM viable?In the study we will discuss these issues.
This paper is organized as follows.In Section 2, the adjoint equation of conservative shallow water equation of Cartesian coordinate was derived, the formula for inversing bottom friction coefficient is given, and the differences between nonlinear shallow water equation and its adjoint equation are analyzed.In Section 3, properties of adjoint equation, such as characteristic structure, hyperbolicity, and homogeneity of flux, are discussed.In Section 4, the forward model and its adjoint model based on finite volume methods and unstructured triangular meshes are described in brief.In Section 4, the correctness of the new adjoint model and code is verified via the gradient check.Then, in Section 5, the effects of the adjoint model by assimilating the two different types of observation are discussed by a series of twin experiments.Concluding remarks are included in Section 6.

Water Shallow Equations and Their Adjoint Equations
where  1 ,  2 ,  3 are weight coefficients; , , V are numerical values; and  obs ,  obs , V obs are observed values in the assimilate window where there are observed data; otherwise, they are equal to numerical value.By means of Lagrange multiplier method, , ,  are introduced into (1), respectively, as Lagrange multiplier.The objective functional is rewritten as There must be ∇   = 0 while the function has extreme value.Utilize (, , V) = 0 on the boundary and initial fields of forward process |  0 = |  0 = V|  0 = 0, and let |  = |  = |  = 0 be the initial field.Thus, the adjoint equation is obtained as And the formula for inversing bottom friction coefficient is where and  is the times of assimilation.
The most significant difference between shallow water equation and its adjoint equation is that shallow water equations are forward integration, while adjoint equations are backward integration, so shallow water equations are often called forward equations and adjoint equations are often called backward equations.In addition, shallow water equations are forced by gravity term, bottom friction term, and other external forces; the variables and terms in shallow water equations have physical meanings.However adjoint equations are derived from forward equations and driven by the difference between observed values and numerical values..The adjoint equations can be rewritten as

Hyperbolicity of Adjoint Equation.
Let  * matrix be the linear combination of Jacobian matrices  * and  * , ],  1 ,  2 are real parameters, and √ 2 1 +  2 2 > 0. Thus, Let  = √, and the eigenvalues of  * are According to the definition of hyperbolic equation, we have the following.
(1) For  ̸ = 0, all eigenvalues of  * are real numbers; that is, there are  real eigenvalues.Thus, the adjoint equations are hyperbolic equations.
(2) When  > 0, that is, when  ̸ = 0, the eigenvalues of  * do not equal each other.Then the adjoint equations are strict hyperbolic equations.

Homogeneity of Flux.
According to the expressions of  * ,  * ,  * ,  * , and  * , we have This means that  * ( * ),  * ( * ) in adjoint equations meet the requirements of homogeneity property, so that the adjoint equation can be used to construct flux vector splitting directly.

Finite Volume Scheme of Adjoint Equation
In the section, the FVM is used to solve the forward equations and their adjoint equations, and the computational control volumes are based on the same unstructured triangular meshes.
As regards the forward equations, a lot of discontinuity capture methods are used to solve the shallow water equations, including the cases with strong discontinuity, such as tidal bore [15].In this study, the numerical flux is evaluated by Roe's scheme [17,18] and the cell variables are reconstructed by the piecewise linear model [19].For the sake of brevity, we do not describe them again.
It has been proved that adjoint equations are hyperbolic equations in Section 3.2.A significant feature of hyperbolic equations is discontinuity of their solution.Therefore, we used the discontinuity capture methods which are used to solve the shallow water equations to solve adjoint equations.
Using FVM based on unstructured meshes, computational domain is performed on triangular element and the variable is set at the center of each element.Let  ⇀  * = ( * ,  * ) be the numerical flux at interfaces of control volumes; the finite volume scheme of adjoint equations can be formulated as follows: where Ω  is the area of  where || =  |Λ| .
According to the derivation in Section 3.   Using the Green formula, the gradient terms in the adjoint equation source terms ( * ) can be approximately computed as where , ,  are the centers of the neighbor cells, respectively (Figure 1), and   is the area of triangle △.For /, V/, V/, /, /, they can be approximately computed similarly to /.

Verification of the Adjoint Model Based on FVM
In this section, we perform data assimilation experiments to verify the correctness and to discuss the calculation efficiency of the adjoint model in tidal bore estuary.The domain is a generalized estuary, the length is 20 km, the width at upstream is 2 km, the width of estuary is 10 km, and the water depth is 10 meters.The flow is calculated on a mesh of unstructured triangle total of 588 nodes and 1064 cells.The minimum side of the grid is about 100 meters (Figure 2).In the numerical experiments, the gravitation constant is taken as 9.8 m/s; the initial velocity and the free surface elevation above the still water level are set to be zero as  The forward model is integrated five period.When flow propagates to the straight channel, the tidal bore is formed.As discussed in [20], the tidal height reaches 3 meters within a short period (tidal bore) 9 kilometers away from the estuary.
For adjoint data assimilation, descent gradient of cost function is obtained by the calculation of ANM.Thus, the correctness of the adjoint model and code should be verified, which can be performed via a gradient check [14,21]; that is, the quantity is considered and it is verified that it is of order 1+().() is the cost function at a general point  (in the study, it is bottom friction coefficient), ∇() is the computed gradient,  is a disturbance quantity, and ℎ is an arbitrary unit vector (such as ℎ = ∇‖∇‖ −1 ).Table 1 shows typical values of Φ() and Log 10 |Φ() − 1| for a sequence of decreasing values of .Clearly, Φ() is 1 + ().Figure 3 is the curve of Φ() and Log 10 |Φ() − 1| varying with disturbance .It can be seen that Φ() has 1order accuracy and the curve of Log 10 |Φ() − 1| varying with  also shows an obvious V-shape, also as expected in literature [14].Thus, it can thus be concluded, at least for this model problem setup, that the new ANM based on FVM and Roe's scheme is correct and feasible.

Bore Estuary
In this section, In order to discuss and display the performance of assimilation of the adjoint model in generalized tidal bore estuary, we perform "twin experiments" to invert Manning's roughness coefficient by assimilating the observation of water level or velocity; the observed values are four kinds: water level or velocity at all the cells, at cells with odd numbers (half of the grids), at cells with the number ending with 3 (1/10 of the grids), and at the four cells (3 km, 0), (8 km, 0), (13 km, 0), and (18 km, 0).The domain and experiment setup are the same as in Section 5. Process for the twin experiments is as follows.
(1) Given true bottom friction coefficient (0.02), calculate forward model 5 periods and then save water level and current velocity of last periods as "observed value." (2) Given the initial or calibrated bottom friction coefficient, calculate forward model 5 periodic and calculate cost function based on the calculated value of the last periodic and "observed value." If the cost function satisfies the specified termination conditions, then the assimilation process ends; otherwise, proceed to the next procedure.
(3) Calculate the ANM and calibrate bottom friction coefficient.Then go back to procedure (2) to continue assimilation.The initial values and inverted results for the four cases are given in Table 2.It can be seen from Table 2 that the assimilation algorithm converges to the true value (0.02) when initial values are given as 0.01 in the four tests.It also can be found that when more observations were assimilated, the assimilation algorithm conversed faster.The closest result to the true value with minimum assimilation times ( 21) is obtained, when assimilating water level at all cells.Satisfactory results are also obtained when assimilating water level data at only four cells, though the assimilation times (51) are greater.When the initial values of 0.02, 0.05, and 0.15 are given, similar results are also obtained.
It should be pointed out that the adjoint assimilation algorithm is quite sensitive to initial values when assimilating water level data at only four cells.In this case, the initial values should near to the true value, or the algorithm does not converge.

Assimilate Different Quantities of Current Velocity
Observed Data.Assimilating the observation of current velocity, if we let the initial guesses drag coefficient equal to the value which used in assimilating the observation of water level, the assimilation process cannot converse.When we let the initial guesses drag coefficient near to the true value, it works.By assimilating the observation of current velocity at all cells, we let the initial guess coefficient be 0.015; after 29 assimilations, we get the inverted result to be 0.019977.At onehalf of all cells, the initial guess is 0.015; after 37 assimilations, the inverted result is 0.00045.At one-tenth of all cells, the initial guess is 0.0015; after 52 assimilations, the inverted result is 0.019964.At four cells, we let it be 0.017; after 73 assimilations, we get the inverted result to be 0.019959.The initial values and inverted results for the four cases are given in Table 3. Figure 4 shows the first fifteen iterations of the cost function with assimilating water level observed value at all cells (real line) and with assimilating current velocity observed value at all cells (dashed), respectively.It can be seen that the performance with assimilating the two types of observed value behaves differently.When assimilating current velocity, the cost function decline fast, but its convergence speed is slow; When water level was assimilated, the cost function convergence to minimum value smoothly.The phenomena may result in the fact that the observation of current velocity is more sensitive than the water level for the coefficient.So when current velocity was assimilated, the cost function responded more quickly; also because of this, it is hard to converge.The results suggest that the two types of observation should be assimilated by a suitable balance.

Conclusions and Outlook
This study gives a detailed analysis and summary of the properties of adjoint equation and analyzes the characteristic structure of adjoint equation, hyperbolicity of equation, equation homogeneity, and the splitting characteristic of Jacobian matrix.It has been proven that the adjoint equation belongs to hyperbolic equation and demonstrated that the attention should be paid to the solution discontinuity of hyperbolic equations when constructing adjoint numerical model.We also demonstrate the homogeneity of the adjoint equation and demonstrate that ANM in flux vector splitting can be constructed directly.Based on the theoretical analysis, we verify the accuracy and assimilation effect of ANM constructed by finite volume scheme.Regarding the disputes arising from two different construction methods of ANMs, we believe that it is associated with the design and selection of numerical algorithm, which significantly influence the accuracy and validity of ANM.Thus, theoretical analysis on adjoint equation should be carried out before developing ANMs.However, this study is only a preliminary research on adjoint equation theory.Further researches on issues related to adjoint equation are required, including the basic properties of adjoint equation solutions, relation of adjoint equation, and certain mathematical and physical equations as well as physical meaning of adjoint equation.The clarification of these problems may provide necessary guidance for the construction of reasonable adjoint numerical model.
The "twin experiments" of generalized tidal bore estuary to invert Manning's roughness coefficient show that the spatial distribution of the bottom friction coefficient will influence the precision of the inversion; plenty of observation data can improve the accuracy of the inversion.The results also indicate that the two types of observation should be assimilated in a suitable balance th control volume, Δ  is the length of the  th side of the  th control volume, and  *  ,  *  is the mean value of  *  ,  *  on the control volume Ω  .The adjoint variables at the two sides of control volume boundary are considered as different values ( *  ,  *  ).The variables are set at the center of grid element, and  *  ,  *  are reconstructed with piecewise constant approximation.The numerical flux of adjoint equations based on Roe's scheme can be written as  * = 1 2 [ ( * , ) +  ( * , ) − || ( * , −  * , )] ,

Figure 4 :
Figure 4: The first fifteen iterations of the cost function.
the initial conditions.The time steps of the forward model and the adjoint model are both 40 seconds.Tidal component  2 at the east open boundary of estuary is given with amplitude of 2 meters.The absorbing extrapolation boundary condition for upstream is used, and the reflection boundary condition is used for solid wall boundary.

Table 2 :
Manning's roughness coefficient () inverse results of assimilating different quantities of water level data.

Table 3 :
Different case of Manning's roughness coefficient inversion.