Effects of eddies on the subduction and movement of water masses reaching the $137^{\circ }\,\hbox {E}$ section using Lagrangian particles in an eddy-resolving OGCM

Abstract

The effects of eddies on the subduction and movement of water masses reaching the 137$^{\circ }$ E section are examined in a nominal 10-km resolution ocean general circulation model using a backward-particle tracking method. Target water masses are the Tropical Water (TW), the Eastern Subtropical Mode Water (ESTMW), the Subtropical Mode Water (STMW), and lighter variety of the Central Mode Water (L-CMW). Each particle is classified as a typical water mass according to its physical properties in the subduction area, and into eddy and non-eddy components based on the Okubo–Weiss parameter. During subduction, each water mass tends to be located in anticyclonic eddies rather than cyclonic eddies. The effects of eddies on the spatial distribution of water mass and the time taken by the water mass to reach the 137$^{\circ }$ E section differ for each water mass. For the TW, the water mass in the mesoscale eddies tends to be distributed along the eastward Subtropical Countercurrent (STCC), which moves eastward from the 137$^{\circ }$ E to Hawaii. The eddy component takes lesser time to reach the 137$^{\circ }$ E section as compared to the non-eddy component. For ESTMW, a similar pattern appears around STCC but its effect is confined to the west of Hawaii. For STMW, subduction and distribution occur predominantly in anticyclonic eddies. A part of L-CMW crosses the Kuroshio Extension when cyclonic eddies are pinched off from the troughs of the Kuroshio Extension, reaching the 137$^{\circ }$ E section within 2 years.

Appendices

Appendix 1: particle-tracking method

The particle-tracking algorithm was originally developed by Matsumura and Ohshima (2015) and incorporated in MRI.COM. This method achieves high-speed processing by reducing unnecessary memory access as much as possible, by considering that the speed of the supercomputer is limited by the memory access speed, as in the current architecture. The algorithm calculates the trajectory of each particle using the fourth-order Runge–Kutta method through the linearly interpolated three-dimensional velocity field within the model grid. In MRI.COM, this algorithm can be operated both online and offline. Although effects such as random walk for representing diffusion can be incorporated, they are not used in this study.

The model applies a convective adjustment scheme, comparing every two vertically adjacent grids, and if unstable, a large vertical diffusion (1.0 $\hbox {m}^2$ $\hbox {s}^{-2}$) is set between the grids. We consider the deepest grid among the unstable grids as the mixed layer in the particle-tracking scheme. We use the criterion of the MLD that the depth at which the $\sigma _0$ changes from the surface is 0.03 kg $\hbox {m}^{-3}$ greater than that of the surface layer only for comparison with the observation in Fig. 4. These two criteria show similar MLD patterns (not shown).

In case of online calculation, particles flow at the same time interval as the physical field is updated each time step. In case of offline calculation, the trajectory of the particles is calculated by time-interpolating the physical field saved at predetermined intervals. However, in principle, the offline calculation can be performed in exactly the same manner as the online calculation, if the outputs of the physical fields are saved at every time step; however, this requires enormous computer resources and is unfeasible.

Here, we examine which interval is adequate for the offline calculation. Five-day and 1-day intervals are useful because they can divide one year (365 days) into regular intervals and are frequently used for the storage intervals of OGCMs. We examine the performance of the 5-day average fields and 1-day average fields, by conducting five experiments. Each experiment comprises a 1-year forward-tracking release of 49 particles per grid from the 165$^{\circ }$ E section. Two offline calculations are conducted with the 5- and 1-day averaged fields, which are used with linear interpolation with time. In these experiments 49 particles are evenly distributed within a grid cell, while the other three online calculations differ only in the initial position of the particles. One initial condition is the same as that of the offline calculations. For the other two initial conditions, the horizontal position of each particle within a grid is randomly distributed.

Figure 18a shows the horizontal distribution of the number of particles after the 1-year release of one of the online calculations. The other four experiments show similar patterns (not shown), but are slightly different. The standard deviation ($\sigma$) of the three online calculations is calculated as the square root of the unbiased dispersion of the three experiments. The experiment using the 1-day average fields is within 1$\sigma$ in most part (Fig. 18b). Even in the strong mesoscale eddy at 150$^{\circ }$ E after the 1-year release, the error is within 2$\sigma$. The performance of the experiment conducted using the 5-day average fields might be acceptable, but the error tends to exceed 1$\sigma$ and sometimes 2$\sigma$. Because we mainly examine the eddy effect in this study, we decide to use the 1-day-average fields.

