Mathematical Problems in Engineering

Volume 2015, Article ID 309525, 11 pages

http://dx.doi.org/10.1155/2015/309525

## Estimation of Oceanic Eddy Viscosity Profile and Wind Stress Drag Coefficient Using Adjoint Method

^{1}Laboratory of Physical Oceanography, Ocean University of China, Qingdao 266100, China^{2}State Key Laboratory of Satellite Ocean Environment Dynamics, Second Institute of Oceanography, SOA, Hangzhou 310012, China

Received 8 March 2015; Revised 29 May 2015; Accepted 15 June 2015

Academic Editor: Carla Faraci

Copyright © 2015 Qilin Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

Adjoint method is used to assimilate pseudoobservations to simultaneously estimate the OEVP and the WSDC in an oceanic Ekman layer model. Five groups of experiments are designed to investigate the influences that the optimization algorithms, step-length, inverse integral time of the adjoint model, prescribed vertical distribution of eddy viscosity, and regularization parameter exert on the inversion results. Experimental results show that the best estimation results are obtained with the GD algorithm; the best estimation results are obtained when the step-length is equal to 1 in Group 2; in Group 3, 8 days of inverse integral time yields the best estimation results, and good assimilation efficiency is achieved by increasing iteration steps when the inverse integral time is reduced; in Group 4, the OEVP can be estimated for some specific distributions; however, when the VEVCs increase along with the depth at the bottom of water, the estimation results are relatively poor. For this problem, we use extrapolation method to deal with the VEVCs in layers in which the estimation results are poor; the regularization method with appropriate regularization parameter can indeed improve the experiment result to some extent. In all experiments in Groups 2-3, the WSDCs are inverted successfully within 100 iterations.

#### 1. Introduction

The development of ocean observation technology provides a wealth of fruits, especially the satellite data as well as acoustic doppler current profilers (ADCP) data. It is a mass of observations that make the application of data assimilation techniques to oceanography possible. In recent years, researchers have developed some data assimilation methods, such as successive corrections [1], optimal interpolation [2, 3], Kalman filtering [4, 5], and variational assimilation method (most notably adjoint method) [6].

During the past few decades, many researches have showed that the adjoint method occupies an increasingly important position in data assimilation [7]. The adjoint method is an advanced method with a combination of variation principle and optimal control theory [8]. The method has an advantage of directly assimilating various observations distributed in time and space into numerical models while maintaining dynamical and physical consistency with the model. The adjoint method can be used to solve numerous different types of practical problems, such as estimating model parameters, inverting initial conditions, and researching steady flows and circulation problems [9, 10] and surface heat flux problems [11]. So far, the adjoint method has been widely used for parameter estimation and initial conditions inversion. Yu and O’Brien [12] estimated the wind stress drag profile and the oceanic eddy viscosity profile (OEVP) based on a modified Ekman layer model. Lardner et al. [13] estimated the bottom drag coefficients for a two-dimensional numerical tidal model of the Arabian Gulf. Chen et al. [14] constructed an adjoint model for simulation of internal tides and simply tested with a series of ideal experiments in which several prescribed spatial distributions of OBCs were successfully inverted by assimilating the model-generated pseudoobservations. Later on, studied the estimation of open boundary conditions and the bottom friction coefficients based on the adjoint model. Cao et al. [15] optimized the open boundary conditions in a 3D internal tidal around Hawaii. Seiler [16] estimated open boundary conditions with the adjoint method based on a quasi-geostrophic ocean model for a midlatitude jet. Yu and O’Brien [17] researched the initial conditions adjustment and parameter estimation simultaneously with the adjoint method. Peng et al. [18, 19] corrected the errors in the initial conditions in storm surge simulation using the adjoint method. Predecessors have done a vast amount of work to research the OEVP. Lentz [20] studied the sensitivity of the inner-shelf circulation to the form of the eddy viscosity profile. Pohlmann [21] calculated the annual cycle of the vertical eddy viscosity in the North Sea with a three-dimensional baroclinic shelf sea circulation model. Tian et al. [22] studied the vertical distribution of diapycnal diffusivity in the South China Sea, the Luzon Strait, and the North Pacific based on the observations.

In this paper, the adjoint method in conjunction with an oceanic Ekman layer model is used to obtain the optimal estimates of the vertical distribution of eddy viscosity and the surface wind drag coefficients from pseudoobservations of wind velocities and current velocity profiles within the water columns. At present, the error of ADCP survey technique could be controlled within 5 meters. Therefore, the Ekman layer (the thickness is 100 meters in this paper) could be divided into 20 layers evenly.

To solve the problem above, good many feasible large-scale optimization algorithms can be found in [23]. However, the studies discussing the performances of various optimization algorithms in the meteorological and oceanographic application are still relatively rare. The gradient descent (GD) algorithm is a fundamental and probably the earliest algorithm for nonlinear optimization so that it attracts many researchers’ attention since it was proposed. It has the negative gradient direction as the search direction and is one of the simplest minimization algorithms that require calculation of derivatives. The GD algorithm is particularly useful when the dimension of the problem is very large. The limited-memory BFGS (L-BFGS) algorithm is an adaptation of the standard BFGS algorithm to large-scale optimization problems. This algorithm has been proved to be globally convergent for convex objective functions and provides a fast rate of linear convergence as well as requiring minimal storage. The conjugate gradient (CG) algorithm is an algorithm that establishes a link between conjugate direction algorithm and steepest descent algorithm, in order to obtain a good algorithm which is effective and has good convergence.

