## Abstract

Based on the Regional Ocean Modeling System (ROMS) and the conditional nonlinear optimal perturbation (CNOP) method, we explore the nonlinear optimal triggering perturbation of the Kuroshio large meander (LM) and its evolution, and reveal the role of nonlinear physical processes in the formation of the LM path. The results show that the large amplitudes of the perturbations are mainly located in the upper 2000 m in the southeastern area of Kyushu (29°–32°N, 131°–134°E), where the eastward propagation of the cold anomaly is vital to the formation of the LM path. By analyzing the depth-integrated vorticity equation of the perturbation, we find that linear advection, namely, the interaction between the perturbation and the reference field, tends to move the cyclonic eddy induced by the optimal triggering perturbation eastward, while the nonlinear advection associated with the interaction of perturbations tends to move the cyclonic eddy westward. The opposing effects of the nonlinear advection and the linear advection slow the eastward movement of the cyclonic eddy so that the eddy has a chance to effectively develop, eventually leading to the formation of the Kuroshio LM path.

## 1. Introduction

The Kuroshio, the well-known western boundary current of the wind-driven subtropical gyre in the North Pacific Ocean, exhibits a remarkable bimodal feature as it flows through the southern region of Japan and forms either a large meander (LM) path or a non-large-meander (NLM) path (or the straight path) (Taft 1972). The NLM path is also sometimes divided into the offshore NLM (oNLM) path and the nearshore NLM (nNLM) path (Fig. 1) (Kawabe 1995). Once the NLM path or the LM path forms, it will persist for a few years to a decade, but the transition between the two types of paths occurs over only a few months (Kawabe 1986, 1987, 1995).

The bimodal feature of Kuroshio south of Japan has important effects on local weather and climate. For example, Xu et al. (2010) noted that the LM path during 2004–05 caused a decrease in local wind speed and precipitation. Nakamura et al. (2012) suggested that the Kuroshio’s bimodal characteristics could significantly affect storm tracks. Hayasaki et al. (2013) indicated that the Kuroshio LM path has an important effect on extratropical cyclone activity. In addition, the Kuroshio path variations also have impacts on fisheries and ship navigation. Hence, research on variations in the Kuroshio’s path is of great significance.

In the past few decades, studies of the Kuroshio path variations south of Japan have mostly focused on the dynamical mechanism. Some previous studies suggested that the particular bottom topography, especially the existence of the Izu Ridge, influences the Kuroshio path variations (Masuda 1982; Chao 1984; Masuda et al. 1999). In addition, some studies have shown that the Kuroshio path selection south of Japan is closely related to the upstream volume transport (Akitomo et al. 1997; Chao and McCreary 1982; Kawabe 1980, 1995; Yoon and Yasuda 1987). However, Qiu and Miao (2000) found that the Kuroshio bimodal feature after 1975 did not clearly correlate with the upstream transport by analyzing long-term observational data. Furthermore, with the aid of a two-layer primitive equation model, they posited that this behavior could be explained as a self-sustained intrinsic oscillation. Additionally, the importance of the baroclinic instability processes in the formation of the Kuroshio LM path has also been emphasized in many studies (Miyazawa et al. 2004; Tsujino et al. 2006; Usui et al. 2008; Fujii et al. 2008).

The above studies show that different factors affecting the formation of the LM have been investigated from the perspective of specific ocean dynamical processes. In fact, the Kuroshio system can be regarded as a nonlinear dynamical system (Schmeits and Dijkstra 2001). In such a system, an initial perturbation may be able to trigger a transition from one steady state to another owing to nonlinear instability, as shown in Mu et al. (2004, 2017) and Kerswell et al. (2014). The initial perturbation that can most easily trigger the transition is called an optimal triggering perturbation. It is therefore interesting to look for the triggering perturbation associated with an LM occurrence because we can easily forecast an LM occurrence if we observe the signal of such an initial perturbation. Furthermore, if the triggering perturbation is found, the roles of some physical processes, especially nonlinear physical processes, in LM formation can be revealed by exploring the nonlinear evolution of the perturbation, which will deepen our understanding of the dynamics of LM path formation.