Fig. 18

a Number of particles 1-year after being released from the 165$^{\circ }$ E section for the online calculation. b The blue (red) lines show the difference in number of particles along 29$^{\circ }$ N between the online and offline calculations using the 1-day (5-day) average fields. The dark gray shade shows the standard error ($\sigma$) estimated from the three online calculations. The light gray shade show 2$\sigma$

Appendix 2: trajectories and propagation speed of eddies using the OW-parameter

To confirm the fidelity of the use of the OW-parameter (Fig. 19), we show that the trajectories and speed of eddies based on the OW-parameter are largely consistent with the previous studies. The following procedures are conducted to calculate the trajectories of eddies using the modeled OW parameters at the surface. Note that this procedure is conducted independently of the Lagrangian particle tracking method to confirm that the validity of the OW parameter and performance of the modeled eddies are realistic.

Fig. 19

The black lines show the contours of $OW = 2.0\times 10^{-12}$ $\hbox {s}^{-1}$ of AVISO on 1 June, 2017. The red and blue circles show the anticyclonic and cyclonic eddies identified in Mesoscale Eddy Trajectory Atlas (2019) on the same day. The radius of each circle is calculated so that each area circle is the same as the eddy area identified in Mesoscale Eddy Trajectory Atlas (2019)

In the area of $OW < \mathrm{{ow}}_\mathrm{{cr}} (-2.0 \times 10^{-12}$ $\hbox {s}^{-2})$, when a maximum (minimum) of $h' = h - {\overline{h}}$ is found, it is considered as the center of an anticyclonic (cyclonic) eddy. In the next 5-day averaged fields, we attempt to find whether such a maximum (minimum) exists in $OW < \mathrm{{ow}}_\mathrm{{cr}}$. We search for such a maximum (minimum) within 1.5$^{\circ }$ longitude and 0.5$^{\circ }$ latitude centered on the anticyclonic (cyclonic) eddy. When we find such a maximum (minimum), we consider that the eddy center moves to this position. When we cannot find such a maximum (minimum), we consider that the eddy has collapsed. Figure 20 shows the trajectories of eddies with a lifespan greater than 120 days. We can see a reduced number of eddies along the axes of KE, predominantly cyclonic north of it, and the predominance of anticyclonic eddies on its poleward side. Various places are scarce of long-lived mesoscale eddies, such as the northern North Pacific, and the zonally banded area between 15$^{\circ }$ and 20$^{\circ }$ N west of Hawaii. This pattern is consistent with the trajectory map of Chelton et al. (2011). We then estimate the propagation speed of eddies by dividing the first and last positions of the eddy trajectories by the propagation time (Fig. 21). We assume that during the propagation time, the eddy moves with the constant velocity needed for calculating the horizontal distribution of the propagation speed (Fig. 22).

Fig. 20

Trajectories of cyclonic (blue lines) and anticyclonic (red lines) eddies during the 2013–2018 in the model for lifetimes > 120 days (about 16 weeks). The starting point of each trajectory is stressed

Fig. 21

The red dots indicate the latitudinal variation of westward zonal propagation speeds estimated by each eddy trajectory whose lifetime $\ge$ 115 days. The blue shaded line is the latitudinal profile of the zonally averaged westward phase speeds of long baroclinic Rossby waves

Fig. 22

The shading shows the estward zonal propagation speeds estimated by each eddy trajectory whose lifetime $\ge$ 115 days. Unit is cm $\hbox {s}^{-1}$. The black contours show the average of SSH between 1993 and 2017. The contour interval is 0.3 m

Appendix 3: post-processing of particle analysis

Each particle has a unique ID that increases every year. On August 11, 2017, the ID is numbered $1, \cdots , n_0$. The particle group released this year is labeled as $E = E_{2017}$. We express $n \in E_{2017}$, when a particle whose ID is n belongs to the particle group $E_{2017}$. Similarly, on August 11, 2016, the ID is numbered $n_0+1, \cdots , 2 n_0$. The particle group released this year is labeled as $E = E_{2016}$.

Then the particle ID, n, is associated with the typical water mass using PV, potential temperature ($\theta$), salinity(S), longitude, latitude, MLD, and velocity when it is ventilated (Table 2). In this study, we consider that the particle does not change its associated transport and water mass category. When a particle whose ID n belongs to one of the four water masses, it is written as $n \in W$, where W is the group corresponding to the water mass.