The organization of this paper is as follows. The numerical model is given briefly in Section 2. In Section 3, five groups of numerical experiments and results analysis are performed. To be specific, in the first group, three optimization algorithms are used to investigate which is more effective; in the rest of the groups, the influences that step-length, inverse integral time of the adjoint model, prescribed vertical distribution of eddy viscosity, and regularization parameter exert on the inversion results are investigated, respectively. A summary of the results and conclusions are given in Section 4.

#### 2. Numerical Model

A horizontally unbounded ocean surface layer with depth is considered. The -axis points vertically upwards, and the sea surface is located at . The Coriolis parameter, , is taken to be a constant. The oceanic Ekman layer model is governed by the following equations:The corresponding boundary conditions areThe initial conditions arewhere (positive to the east) and (positive to the north) are horizontal velocity components in ocean, while and are the wind velocity components. Assume that the wind force is only along the -direction, and so is 0. is the density of water, and is the density of air. is the vertical eddy viscosity coefficients (VEVCs) while is the wind stress drag coefficient (WSDC), which are the parameters to be estimated. The VEVC is a function of depth. The numerical model is formulated using a difference discretization on a grid with spatial increment . Model parameters are set as follows: , , , , , , , and .

A cost function which quantifies the discrepancy between the model results and the observations can be constructed aswhere the caret denotes the observation values. is the weighting matrix and theoretically should be the inverse of observation error covariance matrix [24]. However, determining the “exact” form is far from easy [25]. The general practice is to assume that each observation is not relevant and is equally important. So the covariance matrices can be set to unit matrices. Therefore, is simplified to be unit matrix [26].

Equations (1a), (1b), (1c), and (1d) can be enforced by introducing a set of undetermined Lagrange multipliers, which are essentially the adjoint variables. This leads to the formation of the Lagrange function, which is given bywhere and are Lagrange multipliers for the -component and -component, respectively.

The aim is that the cost function is minimized subject to the initial boundary value problem (1a), (1b), (1c), and (1d). For this purpose, the adjoint equations derived from the Lagrange multiplier method are given as follows:The corresponding boundary conditions are The corresponding initial conditions arewhere ( days) is the total integral time of this oceanic Ekman layer model. The currents are simulated by integrating positive (i.e., oceanic Ekman layer) model, while the control parameters are optimized by integrating inverse (i.e., adjoint) model. In general, the inverse integral time of inverse model is equal to the integral time of positive model. The gradient of the Lagrange function with respect to and yields

As shown in Section 1, there are many feasible large-scale optimization algorithms. Take the GD algorithm [27] for instance, the control parameters can be optimized using the following formulas:where is the cost function which quantifies the discrepancy between the model results and the observation values and denotes the normalized gradient with respect to , andwhere and are step-lengths which are used to adjust the parameters preferably.

#### 3. Numerical Experiments and Results Analysis

In the first place, the performances of GD, CG, and L-BFGS are compared via a group of ideal experiments in Section 3.1. In the CG algorithm, the coefficient of search direction employs FR formula, which is showed by Fletcher and Reeves [28]. The FR formula is given as follows:where is the coefficient of search direction of th-iteration and is the gradient of the cost functions with respect to control parameters of th-iteration.

As we all know, step-length affects the result involved with numerical optimization. If VEVC is estimated successfully when the inverse integral time of the adjoint model is reduced, the workload of calculation will be saved greatly. The actual oceanic eddy viscosity profile cannot be observed and it is unknown. If different prescribed distributions of VEVCs are estimated successfully, it can be inferred that the model based on the adjoint method has a good ability to estimate the VEVCs. Based on the above considerations, the following four groups of experiments are carried out in order to investigate the influences that the three optimization algorithms, step-length, inverse integral time of the adjoint model, and prescribed distribution of the VEVCs exert on the inversion results, respectively. In all four groups, the given value of the WSDC is set to be 0.0012. Therefore, the effect of is not investigated in this paper, and the related step-length, , is set to be a constant. The initial values of VEVC and WSDC are 0.008 and 0.0007, respectively. In the first two groups, the prescribed OEVP refers to Yu and O’Brien [12].

##### 3.1. Group 1

The performances of GD, CG, and L-BFGS algorithms are compared via a series of ideal experiments in this group. Because of simplicities of the WSDCs in this model, the experiments are performed with only VEVCs briefly. The inverse integral time of the adjoint model is 10 days. In the process of parameter optimization, we tried different iteration numbers. And the best result in experiments with different iteration numbers is adopted as the result of the algorithm. Clearly, the step-length is set to 1 as a matter of priority [23].

The prescribed and inverted VEVCs and variation curves of cost functions are shown in Figure 1. The root mean square errors (RMSEs) between VEVCs before and after assimilation are listed in Table 1.