A key problem is how to search for the nonlinear optimal triggering perturbation. In fact, some researchers have made certain progress regarding initial conditions and transient growth. For example, Moore and Mariano (1999) applied a singular vector analysis to the Gulf Stream in an attempt to investigate the dynamics of perturbation growth, although their optimal perturbation calculation was only valid for several days. Fujii et al. (2008) also adopted the singular vector analysis method to search for the initial perturbation affecting the formation of the Kuroshio LM. However, in their selected reference state, the meander path was already present. Therefore, the initial perturbation they investigated does not refer to the triggering perturbation sought in this study. Wang et al. (2013b) and Zhang et al. (2015) used the conditional nonlinear optimal perturbation (CNOP) method proposed by Mu et al. (2003) to preliminarily investigate the nonlinear optimal triggering perturbation for LM formation and decay. However, their model was a barotropic inflow–outflow model and did not consider realistic topography. That model was too simple to effectively capture the effects of some physical processes, especially those related to nonlinear physical processes.

In this study, we look for the optimal triggering perturbation of the Kuroshio LM path formation based on a realistic high-resolution ocean general circulation model together with the CNOP approach and reveal the role of the nonlinear physical processes (mainly referring to the interaction between the perturbation and the perturbation in our study). The rest of the paper is organized as follows. The model configuration and simulation are described in section 2. In section 3, we will briefly describe the CNOP method and its application to the Kuroshio LM. The results, including the nonlinear optimal triggering perturbation of the LM path and its evolution, are presented in section 4. Finally, a summary and discussion are presented in section 5.

## 2. Model configuration and simulation

We adopted the Regional Ocean Modeling System (ROMS) as the original model to simulate the Kuroshio south of Japan. ROMS is a state-of-the-art hydrostatic, free-surface, three-dimensional primitive equation ocean model that employs a vertical terrain-following coordinate system and horizontal generalized orthogonal curvilinear coordinate system (Song and Haidvogel 1994). It is equipped with various options for advection schemes, parameterizations of horizontal diffusion and vertical mixing, and open boundary conditions (Shchepetkin and McWilliams 2003, 2005). Furthermore, the ROMS adjoint model (Moore et al. 2004) was also used to explore the optimal triggering perturbation of the LM path in this paper.

Sakamoto et al. (2005) and Sato et al. (2006) noted that a horizontal grid resolution of at least 10–20 km is required to reproduce oceanic mesoscale eddies and the western boundary currents in midlatitudes. Smith et al. (2000) also showed the importance of a high resolution when modeling the Gulf Stream. However, because of the lack of observational data and limits on computational resources, we adopted a one-way nested numerical simulation in which a 3:2 grid ratio was used to decrease the grid size in the zonal (meridional) direction from (1/8)° (~14 km) in nest 1 to (1/12)° (~8 km) in nest 2 (Fig. 2). The nest-1 domain spans 20°S–60°N and 100°E–70°W, thus covering the entire North Pacific, whereas the nest-2 domain spans 23°–46°N and 122°–162°E, including the region in which the Kuroshio path variations occur. In the two simulation runs, the model domain was split into 32 sigma levels in the vertical direction, with enhanced resolution near the upper boundary (specified by the same surface stretching parameter *θ*_{s} = 6.0 and bottom stretching parameter *θ*_{b} = 0.1). The critical depth controlling the stretching of the *s* coordinate was 10 m. In addition, the baroclinic-mode time steps were 300 s for nest 1 and 200 s for nest 2. The bottom topography was extracted from the ETOPO2 2-min topographic dataset, for which the maximum depth was set to 5000 m and the minimum depth along the coast was 75 m.

Nest 1 was forced by monthly climatology data with an annual cycle, which was derived from the Comprehensive Ocean–Atmosphere Data Set (COADS) (Diaz et al. 2002). The initial velocity and free-surface elevation were set to zero, and the initial temperature and salinity were prescribed by the monthly climatology of the *World Ocean Atlas* (*WOA2009*; Antonov et al. 2010; Locarnini et al. 2010). All the lateral boundaries were open, and the boundary temperature and salinity were obtained from *WOA2009* while the other variables were calculated through the geostrophic relation. The monthly averaged sea surface height (SSH), temperature, salinity, and velocities were stored as outputs. Because of the limited computational resources, nest 1 was integrated for 50 years, and the final 20-yr outputs were used to drive nest 2.