To obtain the horizontal distribution, we coarse-grain the particle position corresponding to the water mass on a $0.5^{\circ }\times 0.5^{\circ }$ grid rather than on the native $1/11^{\circ }\times 1/10^{\circ }$ grid of the model. We choose this $0.5^{\circ }\times 0.5^{\circ }$ resolution to eliminate gaps but retain sufficient details after trial and error.

The subduction rate is calculated by the summation of the transport ($\psi _n$) of the ventilated particles.

$$\begin{aligned} S(I,J,W, E) = \sum _{n\in W, n\in E} \psi _n \qquad \{n \mid \mathbf {x}_n(t_n^\mathrm{{sbd}}) \in G(I,J) \} , \end{aligned}$$

(4)

where subscript n is the particle ID, $t^\mathrm{{sbd}}$ is the time of subduction, $\mathbf {x}(t)$ is the horizontal location of the particle, and G(I, J) is a grid of the $0.5^{\circ }\times 0.5^{\circ }$ resolution whose indices are I, J. By comparing the S in different E, we can examine how differences in S contribute to differences in the 137$^{\circ }$ E section, which will be examined in the future. The climatology of the subduction is calculated by the average for E, $\overline{S{(I,J,W, E)}}^E$.

We investigate how each water mass is selectively created by the anticyclonic eddies, by integrating the assigned transport for each bin of the OW parameter,

$$\begin{aligned} \varPhi _H^\mathrm{{sbd}}(k, W, E)&= \sum _{n\in W, n\in E} \psi _n \qquad \{n\mid \mathrm{{ow}}_n(t_n^\mathrm{{sbd}}) \in \mathrm{{OW}}_k, \zeta _n(t_n^\mathrm{{sbd}}) < 0\}, \end{aligned}$$

(5)

$$\begin{aligned} \varPhi _L^\mathrm{{sbd}}(k, W, E)&= \sum _{n\in W, n\in E} \psi _n \qquad \{n\mid \mathrm{{ow}}_n(t_n^\mathrm{{sbd}}) \in \mathrm{{OW}}_k, \zeta _n(t_n^\mathrm{{sbd}}) > 0\}, \end{aligned}$$

(6)

where $\mathrm{{ow}}_n(t)$ is the value of the OW-parameter for the particle whose ID is n and time is t, $\mathrm{{OW}}_k$ is the interval corresponding to the kth bin of the OW parameter, L and H present cyclonic and anticyclonic, and $\zeta$ is the relative vorticity. The climatology is calculated by the average for E, $\overline{\varPhi _{{\mathrm{H}}}^\mathrm{{sbd}}(k, W, E)}^E, \overline{\varPhi _{{\mathrm{L}}}^\mathrm{{sbd}}(k, W, E)}^E$. Similarly, we investigate how such signals remain along the 137$^{\circ }$ E section by calculating

$$\begin{aligned} \varPhi _H^{137E}(k, W, E)&= \sum _{n\in W, n\in E} \psi _n \nonumber \\&\{n\mid \mathrm{{ow}}_n(t_n^{137E}) \in \mathrm{{OW}}_k, \zeta _n(t_n^{137E}) < 0\}, \end{aligned}$$

(7)

$$\begin{aligned} \varPhi _L^{137E}(k, W, E)&= \sum _{n\in W, n\in E} \psi _n \nonumber \\&\{n\mid \mathrm{{ow}}_n(t_n^{137E}) \in \mathrm{{OW}}_k, \zeta _n(t_n^{137E}) > 0\}, \end{aligned}$$

(8)

where $t^{137E}$ is the time of the release from the 137$^{\circ }$ E section. The climatology is calculated by the average for E, $\overline{\varPhi _H^{137E}(k, W, E)}^E, \overline{\varPhi _L^{137E}(k, W, E)}^E$.

We visualize the horizontal location of each water mass by conducting the coarse-graining separately for the particles inside or outside eddies.

$$\begin{aligned} \varPsi ^\mathrm {eddy}(I,J,t, W, E)&= \sum _{n\in W, n\in E} \psi _n \nonumber \\&\{n\mid \mathrm{{ow}}_n(t) < \mathrm{{ow}}_\mathrm{{cr}}, \mathbf {x}_n(t) \in G(I,J)\}, \end{aligned}$$

(9)

$$\begin{aligned} \varPsi ^\mathrm {noneddy}(I,J,t, W, E)&= \sum _{n\in W, n\in E} \psi _n \nonumber \\&\{n\mid \mathrm{{ow}}_n(t) > \mathrm{{ow}}_\mathrm{{cr}}, \mathbf {x}_n(t) \in G(I,J)\}. \end{aligned}$$