Before investigating the nonlinear optimal triggering perturbation of the LM path, it is necessary to examine the capability of ROMS to simulate the Kuroshio. The simulated climatology SSH data are shown in Fig. 3a. For comparison, the corresponding climatology SSH data from the satellite altimetry (AVISO) are presented in Fig. 3b. It can be seen from the figures that ROMS successfully reproduced the spatial pattern of the Kuroshio south of Japan and the Kuroshio Extension (KE) system, with a reasonable recirculation gyre in the Shikoku Basin and two quasi-stationary meanders east of 140°E. Although the amplitudes of two KE jet meanders were slightly small, the model showed good skill in simulating the Kuroshio separation between 34° and 35°N from the coast of Japan. Therefore, ROMS can simulate the Kuroshio well.

## 3. CNOP method and its application to the Kuroshio LM

In this section, we show how to search for the nonlinear optimal triggering perturbation of the Kuroshio LM path using the CNOP approach. Actually, the CNOP method can be regarded as a natural extension of the linear singular vector method to the nonlinear system, overcoming the limitation of the linear approximation. In the following, we briefly describe the approach.

Let **X** = (*u*, *υ*, *T*, *s*, *h*)^{T} denote the prognostic variables where *u* (*υ*, *T*, *s*, *h*) represents the zonal velocity (meridional velocity, temperature, salinity, and SSH), and the solution of the nonlinear system at time *t* can be expressed as

where **X**_{0} = (*u*_{0}, *υ*_{0}, *T*_{0}, *s*_{0}, *h*_{0})^{T} is the initial state vector, **X**_{t} = (*u*_{t}, *υ*_{t}, *T*_{t}, *s*_{t}, *h*_{t})^{T} is the state vector at time *t*, and **M**_{t} denotes the nonlinear propagator from the initial time to the future time *t*. When we superimpose an initial perturbation **x**_{0} on the initial state **X**_{0}, the solution can be written as

where **x**_{t} denotes the nonlinear evolution of the initial perturbation **x**_{0} at time *t*.

To seek the initial perturbation that brings about the largest nonlinear evolution at time *t*, a nonlinear constraint optimization problem is defined as

where represents the objective function that measures the nonlinear evolution of the initial perturbation **x**_{0} at time *t*. The denotes the constraint of the initial perturbation, and *δ* is the constraint radius. Here, is the constraint norm and is the objection function norm. The solution **x**_{0δ} to the optimization problem in Eq. (3) is called the CNOP.

To calculate the CNOP, the nonlinear optimization system is built based on the ROMS and the spectral projected gradient (SPG) algorithm (Birgin et al. 2000). The detailed steps of the algorithm used to compute the CNOP are shown in the appendix. In the following, we describe the major preparations of calculating the CNOP:

Selection of the reference state. We focus on the perturbation that triggers the LM from the NLM path in this paper; thus the NLM path is selected as the reference state to calculate the CNOP (viz., the nonlinear optimal triggering perturbation of the LM path). To quantify the Kuroshio path characteristics, we introduce the Kuroshio path index, which is defined as the latitude of the southernmost point of the Kuroshio axis (referring to the 16-cm SSH isoline) within the longitude band from 135° to 140°E, as shown in Fig. 4a. According to the observation data, it is found that when the Kuroshio takes the LM path, the kinetic energy over the area (25°–35°N, 135°–140°E) is large, whereas it is small when the Kuroshio travels a straight path. Therefore, the eddy kinetic energy south of Japan (25°–35°N, 135°–140°E) can also be regarded as a characterization of whether the Kuroshio takes the NLM path or the LM path, as shown in Fig. 4b. The two indexes are obviously consistent, and the Kuroshio always travels on a straight path. Then we arbitrarily selected two cases (denoted as case 1 and case 2) as the reference states to calculate the nonlinear optimal triggering perturbation. It needs to be mentioned that ROMS can also successfully capture the Kuroshio bimodal characteristics when the harmonic horizontal viscosity coefficient for nonlinear model (VIS2) and the harmonic horizontal diffusion of tracer (TNU2) are 40.0 m

^{2}s^{−1}. As the model input parameters, VIS2 and TNU2 have a certain effect on the model results according to our experiments. In addition, the simulated LM formation process is similar to the observed result (figure not shown). However, we are more interested in the perturbation that triggers the LM path after the Kuroshio has followed the NLM path for a long time; thus we adopt the simulation result when VIS2 and TNU2 are 50.0 m^{2}s^{−1}.Selection of the optimization time. According to the observational data and our simulated results, the transition time from the straight path to the Kuroshio LM path is approximately 2–4 months, so we have chosen 70 days as the optimization time in this paper. As shown in Fig. 4a, case 1 extended from 1 September of model year 9 to 10 November of model year 9, and case 2 extended from 1 November of model year 16 to 10 January of model year 17.

- Choosing the suitable initial constraint. Similar to Li et al. (2014), we select the sum of the quadratic perturbations normalized by the standard deviation as the initial constraint. In fact, it is the norm of the initial perturbations and contains the kinetic energy perturbations and involves the potential energy perturbations: This equation represents a ball constraint, and
*δ*is the radius of the ball, which is dimensionless. The and represent the perturbation and the standard deviation of the zonal velocity (meridional velocity, temperature, salinity, and SSH) in the upper 5000 m in the model domain (23°–46°N, 122°–162°E). They are calculated using the monthly outputs during model years 15–20 in nest 2. The initial perturbation in Eq. (4), which actually refers to the initial guess perturbation for the nonlinear optimization system, is obtained from the difference between the monthly averages of two adjacent months. Different initial guess perturbations can result in similar CNOP-type perturbations after several iterations, but the corresponding iteration numbers differ. Usually, about 25 iterations are required for a random initial perturbation, requiring approximately 10 days with 40 CPUs (six cores per CPU) in a Dawning server. - Selection of the objective function. As mentioned above, to use the SPG algorithm, the gradient of the objective function must be supplied. However, for the first path index (measured by latitude) shown in Fig. 4a, the gradient information is difficult to obtain. Moreover, the convergence of the SPG algorithm is better for the L2-norm objective function. Therefore, we chose the kinetic energy (KE) of the nonlinear evolution of the initial perturbation over the area (25°–35°N, 135°–140°E), namely, the second path index (measured by eddy kinetic energy, shown in Fig. 4b), as the objective function. The specific formula can be expressed as where and are the zonal and meridional velocity perturbations at time
*t*, and*ρ*_{ref}denotes the reference density. Note that the ROMS adjoint model is used to provide the gradient information of the objective function while calculating the CNOP. Selection of the constraint radius

*δ.*The number of grid points in nest 2 was 479 × 339 × 32 ≈ 5.2 × 10^{6}. When the amplitudes of velocity, temperature, salinity, and SSH anomalies in the initial perturbation were set to approximately 0.01 m s^{−1}, 0.01°C, 0.01 psu, and 0.001 m, respectively, the dimensionless constraint radius according to Eq. (4) is approximately*δ*= 8 × 10^{8}. However, the perturbations that we are searching for are concentrated in a relatively small area. In addition, we also need to consider the stability of the numerical simulation. Too large a constraint radius will easily induce simulation instability. Given the above factors, we set the constraint radius to*δ*= 1 × 10^{8}_{.}Actually, when the change range of the constraint radius is small, the patterns of the calculated optimal triggering perturbations are similar to each other. Thus, we only show the results for*δ*= 1 × 10^{8}in this paper.

## 4. Results

The nonlinear optimal triggering perturbations for case 1 and case 2 were obtained under the above settings. In this section, we describe the structures of the perturbations, explore their evolution, and reveal the role of the nonlinear physical processes during the LM formation.

### a. The nonlinear optimal triggering perturbation