(10)

We also call $\varPsi ^{{\mathrm{eddy}}}$ an eddy component of the water mass, and $\varPsi ^{{\mathrm{noneddy}}}$ a non-eddy component of the water mass. These numbers have the unit of transport [$\hbox {m}^3$ $\hbox {s}^{-1}$]. When spatially integrated, they demonstrate the time series of the water-mass subduction gathering toward the 137$^{\circ }$ E section for the eddy and non-eddy components,

$$\begin{aligned} \varXi ^{\mathrm {eddy}} (t, W, E)&= \sum _{I,J}\varPsi ^\mathrm {eddy}(I,J,t, W,E), \end{aligned}$$

(11)

$$\begin{aligned} \varXi ^{\mathrm {noneddy}} (t, W, E)&= \sum _{I,J}\varPsi ^\mathrm {noneddy} (I,J,t,W,E). \end{aligned}$$

(12)

Because it is rather challenging to understand the meaning of the horizontal distribution of transport, we change it to the probability for each particle by dividing the total transport for each water mass. This is also customary practice in the visualized trajectory atlas (van Sebille et al. 2018).

$$\begin{aligned}&P^\mathrm {eddy}(I,J,t,W,E) \nonumber \\&\quad = \frac{\varPsi ^\mathrm {eddy}(I,J,t,W,E)}{\varXi ^\mathrm{{eddy}} (t, W, E) + \varXi ^\mathrm {noneddy} (t, W, E) }, \end{aligned}$$

(13)

$$\begin{aligned}&P^\mathrm {noneddy}(I,J,t,{{\mathcal {W}}},E) \nonumber \\&\quad = \frac{\varPsi ^\mathrm {noneddy}(I,J,t, W,E)}{\varXi ^\mathrm{{noneddy}} (t, W, E) + \varXi ^\mathrm {noneddy} (t, W, E) }. \end{aligned}$$

(14)

Then by averaging for t and E, we obtain the climatological probability distribution for each water mass for the eddy component ($\overline{P^\mathrm {eddy}(I,J,t, W,E)}^{t,E}$) and the non-eddy component ($\overline{P^\mathrm {noneddy}(I,J,t,W,E)}^{t,E}$).

We also calculate the time required to reach the 137$^{\circ }$ E section ($T^\mathrm {eddy}$, $T^\mathrm {noneddy}$) using the recorded elapsed time from the release position every five days $(\tau _{n}(t))$. By weighting the transport assigned to each particle, we obtain

$$\begin{aligned} T^\mathrm {eddy}(I,J,t,W,E)&= \frac{\tau _n(t)~\varPsi ^\mathrm {eddy}(I,J,t,W,E)}{\varXi ^\mathrm{{eddy}} (t, W, E) + \varXi ^\mathrm {noneddy} (t, W, E) }, \end{aligned}$$

(15)

$$\begin{aligned} T^\mathrm {noneddy}(I,J,t,W,E)&= \frac{\tau _n(t)~\varPsi ^\mathrm {noneddy} (I,J,t,W,E)}{\varXi ^\mathrm{{eddy}} (t, W, E) + \varXi ^\mathrm {noneddy} (t, W, E) }. \end{aligned}$$

(16)

By averaging for t and E, we obtain the average time required to reach the 137$^{\circ }$ E section for the eddy component ($\overline{T^\mathrm {eddy}(I,J,t,W,E)}^{t,E}$) and the non-eddy component ($\overline{T^\mathrm {noneddy}(I,J,t,W,E)}^{t,E}$).

For statistical testing between $\overline{T^\mathrm {eddy}}(i,j,W)$ and $\overline{T^\mathrm {noneddy}}(i,j,W)$, we first calculate the annual-mean members for ${T^\mathrm {eddy}}(i,j,W)$ and ${T^\mathrm {noneddy}}(i,j,W)$ between 1994 and 2017. We consider that the annual mean in each year between ${T^\mathrm {eddy}}(i,j,W)$ and ${T^\mathrm {noneddy}}(i,j,W)$ are related; therefore, we conduct dependent t-tests for the null hypothesis ${T^\mathrm {eddy}}(i,j,W) = {T^\mathrm {noneddy}}(i,j,W)$ for a 90% significance level. Considering the significant temporal variations in subduction in each year (Fig. 9), we set the degree of freedom as $2017 - 1994 + 1$.