According to Eq. (4), the norms of the calculated CNOPs for case 1 and case 2 are spatially shown in Fig. 5. Note that they are both dimensionless and vertically integrated from the sea surface to the maximum depth. It can be seen that the large amplitudes of the perturbations are mainly located in the upstream region of the Kuroshio path variations, that is, southeast of Kyushu (29°–32°N, 131°–134°E), although there are small differences in the spatial distributions and magnitudes between the two cases. In addition, because the ROMS is a more realistic ocean model and contains more physical processes than the 1.5-layer shallow-water model, the calculated CNOPs not only possess three-dimensional structures but also differ somewhat in location from those of Wang et al. (2013a), even though they are both in the upstream region. To further investigate the vertical structure of the perturbations, Figs. 6 and 7 show the density (shaded; calculated from the relation between temperature and salinity) and velocity (vector) components in different vertical layers in the perturbation field for case 1 and case 2, respectively. We can see from the figures that the density and velocity perturbations are mainly located in the upper 2000 m, and in particular, the perturbations around (29°–32°N, 131°–134°E) in 800–2000 m are more important in the LM formation process. This is almost consistent with that of Fujii et al. (2008), but our focus in this paper differs from theirs. Specifically, they emphasized the effect of the initial perturbation on LM formation, whereas we are focusing on the perturbation that triggers the LM path from the NLM path. In our study, it is worth noting that the effect of velocity perturbation is minor relative to the temperature and salinity perturbation. Furthermore, salinity makes a larger contribution to density than temperature at the initial time, but the temperature contribution becomes more important over time.

To examine whether the perturbations can trigger the transition from the NLM path to the Kuroshio LM path, we superimposed them on the reference states and integrated the nonlinear model for 70 days. As we expected, the perturbations for case 1 and case 2 caused the transition from the NLM path to the LM path on day 70 (figure not shown). In addition, to demonstrate that nonlinearity is essential in the model integration, we ran another experiment in which the optimal triggering perturbation was halved in amplitude and imposed on the initial field. Taking case 1 as an example, Fig. 8a presents the evolution of the SSH perturbation on day 70, which has been multiplied by a factor of 2, whereas Fig. 8b shows the SSH perturbation evolution on day 70 for the optimal triggering perturbation, in which the black line denotes the Kuroshio axis. Obviously, the amplitude of the meander caused by the halved optimal triggering perturbation is too small to be called the LM path despite being multiplied by a factor of 2. A zsimilar result can also be obtained for case 2. These results directly verify that the model response is indeed nonlinear.

To further investigate the physical processes for the LM formation caused by the nonlinear optimal triggering perturbation, Figs. 9 and 10 separately show the time series of SSH in the perturbation fields for case 1 and case 2, and the thick black line represents the Kuroshio axis. Note that the shading scales are different in the different figures. Here we take case 1 as an example for analysis. On day 10, a negative cyclonic anomaly is located around (31.5°N, 133.5°E) accompanied by two positive anticyclonic anomalies. With the eastward propagation of the negative cyclonic anomaly along the coast of Japan, the positive anomaly in the upstream region gradually weakens, as shown on day 35. Then the further development of the cyclonic eddy causes the meander to form south of Japan, during which the perturbation tends to enlarge the negative anomaly of the meander southward. Finally, on day 70, the Kuroshio LM path occurs after the negative anomaly is well developed. In addition, although the positive anomaly in the downstream region is not large compared to the cold anomaly, it cannot be ignored. With the formation of the meander, the positive anomaly remains in the same position on the downstream side of the meandering path. Actually, the positive and negative anomalies positioned on the upstream and downstream sides of the meander tend to slow its propagation by displacing the meander backward. As noted in Fujii et al. (2008), when the propagation is slow, the meander tends to grow larger. For case 2, we obtain similar findings. Besides, we find that the positive anomaly downstream of the meander in case 1 is slightly weaker than that in case 2, which may explain why the LM in case 1 is farther to the east than its counterpart in case 2.

### b. The role of nonlinear physical processes

As stated above, the nonlinear evolution of the optimal triggering perturbation successfully leads to the transition from the NLM path to the LM path. The following question then arises: What role do the nonlinear physical processes play in the transition?

To address the above question, we adopt the following momentum equations for the perturbations, by means of which the nonlinear terms are presented:

where **V** = (*u*, *υ*) represents the horizontal velocity vector of the reference state, and the prime denotes the corresponding perturbation field. The is the nonlinear advection vector of the perturbation by the perturbation, which is exactly the nonlinear physical processes we are referring to. Because we are concerned about the role that the nonlinear physical processes play in the LM formation in this paper, we do not discuss the other terms apart from the nonlinear advection term in Eqs. (6) and (7). As shown in Fig. 11, the vector diagrams of the nonlinear term in the upper 2000 m for case 1 and case 2 were averaged for the last 30 days. The blue lines denote the Kuroshio axis in the reference state, whereas the other lines represent the Kuroshio axis caused by the optimal triggering perturbation at different times in the last 30 days. The figures show that the amplitude of the meander path indeed increases over time. The strong anticyclonic eddy on the western side implies the eastward propagation of the meander path, while the anticyclonic eddy on the eastern side implies the narrowing of the cold cyclonic eddy, thereby preventing the meander from moving downstream. However, the narrowing of the cold cyclonic eddy corresponding to the core of the LM only occurs in the east–west direction, while the meander path gradually expands toward the southeastern direction. That is, the amplitude of the meander is increasing, and the nonlinear advection associated with is associated with southeastward enlargement the Kuroshio path south of Japan. In fact, the vector diagrams just intuitively show the direction of the effect of the nonlinear physical processes in the formation of the LM path. The dynamical mechanism of the nonlinear processes requires further investigation.

Some theories regarding the nonlinear advection effect have been proposed. For example, Cushman-Roisin et al. (1993) confirmed two opposing mechanisms governing the temporal evolution of meanders for the Gulf Stream. They stated that the advection term and the beta term of the depth-integrated vorticity equation represent the major balance for the LM. In addition, Tsujino et al. (2006) also reached a similar conclusion based on a high-resolution general circulation model. In their studies, the stated advection effect was associated with the interaction between the perturbation and the mean current field or the combined role of perturbation–mean current field interaction and perturbation–perturbation interaction. Actually, from the perspective of the perturbation evolution equation, the advection term can be divided into linear advection (the perturbation–current interaction) and nonlinear advection (the perturbation–perturbation interaction). In our study, nonlinear advection refers to the perturbation–perturbation interaction, which is a part of the advection effect mentioned in the previous studies. Therefore, we also introduce the depth-integrated relative vorticity equation to explore the dynamical mechanism of the nonlinear physical processes. For diagnostic analysis, the simplified depth-integrated vorticity equation still explains the main balanced forces influencing the LM path:

where is the vorticity of the reference field. The *A*_{M} is the horizontal viscosity coefficient, *f* is the Coriolis parameter, and *β* is the meridional derivative of *f*. The prime denotes the corresponding perturbation caused by the nonlinear evolution of the optimal triggering perturbation. The terms on the right-hand side in Eq. (8) separately represent the advection terms, the convergence and divergence terms, the beta effect, and the viscosity term. The advection terms include the linear advection, that is, the perturbation–reference flow interaction *J*_{1} and the nonlinear advection term of the perturbation by the perturbation *J*_{2}.

To verify the feasibility of using the simplified vorticity equation to identify the main balanced forces influencing the LM path, we compared the right-hand tendency and the left-hand tendency based on the model monthly averaged outputs and found that their spatial patterns are similar to each other (figure not shown). Therefore, the right-hand terms in Eq. (8) basically represent the vorticity perturbation budget. The next step is to clarify the dynamical mechanism of the nonlinear processes by analyzing each term on the right-hand side of Eq. (8).

As shown in Figs. 12 and 13, the contours denote the actual anomalous vorticity field, and the colors represent the anomaly terms in the vorticity equation. Each term in Eq. (8) was integrated from 2000 m to the sea surface based on the 10-day average fields during the LM formation. Considering the growth rate of the meander for case 1 and case 2 is a little different; each term of the vorticity equation for case 1 has been averaged from day 31 to day 40, whereas those for case 2 have been averaged from day 41 to day 50. The figures clearly show that the anomalous vorticity in the two cases is positive in the core of the meander but negative on its eastern flank, and the tendency in the upper-left panel is negative on the western side of the meander but positive on the eastern side, indicating that the meander is propagating eastward along the coast of Japan. Moreover, comparison of the right-hand terms in Eq. (8) reveal that the linear and nonlinear advection terms are dominant, representing the main balance of forces in the vorticity equation during the LM formation. In the growth of the cyclonic eddy corresponding to the core of the LM, the linear advection term plays the dominant role, although the positive divergence terms off the southern tip of the meander in Figs. 12 and 13 also partly explain the local growth of the cyclonic eddy. However, during the propagation of the cyclonic eddy, the role of the nonlinear advection also becomes very important. As indicated by the black arrows in Figs. 12 and 13, the spatial patterns of the linear and nonlinear advection terms are almost opposite. The linear advection is negative in the western region of the meander, which means a weakening of the local cyclonic eddy, whereas it is positive in the eastern region of the meander, which means a strengthening of the local cyclonic eddy. In contrast, the nonlinear advection is positive in the western region of the meander but negative in the eastern region of the meander. Based on these findings and the anomalous vorticity, we conclude that the linear advection attempts to make the cyclonic eddy propagate eastward, while the nonlinear processes attempts to make the cyclonic eddy propagate westward. The counteracting effects of the nonlinear processes and the linear advection slow down the eastward movement of the cyclonic eddy, allowing it to effectively develop. Eventually, sufficient evolution of the cyclonic eddy leads to the formation of the Kuroshio LM path. Therefore, the nonlinear physical processes indeed play an important role in promoting the formation process of the LM path. It needs to be mentioned that the spatial pattern of the beta effect is similar to that of the nonlinear advection, also contributing to the westward propagation of the meander, although the beta term is very small.

## 5. Summary and discussion

Based on the ROMS and the CNOP method, the nonlinear optimal triggering perturbation of the Kuroshio LM path and its evolution have been investigated. Prior to this analysis, we first checked the ability of ROMS to simulate the Kuroshio. ROMS can effectively reproduce the spatial pattern of the Kuroshio south of Japan and the KE system, with a reasonable recirculation gyre in the Shikoku Basin and two quasi-stationary meanders east of 140°E. In addition, ROMS can also successfully capture the bimodal feature of the Kuroshio south of Japan using a fine-resolution model nested in a coarse-resolution North Pacific Ocean model. In addition, the formation process of the Kuroshio LM path in our simulation is similar to that based on the observational data.

Furthermore, we focused on the perturbation that triggers the LM path from the NLM path; thus we chose case 1 and case 2 in the straight path as the reference field to calculate the CNOPs. By investigating the spatial and vertical patterns of the optimal triggering perturbations for case 1 and case 2, we found that the large amplitudes of the perturbations are mainly located in the upper 2000 m in the southeastern region (29°–32°N, 131°–134°E) of Kyushu. In addition, the time series of SSH in the perturbation field show that the propagation along the Kuroshio axis downstream of the cold cyclonic eddy east of Kyushu causes the meander to further develop, eventually leading to the LM formation on day 70. This process is similar to that in Fujii et al. (2008) and Wang et al. (2013a), although the reference state selected for our study differs from theirs.

The CNOP method is fully based on a nonlinear model, unlike linear methods such as the first singular vector (FSV), which represents the fastest linearly growing initial perturbation; so the role of the nonlinear processes in the formation of the LM path was revealed in this paper. Above all, the importance of the nonlinear physical processes is stated through the nonlinear advection vector diagrams of the momentum equations. To further analyze the physical mechanism of the nonlinear processes, we introduced the perturbation vorticity budget equation, which could represent the main balance of forces influencing the Kuroshio LM path. The results showed that the nonlinear term slows the eastward propagation caused by the linear advection of the interaction between the perturbation and the reference flow. This causes the cyclonic eddy to sufficiently develop, which prompts the transition from the NLM path to the LM path. Therefore, we concluded that the nonlinear processes act as an important positive role in the LM formation. It needs to be emphasized that our conclusion is not contradictory to that of Fujii et al. (2008), who used a linear analysis, because in their study, the perturbation evolution was still obtained by integrating a nonlinear model although the initial perturbations were obtained by means of a linear method. In other words, the nonlinear effect was already present in the LM formation. In addition, we also investigated the role of the baroclinic and barotropic instability in the LM formation from the perspective of energetics analysis and found that the baroclinic instability indeed imposes a greater effect on the growth of the cyclonic eddy; however, this was not our focus, so we presented no further discussion in this paper.

In summary, this paper not only presents the three-dimensional optimal triggering perturbation of the LM path but also reveals the important role of the nonlinear physical processes (the interaction of the perturbations) in the LM formation based on a realistic ocean circulation model for the first time. This further deepens our understanding of the transition process from the NLM path to the LM path. However, considering the computational limitations, we present the results for only two cases here. In the future, more cases will be chosen and discussed. Furthermore, we found that the patterns of the obtained nonlinear optimal triggering perturbations for case 1 and case 2 are localized, as noted by Wang et al. (2013a), and this finding will help us effectively identify the areas of sensitivity for targeted observations of LM path predictions. In addition, the vertical distribution of the perturbations, which cannot be obtained based on the simple model by Wang et al. (2013a,b), also offers us useful guidance for implementing targeted observations. These factors will be our focus in a later study.

## Acknowledgments

This study was supported by the National Natural Scientific Foundation of China (41230420, 41576015), the Qingdao National Laboratory for Marine Science and Technology (QNLM2016ORP0107), the NSFC Innovative Group Grant (41421005), the NSFC–Shandong Joint Fund for Marine Science Research Centers (U1606402), and the National Programme on Global Change and Air-Sea Interaction (GASI-IPOVAI-06).

### APPENDIX

#### Computation of the CNOP by the Spectral Projected Gradient Algorithm

The SPG algorithm was used to solve the minimum problem of a nonlinear function subject to a set of constraints on the variables. Actually, the derivations of the algorithm have been presented by Birgin et al. (2000). However, to more clearly describe how the SPG algorithm is used to obtain the CNOP, we repeat the procedure. The main steps are as follows (Wang et al. 2012):

Step 1: Given an integral , set the parameters and ; set a constant parameter and ; and is arbitrary. Given a starting iterate point , if , replace by . The is an projection operator and defined as . Set , where

*k*indicates iteration step.Step 2: If , the iteration stops and go to step 8.

Step 3: Compute .

Step 4: Determine

*λ*_{k}using the following line search algorithm.Step 4.1: Compute , , , where 〈〉 denotes the inner product of two variables; set .

Step 4.2: If , the line search stops and go to step 4.6.

Step 4.3: Compute .

Step 4.4: If ( and ), set ; else set

*λ*=*λ*/2.Step 4.5: Compute ; return to step 4.2.

Step 4.6: Set

*λ*_{k}=*λ.*Step 5: Compute , , , and .

Step 6: If , set ; else .

Step 7: Set

*k*=*k*+ 1; go to step 2.Step 8: Set , where

*x*_{*}represents the optimal solution.

Therefore, to compute the CNOP in Eq. (3) using the SPG algorithm, we need to rewrite Eq. (3) as ; namely, . In our study, **J**_{1} denotes the objection function *f* in the SPG algorithm, and the defined constraint condition in Eq. (4) is equivalent to the projection operator *P*_{Ω}. In addition, the parameters in step 1 of the algorithm are set to , while in step 2 of the algorithm, the parameter *ε* is set to 5 × 10^{−3}. In summary, three subroutines that define the objective function, the constraint conditions, and the gradient of objective function need to be supplied. Here the ROMS adjoint model is used to provide the gradient information of objective function respect to the initial perturbations. Then the initial perturbations will gradually approach the optimal solution along the gradient direction during the iteration.

## REFERENCES

*Salinity*. Vol. 2,

*World Ocean Atlas 2009*, NOAA Atlas NESDIS 69, 184 pp.

*Temperature*. Vol. 1,

*World Ocean Atlas 2009*, NOAA Atlas NESDIS 68, 184 pp.

*Kuroshio: Its Physical Aspects*, H. Stommel and K. Yoshida, Eds., University of Tokyo Press, 165–216.

## Footnotes

This article is included in the Ocean Turbulence Special Collection.

© 2018 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).