Journal of Geophysical Research: Oceans
RESEARCH ARTICLE
10.1002/2017JC013359
Key Points:
 Oceanic dispersal of dissolved
radionuclides is investigated using ahigh-resolution, coupledcirculation-ofﬂine Eulerian transportmodel
 Submesoscale eddy-induced
ageostrophic secondary circulations(ASCs) notably promoted verticalmixing shortly after the nuclear
accident ASCs occur primarily because ofshear instability, with baroclinic
instability under late winterconditions when the initial dilution
occurred
Supporting Information:
 Data Set S1Correspondence to:
Y. Kamidaira,kamidaira.yuki@jaea.go.jp
Citation:
Kamidaira, Y., Uchiyama, Y.,
Kawamura, H., Kobayashi, T., &Furuno, A. (2018). Submesoscale
mixing on initial dilution ofradionuclides released from the
Fukushima Daiichi Nuclear Power
Plant. Journal of Geophysical Research:
Oceans, 123, 2808–2828. https://doi.
org/10.1002/2017JC013359
Received 15 AUG 2017
Accepted 16 MAR 2018
Accepted article online 26 MAR 2018
Published online 24 APR 2018

C 2018. American Geophysical Union.
V

All Rights Reserved.

KAMIDAIRA ET AL.

Submesoscale Mixing on Initial Dilution of Radionuclides
Released From the Fukushima Daiichi Nuclear Power Plant
Yuki Kamidaira1,2
Akiko Furuno1, Yusuke Uchiyama2,3, Hideyuki Kawamura1, Takuya Kobayashi1, and1

Nuclear Science and Engineering Center, Japan Atomic Energy Agency, Tokai, Japan, 2Department of Civil Engineering,
Kobe University, Kobe, Japan, 3Coastal and Estuarine Environment Department, Port and Airport Research Institute,
Yokosuka, Japan

Abstract This study developed a submesoscale eddy-resolving oceanic dispersal modeling system comprising a double-nested oceanic downscaling model and an ofﬂine oceanic radionuclide dispersion model.
This was used to investigate the inﬂuences of submesoscale coherent structures (SCSs) and associatedageostrophic secondary circulations (ASCs) on the three-dimensional (3-D) dispersal of dissolved cesium137 (137Cs) released from the Fukushima Daiichi Nuclear Power Plant (FNPP1). Extensive model-data comparison demonstrated that the innermost high-resolution model, with a lateral grid resolution of 1 km,
could successfully reproduce transient mesoscale oceanic structures, the Kuroshio path and stratiﬁcation,and spatiotemporal variations of 137Cs concentrations. Using an accompanying mesoscale eddy-resolving
model (grid resolution: 10 km) as a guide, we showed that submesoscale dynamics are important for
improved representation of both the eddy ﬁeld and the resultant 3-D dispersal of 137Cs, with the temporalvariability of surface 137Cs near the FNPP1 being equivalent to that in the coarse-resolution model. According to energy conversion and spectral analyses, SCSs and ASCs occur most intensively on the submesoscale,primarily because of shear instability. However, baroclinic instability serves as a secondary mechanism. SCSs
have prominent seasonality, reﬂected by intensiﬁcation in the colder months, which is when the FNPP1
accident occurred. Analysis of the vertical ﬂux of 137Cs was performed by decomposition of the variables
into eddy, mesoscale, and submesoscale components using frequency and wave number ﬁlters. It revealed
that 42.7% of the FNPP1-derived 137Cs was transported downward below the mixed layer by eddies withthe major contribution being from ASCs induced by submesoscale eddies.

1. Introduction
At 14:46 JST on 11 March 2011, a magnitude 9.0 earthquake occurred off the Paciﬁc coast of Tohoku, Japan,
generating catastrophic tsunamis that struck the Tohoku coast (Mori et al., 2011). Inundation and runup
height of the tsunamis reached 15 m at the Fukushima Daiichi Nuclear Power Plant (FNPP1). Subsequently,
the cooling system for the reactors was severely damaged, which led to three nuclear meltdowns andhydrogen-air explosions. Consequently, substantial amounts of radioactive substances were released accidentally from the FNPP1 to both the atmosphere and the ocean (Nuclear Emergency Response Headquarters, Government of Japan, 2011).
Several modeling studies on source term estimation, in relation to the accidental release, have been conducted because the acquisition of suitable measurements was extremely difﬁcult to accomplish shortly
after the accident. Among others, the Japan Atomic Energy Agency (JAEA) estimated the temporal variations of leaked radionuclide doses during the accident using atmospheric and oceanic transport, dispersion, and deposition models. Kawamura et al. (2011) reported that the total amount of oceanic direct
release of cesium-137 (137Cs; half-life 5 30.17 years, e.g., National Institute of Standards and Technology,
https://www.nist.gov/pml/radionuclide-half-life-measurements-data) from 21 March to 30 April 2011, wasestimated to be 4 PBq, which agrees well with the model-based estimate by Tsumune et al. (2012). Katata
et al. (2015) estimated the atmospheric release of 137Cs from the FNPP1 to be 14.5 PBq from 12 March to 1
May 2011, which is consistent with estimates reported in several other studies (e.g., Aoyama et al., 2016a;
Inomata et al., 2016; Tsubono et al., 2016).

2808

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Model assessments have also been undertaken to investigate the oceanic dispersal of the radionuclides
that leaked from the FNPP1 on and around the continental shelf off Fukushima, with primary focus on thesurface transport of dissolved 137Cs. For example, Miyazawa et al. (2012) conducted a numerical hindcast for
dissolved 137Cs, released directly from the FNPP1, using an Eulerian passive tracer transport model coupledwith a regional circulation model (horizontal resolution: 1/368, i.e., 3 km) based on the hydrodynamic kernel of the Japan Coastal Ocean Predictability Experiments (JCOPE2) system (Miyazawa et al., 2009). They
performed a series of experiments to determine the sensitivity of oceanic 137Cs dispersal to surface winds,
tides, and river discharges, and they highlighted that the role played by wind facilitates the meridional
extension of the distribution of surface 137Cs on the shelf. Estournel et al. (2012) modeled oceanic 137Cs dispersal by considering both direct release and atmospheric deposition using the SYMPHONIE model (Marsaleix et al., 2009, 2012) with variable horizontal grid sizes that increased linearly from 600 m to 5.5 km with
distance from the FNPP1. They identiﬁed that in addition to the importance of wind, riverine freshwater discharges promoted the offshore transport of 137Cs via surface-concentrated buoyant plumes. Masumoto
et al. (2012) conducted an intermodel comparison of oceanic 137Cs dispersal using ﬁve different regional
oceanic models (including the preceding two models) at horizontal resolutions ranging from 600 m to
3.0 km. The ﬁve models were found to show overall qualitative similarities, despite substantial differences inthe detail of their processes. Modeled surface 137Cs dispersal depends considerably on the lateral grid resolution of the model, where coarser resolution (e.g., grid spacing  1 km) models tend to overestimate thecross-shelf transport. In contrast, ﬁner-resolution models are more capable of reproducing the prevailingcoastally trapped alongshelf transport that acts to retard immediate offshore transport and subsequent
dilution via the synoptic Kuroshio-Oyashio current system. In summary, previous studies have identiﬁed
four key mechanisms that could have determined the initial oceanic dispersal of the FNPP1-derived dissolved 137Cs shortly after the accident occurred: (1) wind-driven, coastally trapped alongshelf currents; (2)buoyancy inputs from adjacent rivers; (3) a mesoscale clockwise lateral circulation that developed southwest of the FNPP1 in April 2011 associated with the migration of a warm-core eddy pinched off from theKuroshio Extension Jet, which interacted with (1); and (4) the Kuroshio-Oyashio system, which was responsible for subsequent offshore dispersal and dilution.
Compared with surface processes, much less attention has been paid to the subsurface three-dimensional(3-D) transport. The numerical hindcast by Tsumune et al. (2012, 2013) could be regarded as one of the mostsuccessful high-resolution (1 km lateral resolution) modeling studies because it has been used as a benchmark
for other models (e.g., Estournel et al., 2012; Kawamura et al., 2014). Its results are generally found consistent
with in situ 137Cs measurements, temporally, spatially, and in terms of vertical distributions, despite the deﬁciencies of their modeling strategy, e.g., the large nesting grid-size ratio (viz., the parent 10 km resolution isdownscaled directly to the child 1 km resolution), lack of river inputs, and the neglect of topography deeper
than 1,000 m. In contrast, the mesoscale eddy-resolving oceanic cesium-134 (134Cs; half-life 5 2.07 years) and137
Cs dispersal modeling (horizontal grid resolution: 1/108) conducted by Kawamura et al. (2014) underestimated subsurface cesium concentrations, particularly in intermediate (depth: 100–500 m) layers. These results
suggest that high-resolution models are required for satisfactory reproduction of both nearshore-conﬁned
alongshelf transport and subsurface 3-D dispersal of dissolved radionuclides in offshore areas.
A lateral grid resolution of 1 km is widely accepted as the threshold maximum model resolution for adequate representation of submesoscale processes (e.g., Capet et al., 2008a), which is crucial for both the generation of substantially enhanced vertical mixing and the resultant material dispersion in the upper ocean(e.g., Badin et al., 2011; Boccaletti et al., 2007; Callies et al., 2015; Kamidaira et al., 2017; Kunze et al., 2015).
The intensities of submesoscale eddies, fronts, and ﬁlaments, collectively referred to as submesoscale coherent structures (SCSs) in this paper, vary seasonally with remarkable ampliﬁcation in the colder season. Thisis because of mixed-layer deepening and surface cooling that precondition baroclinic instability and symmetric instability near the surface (e.g., Sasaki et al., 2014; Thomas & Taylor, 2010). Because the FNPP1 accident occurred in late boreal winter (March), it could be expected that the SCSs would have inﬂuenced thesubsequent 3-D mixing and dispersal of the dissolved radionuclides. Note that the approximate scale ranges
for SCSs reported in the literature are 0.1–10.0 km horizontally, 0.01–1.00 km vertically, and hours-days temporally (e.g., McWilliams, 2016).
Multiple-nested, very high-resolution (grid spacing: 75 m) downscaling ocean modeling studies (Romero
et al., 2013; Uchiyama et al., 2014) have demonstrated that SCSs promote anisotropic, enhanced alongshelf

KAMIDAIRA ET AL.

2809

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

dispersal of Eulerian passive tracer concentrations and buoyantly neutral Lagrangian particles with minimal
cross-shelf dispersal in innershelf and nearshore areas. Both observational and computational studies haveshown that SCSs are accompanied by energetic vertical motions, referred to as ageostrophic secondary circulations (ASCs) (e.g., Capet et al., 2008a; Levy et al., 2001; Mahadevan & Archer, 2000; Nagai et al., 2006,
2012; Perez et al., 2003; Pollard & Regier, 1992). The kinetic energies of both the mesoscale and the submesoscale are dominated by geostrophy (e.g., Sasaki et al., 2014). However, according to quasi-geostrophic
theory, vertical motions depend on the ageostrophic components of the ﬂows rather than on geostrophy.
The thermal wind relation also assumes geostrophy in the basic cross-frontal momentum balance, whichfails in representing the mixed-layer structure where submesoscale ﬁlaments or fronts develop. Gula et al.
(2014) proposed that the turbulent thermal wind balance, considering geostrophy and vertical momentum
mixing, adequately describes cross-frontal ﬂows in a mixed layer that incorporates overturning circulations(i.e., ASCs). Thus, SCSs comprise along-frontal (or along-ﬁlamentary) geostrophic ﬂows and cross-frontal (orcross-ﬁlamentary) ASCs (e.g., McWilliams, 2016), and the latter is important to the enhancement of verticalmaterial transport and mixing. ASCs are formed via strong downwelling from zones of surface convergence
due to dense (cold) fronts and ﬁlaments with peripheral upwelling (McWilliams et al., 2015). Mahadevan
and Tandon (2006) found that ASCs in the surface mixed layer occur through frontogenesis, which invokes
a temporal increase of the horizontal buoyancy gradient under mesoscale external strain mostly on theﬂanks of mesoscale eddies. ASCs are considered to interact with mesoscale vertical ﬂows in the surface
boundary layer and thus, they enhance vertical mixing. Therefore, given that Tsumune et al. (2012, 2013)
were not explicit about the inﬂuences of SCSs, we consider that submesoscale eddy-resolving 3-D modelingis important for successful 3-D dispersal realization, particularly in coastal shelf regions such as near the
FNPP1. However, the submesoscale inﬂuences on the 3-D dispersal of radionuclides, released by the FNPP1
accident, in the coastal shelf area off Fukushima have never been investigated in detail, either observationally or numerically.
In the present study, we developed a submesoscale eddy-resolving oceanic dispersal modeling system forthe coastal margin off Fukushima in order to conduct comparative hindcasts of the initial dispersal and subsequent dilution of dissolved 137Cs after the FNPP1 accident, paying particular attention to the inﬂuences of
SCSs. Our analysis is based on a double-nested ocean downscaling 3-D circulation modeling using theRegional Oceanic Modeling System (ROMS; Shchepetkin & McWilliams, 2005, 2008) embedded within an
assimilative oceanic reanalysis, Multivariate Ocean Variational Estimation for the Western North Paciﬁc(MOVE-WNP; Usui et al., 2006) product with lateral grid resolution of 1/108. Surface wind forcing was provided by assimilative atmospheric reanalysis products, the Grid Point Value of the Global Spectral Model(GPV-GSM; e.g., Roads, 2004) for the parent ROMS (L1) model, and the Mesoscale Model (GPV-MSM; e.g., Isoguchi et al., 2010) for the child ROMS (L2) model. The ROMS-L2 model, which constitutes our principal tool,has 1 km horizontal grid spacing that is considered sufﬁcient for full representation of SCSs. Oceanic 3-D
dispersal of the radioactive passive tracer released from the FNPP1 was simulated using the SEA-GEARNFDM (SGF) Eulerian advection-diffusion model, which was developed at JAEA and designed speciﬁcally for
oceanic radioactive tracers (Kawamura et al., 2014). Both the direct release and the atmospheric depositionof radionuclides (137Cs) from the FNPP1 are taken into account as point and area sources in the SGF model.
Note that the SGF model is a numerical ofﬂine Eulerian tracer-tracking model that exploits spatiotemporallyinterpolated 3-D ﬂow from a circulation model a posteriori. By taking advantage of this functionality, we ran
the SGF model to compute 137Cs dispersal with the SCS-resolving ROMS-L2 model and the under-resolving
MOVE-WNP model results used for comparative examination of submesoscale effects on 137Cs dispersal. Following the FNPP1 accident, oceanic dissolved 137Cs has been observed intensively (e.g., Aoyama et al.,2016a, 2016b; Buesseler et al., 2017; Kaeriyama et al., 2016; Kumamoto et al., 2017). Thus, we beneﬁted fromthe availability of abundant 137Cs data for quantitative validation of the coupled circulation-dispersal
model.
The remainder of this paper is organized as follows. The modeling framework used for the hindcast experiment for the 4 month period following the FNPP1 accident is described in section 2. Section 3 presents anextensive comparison between the model results, ﬁeld observations, and satellite altimetry data to validate
the model’s capability in reproducing oceanic ﬂow ﬁelds and radionuclide concentrations. Seasonal variabilities of the submesoscale phenomena and their mechanisms are examined with consideration of eddy
energy conversion and spectral analyses in section 4. Section 5 discusses the inﬂuences of ASCs associated

KAMIDAIRA ET AL.

2810

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

with submesoscale eddies in late winter to spring on vertical mixing of 137Cs. To conclude, the summaryand discussion are presented in section 6.

2. The Models
2.1. Double-Nested ROMS Model
We used the state-of-the-art ROMS to simulate the circulation off the northeastern Paciﬁc coast of Japan ina double-nested conﬁguration. The ROMS is a free-surface, terrain-following, Primitive-equation 3-D
regional hydrodynamic model that incorporates the hydrostatic assumption and the Boussinesq approximation. A versatile ﬁrst-order boundary turbulence model, the K-Proﬁle Parameterization (KPP) model (Largeet al., 1994), was employed to determine the mixed-layer depth (MLD), vertical eddy viscosity, and verticaleddy diffusivity, Ks. Figure 1 shows the hierarchy of the double-nested ROMS model domains (i.e., the intermediate, parent ROMS-L1 and the innermost, child ROMS-L2) embedded within the outermost MOVE-WNP
domain. The MOVE-WNP product is an operational numerical reanalysis/forecast product of the Japan Meteorological Agency (JMA) for the northwestern Paciﬁc. It is based on the Meteorological Research InstituteCommunity Ocean Model (MRI.COM) (Tsujino et al., 2010) and it is assimilated with considerable amounts ofsatellite and in situ data, including Argo proﬁling ﬂoat data, using a 3-D variational method. The MOVE-WNP
output comprises daily averaged data of sea surface height (SSH), temperature (T), and salinity (S), andmeridional and zonal horizontal current velocities. We adapted a one-way ofﬂine nesting technique,described in Mason et al. (2010) and in Uchiyama et al. (2014, 2017a, 2017b), to successively reduce the horizontal grid size from 10 km (MOVE-WNP), to 3 km (ROMS-L1), and then to 1 km (ROMS-L2). The extent of
the ROMS-L1 (ROMS-L2) domain was 1,200 3 960 km (704 3 512 km) horizontally with uniform square grid
spacing of 3 km (1 km) and 32 vertically stretched layers with reﬁnement near the surface, to encompass
the Kuroshio Extension region off Japan. The monthly climatology of freshwater inputs from 14 (8) majorrivers in ROMS-L1 (ROMS-L2) was taken into account.
The open boundary and initial conditions of SSH, T, S, and zonal and meridional horizontal current velocities
for the ROMS-L1 domain comprised the daily averaged MOVE-WNP data. Model topography was takenfrom the Shuttle Radar Topography Mission product (SRTM30_Plus) (Becker et al., 2009; Rodriguez et al.,2005), which covers the global ocean at 30 geographic arcseconds (1 km), supplemented with the
J-EGG500 (http://www.jodc.go.jp/data_set/jodc/jegg_intro.html) with 500 m resolution for coastal regions.
The GPV-GSM operational atmospheric reanalysis/forecast product of the JMA (horizontal resolution: 0.28 3
0.258) was applied as daily surface wind stresses in the L1 model. Surface heat, freshwater, and radiation
ﬂuxes were given by the monthly climatology of the Comprehensive Ocean-Atmosphere Data Set (COADS)

Figure 1. Double-nested ROMS model domains with bathymetry in color for the parent ROMS-L1 (dotted box) and child
ROMS-L2 (black solid box) models, embedded within the MOVE-WNP domain (not shown, much wider than the displayedarea). Red star indicates the location of the Fukushima Daiichi Nuclear Power Plant (FNPP1). Bold black line identiﬁes thetransect for the cross-sectional plots in Figure 4. Gray box is the area used for the spectral analyses in Figure 12. Thin red
rectangle corresponds to a zoomed-in region of the SGF-L2 domain (right plot) with the locations of the near-surface
monitoring stations (magenta circles) used in Figures 10a–10c. Model domain for SGF-M extends from 108N–608N, 1208E–
1208W (not shown).

KAMIDAIRA ET AL.

2811

Journal of Geophysical Research: Oceans

10.1002/2017JC013359(Woodruff et al., 1987). The 20 day averaged SST and sea surface salinity (SSS) data, obtained from the
MOVE-WNP product, were applied for surface ﬂux restoration to avoid long-term drift and bias caused byerroneous surface ﬂuxes from the monthly COADS climatology. The monthly climatology of the major riverdischarges, compiled by the Japan River Association (http://www.japanriver.or.jp/publish/book/nenpyou_
dvd.htm, in Japanese), were applied for all major rivers that discharge in the L1 domain. A four-dimensional
TS nudging (a.k.a. robust diagnostic) (e.g., Marchesiello et al., 2003; Uchiyama et al., 2017a, 2017b; Wrightet al., 2006) with a weak nudging inverse time scale of 1/20 d21 was introduced toward the 10 day averaged MOVE-WNP T and S ﬁelds to improve the reproducibility of the Kuroshio path in the L1 domain. The
L1 model ran from 1 January 2011 to 1 September 2011.
The innermost ROMS-L2 domain was forced at the lateral open boundaries and it was initialized using the
spatiotemporally interpolated 12 hourly averaged L1 output. The 3 hourly surface wind stresses from the
JMA GPV-MSM reanalysis were exploited, encompassing the entire L2 domain with horizontal resolution of
0.058 3 0.06258. Similar to the L1 model, the SST and SSS were restored toward the MOVE-WNP SST and
SSS as the surface ﬂux correction. The other numerical conditions were the same as the L1 model, although
we did not include TS nudging in the L2 model in order to avoid interference with the spontaneous development and decay of intrinsic variability associated with eddies. Surface ﬂuxes other than wind stresses inthe L2 domain were also provided by the monthly climatology of the COADS. This choice has been
accepted widely in both submesoscale eddy-resolving and submesoscale eddy-permitting model studies
with ROMS (e.g., Capet et al., 2008a; Gula et al., 2014; Uchiyama et al., 2017a), whereas more realistic MLD
responses might be attained with more elaborate heat ﬂux estimates, such as the bulk formula (e.g., Sasaki
et al., 2014). The L2 computational period ran from 11 January 2011 to 1 September 2011. The ﬁrst 60 days
were considered the model spin-up period during which initial adjustment occurred; thus, these data were
excluded from subsequent analyses. All statistical analyses conducted in the present study utilized themodel results for the same 4 month period (1 March to 30 June 2011) unless stated otherwise. It is worth
noting that we did not include tidal forcing because in situ monitoring around the FNPP1 has shown the
tidal currents to be weak (Kubota et al., 1981) and insigniﬁcant in terms of dispersal (Miyazawa et al., 2012).
Furthermore, dispersal and mixing on continental shelves and nearshore areas are affected predominantly
by submesoscale stirring rather than by tides (e.g., Romero et al., 2013).
2.2. SEA-GEARN-FDM Model
Oceanic dispersal of the dissolved 137Cs released from the FNPP1 was investigated using the SGF model(Kawamura et al., 2014). The SGF model is an ofﬂine Eulerian 3-D advection-dispersion model of arbitrarypassive tracers. The source term estimation representing the direct release and the atmospheric depositionof the 137Cs during the accident was implemented as described below. Table 1 lists the numerical conﬁguration of the SGF model. To evaluate the differences in 137Cs dispersal due to the spatial resolution of the ﬂow
ﬁeld, we conducted two SGF runs (i.e., SGF-L2 and SGF-M) by imposing the daily averaged ROMS-L2 (lateral

Table 1
Numerical Conﬁgurations of the Two SEA-GEARN-FDM Model Experiments
ModelsComputational period
Time step
Model domain
Horizontal
grid cells
Horizontal
resolution
Vertical level(z-coordinate)
Horizontal diffusivity coefﬁcient
Vertical diffusivity coefﬁcient
Duration of the oceanic release
Duration of the atmospheric deposition

KAMIDAIRA ET AL.

SGF-L2SGF-M
12 Mar 2011 to 30 Jun 2011, UTC
120 s
360 s
140.58E–1458E,
1208E–1208W,
35.68N–39.28N
108N–608N
451 3 361
1,200 3 500
0.018

0.18

79 levels

54 levels

5 m2 s21
50 m2 s21
KPP
5 3 1024 m2 s21
26 Mar 2011 to 30 Jun 2011, UTC
12 Mar 2011 to 31 May 2011, UTC2812

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 2. Time series plots of the imposed release rates of 137Cs from 12 March to 30 June 2011. Blue line is the oceanic
direct release estimated by Kawamura et al. (2011); red line is the atmospheric release estimated by Katata et al. (2015).

resolution: 1 km) and MOVE-WNP (10 km) outcomes as input data. The ROMS-L2 (MOVE-WNP) ﬂow ﬁelds
were interpolated linearly at the 3-D grid points of the SGF-L2 (SGF-M) model with horizontal resolutions of
0.018  1 km (0.18  11 km) in the latitudinal and longitudinal directions, and 79 (54) vertical z-levels withnear-surface reﬁnement ranging from 2 m (1 m) near the surface to 500 m (900 m) near the bed. The daily
averaged input ﬁeld was further interpolated linearly in time to obtain any arbitrary intermediate input
data. The KPP-derived tracer vertical eddy diffusivity, Ks, was exploited for the vertical eddy diffusivity of the
dissolved 137Cs concentration in the SGF-L2 model. In contrast, a constant value of Ks 5 5 3 1024 m2 s21
was assumed in the SGF-M model, because the MOVE-WNP product does not include its own Ks. This valueof Ks was considered optimal for the present study, according to our preliminary sensitivity test with the
SGF-M. We observed that a smaller value of Ks led to degraded 137Cs reproducibility (not shown).
The 137Cs point source introduced as the direct release from the FNPP1 into the ocean was based on in situ
137
Cs measurements (Kawamura et al., 2011). We also incorporated atmospheric deposition at the sea surface as the area source of 137Cs based on the Worldwide version of the System for Prediction of Environmental Emergency Dose Information version II (WSPEEDI-II) developed at JAEA (Terada & Chino, 2008). The
amount of 137Cs released from 12 March to 30 April 2011, estimated by Katata et al. (2015), was applied as
the atmospheric release rate for WSPEEDI-II (Figure 2). The WSPEEDI-II model domain encompassed the
SGF-L2 (SGF-M) domain with horizontal resolution of 2 km (10 km) in the latitudinal and longitudinal directions. The 3 hourly (daily) deposition amount of airborne radioactive cesium predicted by WSPEEDI-II was
deposited at the sea surface in SGF-L2 (SGF-M) from 12 March to 31 May 2011. Only dissolved 137Cs was
taken into account in the present study because it dominated the initial dose inventory, although the SGF
model could account for suspended matter such as mineral and biogenic particles that exchange 137Cs withwater by absorption and desorption. Note that a truncated half-life of 30 years was used for 137Cs in the
SGF model runs because the simulations lasted for only 4 months.

3. 3-D Dynamics and 137Cs Dispersal
3.1. Mean 3-D Dynamics
We compared the double-nested ROMS model results with both satellite data and ﬁeld observations to validate the reproducibility of the surface dynamics, stratiﬁcation, and eddy intensity. The AVISO delayed-time
mapped absolute dynamic topography from multiple satellite altimetry data (viz., gridded daily composite
at 1/48 resolution) (e.g., Le Traon et al., 1998) was used for the mean surface currents, including the Kuroshio. The time-averaged surface current velocity magnitude from the ROMS-L2 model, overlaid on thatfrom the ROMS-L1 model, shows reasonable quantitative agreement with the satellite observations in terms

KAMIDAIRA ET AL.

2813

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

of both magnitude and spatial patterns, including the meandering Kuroshio path (Figures 3a and 3b). The
SSH variance is shown in Figures 3c and 3d as a proxy for the temporal variability of the synoptic and mesoscale signals. An overall agreement is exhibited between the modeled SSH variance and that of AVISO. For
instance, both variances are intensiﬁed considerably around the Kuroshio and its Extension Jet, corresponding to the meandering path of the Kuroshio.
The modeled stratiﬁcation was compared with the XBT-DBT observations acquired by the JMA’s survey
ship, the R/V Ryofu Maru (http://www.data.jma.go.jp/gmd/kaiyou/db/vessel_obs/data-report/html/ship/cruisedata_e.php?id5RF1107), along the transect off Fukushima that is represented by the bold black line in
Figure 1. Snapshots of temperature in the vertical cross section on 6 June 2011, from the L2 model and the
JMA observations are consistent, indicating a tilted thermocline around the Kuroshio from 36.28N to 36.88N(Figures 4a and 4b). The time-averaged Kuroshio main body exists from 36.28N to 37.28N with mean streamwise velocities of <0.3 m s21, even at the depth of 400 m, consistent with previous studies (e.g., Qiu, 2001)(Figure 4c). Hence, the deep thermocline at around 36.28N reﬂects the subtropical warm water brought by
the Kuroshio. The temporal standard deviation (RMS) of vertical velocity w is magniﬁed around the Kuroshio
main body and near the base of the time-averaged mixed layer on the northern side of the Kuroshio(Figure 4c). In turn, the time-averaged Ks estimated by the KPP model is large mostly on the northern sideof the Kuroshio at O (1023 m2 s21) in the upper ocean, which is three orders of magnitude larger than the
background diffusivity of O (1026 m2 s21) in the ocean interior (Figure 4d). These results demonstrate thatthe vertical advective transport and vertical eddy diffusion through Ks are concurrently responsible for vertical mixing in the upper ocean, whereas the former is anticipated to prevail in the vertical transport into theocean interior below the mixed layer and near the Kuroshio. Apparently, the constant value of Ks used in
the SGF model (section 2.2) is ad hoc, which introduces uncertainty in the dispersal modeling even though
our choice of Ks provides reasonable reproducibility for 137Cs (section 3.3). The introduction of an ofﬂine turbulent closure model could provide a more plausible estimation of Ks for use in the SGF model. However,this is beyond the scope of the present study and we leave this issue for future work.

Figure 3. Plan view plots of (a and b) time-averaged surface velocity magnitude, and (c and d) sea surface height (SSH)
variances. Left plots (a and c) are the ROMS-L2 result superposed on the L1 result, while right plots (b and d) are based on
AVISO satellite altimetry.

KAMIDAIRA ET AL.

2814

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 4. Cross-sectional plots of (a and b) snapshots of temperature from the ROMS-L2 model and JMA observation,(c) mean streamwise velocity at intervals of 0.15 m s21 normal to the transect (red contours) and temporal standard deviation (RMS) of vertical velocity, w (color) from the ROMS-L2, and (d) time-averaged vertical eddy diffusivity, Ks, estimatedby the KPP model, along the transect shown by the bold black line in Figure 1. White curves in the bottom plots are timeaveraged mixed-layer depths (MLDs) estimated by the KPP model.

3.2. Eddies
To validate the time-dependent reproducibility of the model, we examined the surface eddy kinetic energy(EKE):
1
Ke 5 u0 2 1v 02 ;(1)
2
where (u, v) is the horizontal (zonal and meridional) velocity, and the overbar represents an ensembleaveraging operator. The variables assigned with the prime are the ﬂuctuating eddy components obtainedby removing the seasonal variations with a low-pass Butterworth ﬁlter in the frequency domain with a 90
day cutoff period. Figures 5a and 5b indicate that Ke from the ROMS and MOVE-WNP are similar in theiroverall spatial patterns, where intensiﬁcation of Ke occurs around the Kuroshio path; however, the magnitudes are somewhat different. In the ROMS, Ke is increased because the lateral grid resolution is increased
for better resolution of the buoyancy gradients and the associated velocity shear, which results in higher
EKE, particularly around the Kuroshio front.
The differences in the Ke magnitude are clearly illustrated in snapshots of the daily averaged surface relative
vorticity, normalized by the local Coriolis frequency f, as contrasting emergence of mesoscale and submesoscale eddies on 2 April 2011, from the ROMS and MOVE-WNP models (Figures 5c and 5d). The two models
commonly yield positive and negative vortex streamers on the mesoscale at O (100 km), corresponding to
the Kuroshio. However, submesoscale eddies at O (1–10 km) develop much more prominently in ROMS-L2
than MOVE-WNP, which shows that the higher-resolution L2 model generates much ﬁner submesoscaleeddies. The eddies and associated SCSs in the L2 model have greater magnitudes, are distributed more
extensively, and are more diversiﬁed and strengthened than those in MOVE-WNP, all of which coincide withthe resultant Ke distributions (Figures 5a and 5b). As shown in Figure 6, the SCSs in the L2 model are accompanied by intensive vertical motions, as shown in the cross-sectional snapshot of the daily averaged verticalvelocity w along the transect deﬁned offshore from the FNPP1 indicated by the bold blue line in Figure 5b.
In turn, vertical ﬂows are extremely small in the MOVE-WNP model. The intensiﬁed w is produced in the surface mixed layer as ASCs due to SCSs (e.g., Capet et al., 2008a). The SCSs and ASCs develop more intensively
in the colder seasons from winter to spring with the MLD deepening (e.g., Sasaki et al., 2014), as is examined

KAMIDAIRA ET AL.

2815

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 5. (top) surface eddy kinetic energy (EKE), Ke, from (a) the ROMS-L2 on L1, and (b) MOVE-WNP reanalysis. (bottom)
Daily averaged spatial distributions of surface vertical relative vorticity normalized by planetary vorticity (dimensionless)on 2 April 2011 from (c) the ROMS-L2 on L1, and (d) MOVE-WNP models. Blue thick line in Figure 5b identiﬁes the transect
used for the cross-sectional plots in Figure 6.

in detail in section 4.3. In summary, the present double-nested ROMS model is shown to have satisfactory
reproducibility in terms of the synoptic and mesoscale dynamics, including the Kuroshio with the prominent ampliﬁcation of the submesoscale variability overlaid.
3.3. Dissolved Radionuclides
We veriﬁed the modeled 137Cs concentrations against a variety of observed data based on in situ water
sampling. Because the present study focuses on the initial dispersal and the dilution of 137Cs after the
FNPP1 accident, we used the ‘‘best available’’ offshore 137Cs dataset for June 2011 (Buesseler et al., 2012) for
the 3-D spatial distribution. To consider the temporal variation, we used nearshore surface 137Cs data frommid-March to the end of June 2011, collected intensively near the FNPP1 by the Tokyo Electric Power Company (TEPCO) (http://www.tepco.co.jp/en/nu/fukushima-np/f1/index2-e.html). The horizontal 137Cs

Figure 6. Cross-sectional plots of the daily averaged vertical upward velocity, w, from the (left) ROMS-L2 and (right)
MOVE-WNP models. Black line in the left plot is the mixed-layer depth (MLD) estimated by the KPP model. The corresponding transect is shown by the blue line in Figure 5b.

KAMIDAIRA ET AL.

2816

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 7. Modeled 137Cs concentrations (color) from the (top plots) SGF-L2 and (bottom plots) SGF-M models at (a and d)
the surface, (b and e) 100 m depth, and (c and f) 200 m depth. Colored circles indicate the instantaneous in situ concentrations (Buesseler et al., 2012) using the same color scale.

distributions are compared at the surface and at depths of 100 and 200 m for the two model results (i.e.,
SGF-L2 and SGF-M) with the measured data superposed as colored circles in Figure 7. Here we show only
the measured 137Cs considered derived from the FNPP1, which can be determined from the detection levelof another isotope, 134Cs. The number of available observed data in Buesseler et al. (2012) was only 87.
Overall, the two models reproduced the observations reasonably well, presenting high concentrations near
the shore with alongshelf spread and downward decreasing trends; however, they are not perfect. The
SGF-L2 generally outperforms SGF-M with an improvement of the model skill score (Wilmott, 1981) from
0.2917 to 0.3494 (n.b., model skill score is 1 for perfect agreement and 0 for complete disagreement). Even
though the bulk of the hydrodynamics evaluated by ROMS-L2 and MOVE-WNP are closely related, as discussed in sections 3.1 and 3.2, the resultant 137Cs levels are distributed quite differently. The surface 137Csfrom SGF-L2 is more comparable with the observed data, especially near the coast and the Kuroshio path
(Figure 7a). In contrast, the offshore transport manifest by SGF-M is more prominent, and is subsequently
entrained by the Kuroshio. The slight deviation of the predicted location of the Kuroshio, i.e., a northward
shift of about 50 km (0.58) relative to that of the SGF-L2, promotes intense dilution in the coastal area off
the FNPP1 and remarkable underestimation south of 36.58N (Figure 7b). Similarly, the subsurface 137Csfrom SGF-L2 shows reasonable agreement (Figures 7b and 7c), whereas SGF-M yields unrealistically highlevels of 137Cs near the coast, in particular at 100 m depth (Figures 7e and 7f). According to an observational
study (Aoyama et al., 2012), an anticyclonic mesoscale eddy centered at 36.58N, 141.48E in mid-May 2011
might have prevented southward alongshelf transport of 137Cs. Consistent with Aoyama et al. (2012),
SGF-L2 manifests southward coastally trapped alongshelf transport of 137Cs after the end of May 2011,
when the anticyclonic mesoscale eddy was suppressed by a tropical cyclone (not shown).
A similar conclusion can be drawn from the vertical cross-sectional views of 137Cs along the offshore transect from the FNPP1 (blue line in Figure 5b) in Figure 8. Although both models reproduce the 3-D 137Cs distribution moderately well for June 2011 (Figures 8a and 8b), the SGF-M model overestimates the downward
137
Cs transport. This results in higher concentrations at depth, partly because of our choice of the vertical
diffusivity Ks. The SGF-M also indicates an almost uniform pattern of 137Cs dilution in the offshore and

KAMIDAIRA ET AL.

2817

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 8. Cross-sectional plots of 137Cs concentrations from the (left) SGF-L2, and (right) SGF-M. (top) Time-averaged
137
Cs from 9 to 15 June, and (bottom) the daily averages for 2 April. Colored circles in top plots indicate instantaneous insitu concentrations measured in June 2011 (Buesseler et al., 2012). Black curve in Figure 8c is the mixed-layer depth estimated by the KPP model. The corresponding transect is shown by the blue line in Figure 5b.

downward directions due to intense lateral (isopycnal) advection and large vertical (diapycnal) diffusion. In
contrast, the SGF-L2 has greater patchiness, as evidenced by the near-surface peak concentrations at 1428E
and 1438E. These results suggest the coarse-resolution SGF-M model is rather ‘‘diffusive’’ with transport thatis spatially more homogeneous because of the poor representation ofw, which is compensated by strong vertical diffusion via the imposed
Ks. The scatter diagram of modeled versus measured 137Cs, shown in
Figure 9, illustrates that both the SGF-L2 and the SGF-M slightly
underestimate the observed concentrations. The results of a supplemental numerical experiment without atmospheric release suggested
that this underestimation could be caused by errors in the direct oceanic release rate estimation (not shown), in conjunction with the omission of the outgoing 137Cs returning to the L2 domain in SGF-L2,
rather than possible errors in the atmospheric source term estimation.
However, a large proportion of the modeled 137Cs from SGF-L2 iswithin the range of a factor of 10, indicated by the dotted lines relative to the measured 137Cs, with an apparent improvement for higher
concentrations of 137Cs > 1 Bq m23, compared with SGF-M, which
underestimates the measured 137Cs considerably. The improvementwith SGF-L2 is quantiﬁed by the rate of the modeled 137Cs within therange of a factor of 10 of the measured 137Cs. This rate computedwith all 87 data improved from 67.8% (SGF-M) to 71.3% (SGF-L2),
while that at depths ⭌100 m (with only 19 data) improved from
68.4% (SGF-M) to 73.7% (SGF-L2).
Figure 9. Scatter diagram of the 137Cs concentrations from the SGF-M (whitecircles) and SGF-L2 (red circles) models against in situ sampled data (Buesseler
et al., 2012). Error bars are quite small as they are within the size of the symbols,
and thus not shown here. The solid line represents the 1:1 correlation and it isbounded by two dashed lines indicating the deviations with factors of 1/10
and 10.

KAMIDAIRA ET AL.

These results demonstrate that the high-resolution SGF-L2 outperforms the coarse-resolution SGF-M in the overall reproducibility of137
Cs. A plausible reason for a more successful hindcast using SGF-L2
is the usage of a ﬁner lateral grid resolution, which allows more realistic representation of SCSs and ASCs (sections 3.1 and 3.2). We might

2818

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

be able to improve the SGF-M model by tweaking the background Ks and avoiding laborious highresolution modeling (see also section 3.1). To investigate the inﬂuence of SCSs, we examined the crosssectional views of 137Cs in April 2011, i.e., when SCSs are most pronounced. As shown in Figures 8c and 8d,the SGF-L2 results in more intensiﬁed levels of 137Cs near the surface. These extend offshore from the coastwith notable cross-shelf variability at cross-shelf wavelengths of 5–20 km, accompanied by downward penetration beyond the base of the mixed layer via intense ASCs. In contrast, the SGF-M yields a diffusive distribution of 137Cs with high concentrations remaining near the coast. Hence, it is conﬁrmed that SCSs and
ACSs have considerable effect on the 3D dispersal of 137Cs, particularly in the colder months (i.e., March andApril 2011) that correspond to the initial period shortly after the FNPP1 accident. Obviously, the degraded
estimates using SGF-M can be attributed to the coarse-resolution MOVE-WNP, which does not accountproperly for the SCSs but instead requires a tweak of Ks that introduces large uncertainty. However, the
SGF-L2 does have a downside in that some areas are assigned zero 137Cs near the northern and southern
open boundaries (Figure 7). This is because we did not consider outgoing 137Cs returning to the L2 domain.
However, this deﬁciency could be improved by introducing a nested SGF conﬁguration.

3.4. Temporal Variability
Figures 10a–10c compares temporal variations of near-surface 137Cs from the two models and measurements conducted by TEPCO (http://www.tepco.co.jp/en/nu/fukushima-np/f1/index2-e.html) at three stations near the FNPP1, as indicated in the inset in Figure 1. In general, both models provide reasonableagreement with the measurements. Model skill scores between the model results and the measurements
conducted by TEPCO, including those at the three stations, from 13 March to 30 June 2011, are improved
considerably from 0.1244 (SGF-M) to 0.6930 (SGF-L2). The SGF-L2 produced a greater number of apparent
temporal ﬂuctuations for periods of several days, most likely coinciding with eddying ﬂows that promoteintense downward transport (Figures 6 and 8). Temporal variations of the modeled 137Cs inventory in the
SGF-L2 model domain were then estimated using the volume integral of 137Cs, normalized by the total
amount of FNPP1-derived 137Cs for the 4 month period (Figure 10d). Because the inventory integrated fromthe surface down to the depth of 100 m is about the same as that down to the depth of 1000 m, thereleased 137Cs exists mostly in the upper ocean. The SGF-L2 result indicates that the surface mixed layer
stores 53% of the released 137Cs in the ﬁrst several weeks. The 137Cs inventory evaluated from the SGF-M
model is generally smaller than that from SGF-L2 until the end of April. This is because the SGF-M gives rise
to more prominent offshore transport, which is subsequently entrained by the Kuroshio (section 3.2). However, after 30 April, the initial dilution is more evident in SGF-L2 because the SGF-M model results retain
high 137Cs concentrations near the coast, whereas SGF-L2 might overestimate the dilution, partially becauseof the lack of reentering 137Cs. The area-averaged KPP-estimated MLD represents a nearly monotonic
decrease with time, where late winter conditions are attributed to a deeper MLD (Figure 10e). Figure 10f
shows a time series of the surface Ke (EKE), volume-averaged laterally over the SGF-L2 domain and from thesurface to a depth of 400 m vertically. The Ke ﬂuctuates similarly in the two circulation models with overall
intensiﬁcation in the high-resolution ROMS-L2 model compared with MOVE-WNP, particularly in the ﬁrst
month and in June. Because MOVE-WNP is exploited to force the ROMS-L1 and L2 models, the intensiﬁed Ke
is generated intrinsically because of the grid reﬁnement.
These results suggest that near-surface model-data comparison would be insufﬁcient for model validationand for source term estimation. This is because the SGF-L2 result indicates that signiﬁcant vertical transportof 137Cs due to eddies and associated vertical ﬂows are anticipated to alter markedly the surface distributionof 137Cs. Therefore, a higher-resolution model is needed for a more accurate representation of the verticalmixing of 137Cs, which is sensitive to eddying and transient ambient ﬂows that are largely affected by theseasonally varying MLD. In summary, we conﬁrmed that the double-nested ROMS coupled with the SGF(viz., SGF-L2) model was satisfactory in reproducing the 3D dispersal and initial dilution of dissolved 137Cs
originating from the FNPP1 accident. The open boundary artifact causing the long-term dilution of 137Cs in
the SGF-L2 could be improved if desired. We should emphasize that the SGF-M model is better suited thanSGF-L2 to basin-scale and long-term analyses, whereas the high-resolution model deﬁnitely outperforms
the former near the center of the accident and in coastal shelf regions, producing higher accuracy in the 3D representations of short-term dispersal.

KAMIDAIRA ET AL.

2819

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 10. Time series plots of various quantities. Near-surface 137Cs concentrations from the SGF-M (black) and SGF-L2(red) models, and measured data (white circles) around the FNPP1 at (a) T1, (b) T2, and (c) T3. The uncertainties of theobservations by TEPCO are not reported. Locations of the measurements are shown in Figure 1. (d) Inventory (the total
amount) of 137Cs in the SGF-L2 model domain, integrated vertically from the surface down to the labeled depths(z 5 2100 m, 21,000 m, and the base of the surface mixed layer denoted by MLD) normalized by the total 137Cs load
derived from the FNPP1 as the sum of the oceanic direct release and atmospheric deposition. (e) Spatial MLD estimatedby the KPP model of the ROMS-L2. (f) Volume-averaged surface (z > 2400 m) eddy kinetic energy, EKE (Ke), for the SGF-L2
domain.

4. Submesoscale Dynamics
4.1. Extraction of Submesoscale
SCSs in the North Paciﬁc are known to have strong seasonality, with much greater energetic variability incolder seasons than in warmer seasons, and mixed-layer deepening (Figure 10e) plays an important role inthe development of SCSs and ASCs (e.g., Sasaki et al., 2014). To separate submesoscale variability from theeddy component of any modeled prognostic ﬁeld, we introduced a horizontal spatial ﬁlter. This incorporated a convolution of the 2-D Gaussian smoothing kernel function, computed with a diameter of 20 km(viz., discrete 21 3 21 grid points for the ROMS-L2 and SGF-L2 models) relative to each model grid point
and a standard deviation of 4.5 km. This means that length scales shorter than 20 km constitute our deﬁnition of the submesoscale structures to be isolated. The variable V is thus decomposed as:

KAMIDAIRA ET AL.

2820

Journal of Geophysical Research: Oceans
V5V 1V 0 ; V 0 5Vm0 1Vs0 ;10.1002/2017JC013359(2)

where V and V’ are the seasonal and eddy components, respectively, as deﬁned in section 3.1, and subscripts m and s represent the mesoscale and submesoscale components, respectively, decomposed by the
2-D Gaussian ﬁlter.
4.2. Seasonal Variability
The monthly averaged surface submesoscale EKE from ROMS-L2 (Figures 11a–11c) shows that the EKE isgreatest in April 2011, not only around the mean Kuroshio path indicated by the blue contours for velocity
magnitude, but also in an extensive area to the north of the Kuroshio that includes the coastal area. The
associated vertical velocity variance w02 at the MLD (i.e., base of the surface mixed layer, Figures 11d–11f)
varies consistently with the EKE, where w02 is also largest in April, which is attributed to ASCs associatedwith SCSs. We then estimated the azimuthally averaged wave number spectra of the surface horizontal
velocity (i.e., kinetic energy spectra; Figure 12a) and the vertical velocity at the MLD (Figure 12b) for the subdomain deﬁned by the gray square in Figure 1 for each of the 3 months from April to June 2011. The spectral computations were performed ofﬂine using daily averaged values from ROMS-L2, with subtracting the
spatial mean values over the area used for the spectral analyses (the gray square in Figure 1) to deﬁne ﬂuctuating components. Two-dimensional wave number spectra were determined daily by a 2-D FFT with theazimuthal averaging operation posteriori, and then temporally averaged for each month. The maximum
resolved wavelength is 450 km, while the Nyquist wavelength is 2 km.
The surface kinetic energy spectra have similar shapes for the 3 months; however, the spectral slope is gentle in April with temporal steepening occurring until June (Figure 12a). In April, the spectral power is highest
at higher wave numbers and lowest at lower wave numbers. This indicates that kinetic energy in the submesoscale range is higher in April than in other months. The spectral slope is close to k23 rather than k22
for higher wave numbers, where the latter has been observed in fully developed submesoscale turbulent
regimes (e.g., Capet et al., 2008a). Figure 12b shows more striking submesoscale dominance in the upper
ocean. Similar to Figure 12a, the w spectrum in April is notably enhanced for a notably extensive wave number range, peaking at wavelengths of about 10–20 km (viz., submesoscale), even with the deepest mixed
layer among the ﬁrst 3 months after the accident (Figure 10e). These results are consistent with the spatial
distributions of the submesoscale EKE and w02 shown in Figure 11. Therefore, submesoscale variability was

Figure 11. (top) Monthly averaged surface submesoscale EKE from the ROMS-L2 for (a) April, (b) May, and (c) June 2011.
(bottom) Submesoscale vertical velocity variance w02 at the mixed-layer depth (MLD) from the ROMS-L2 for (d) April,(e) May, and (f) June 2011. Contours represent monthly averaged surface velocity magnitude at intervals of 0.25 m s21,
depicting the monthly averaged Kuroshio path.

KAMIDAIRA ET AL.

2821

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 12. Wave number spectra of (a) kinetic energy at the surface, (b) vertical velocity at the mixed-layer depth (MLD),(c) barotropic conversion rate, KmKe, averaged vertically over the surface mixed layer, and (d) baroclinic conversion rate,
PeKe, averaged vertically over the surface mixed layer. Curves are for each of the 3 months from April to June 2011.
Straight lines indicate the spectral slopes of k25/3 (black), k22 (dashed), and k23 (gray), where k is the wave number.

still evident with reasonably intense ASCs in April 2011, although it was slightly weakened from the fully
developed stage that perhaps occurred in January and February. These results also conﬁrm that the combination of the COADS climatological heat, freshwater, and radiative surface ﬂuxes with the realistic GPV-MSM
wind stresses may be able to reproduce the development of SCSs and ASCs to some extent (see also section
2.1). In addition, the w spectra have a secondary peak near the wavelength of 100–200 km (viz., mesoscale),which conﬁrms that both mesoscale and submesoscale variability concurrently inﬂuence vertical mixing in
the upper ocean.
4.3. Eddy Generation Mechanisms
To investigate a possible cause of the instability and the eddy-mean interaction mechanisms, we examinedthe energy conversion rates in the Ke (EKE) conservation equation (e.g., Dong et al., 2007; Kamidaira et al.,2017; Klein et al., 2008; Marchesiello et al., 2003;). In this study, we focus on the two most important quantities: the release of eddy available potential energy (Pe) to feed Ke (PeKe), and the work done by the Reynolds
stress on the mean velocity shear (KmKe), which represents the energy exchange between the mean kinetic
energy (Km) and Ke (EKE):


@u
@u
@u
@v
@v
@v
Km Ke 52 u0 u0
1u0 v 0
1u0 w 01v 0 u0
1v 0 v 0
1v 0 w 0
;
(3)
@x
@y
@z
@x
@y
@z
Pe Ke 52

g 0 0
qw ;
q0

(4)

where x, y, and z are the horizontal and vertical coordinates, q is the density of sea water, q0 5 1027.5 kg m23
is the Boussinesq reference density, and g is gravitational acceleration. If the barotropic conversion (or shear
production) rate KmKe > 0, then barotropic and vertical shear instability occurs, destabilizing the mean ﬂow
to generate eddies; thus, Ke increases. In contrast, when the baroclinic conversion rate PeKe > 0, baroclinic

KAMIDAIRA ET AL.

2822

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

instability takes place to magnify Ke. On the other hand, negative KmKe (PeKe) corresponds to coalescence ofeddies (restratiﬁcation) that leads to ﬂow stabilization.
Following the analysis of Sasaki et al. (2014), we examined the wave number spectra of KmKe and PeKe foreach of the ﬁrst 3 months averaged vertically over the surface mixed layer to investigate the scales involvedin the Ke source in the upper ocean (Figures 12c and 12d). Evident similarities are observed in both thespectra as well as the w spectra (Figure 12b). The spectral power is highest in April with a decreasing trendwith time. A distinctive primary peak appears at the wavelengths of about 10–20 km (i.e., submesoscale),particularly in KmKe in April. Another spectral peak is centered on the wavelengths of about 100–200 km(i.e., mesoscale), which largely dominates over the submesoscale spectral peak in the 2 warmer months.
Based on a submesoscale eddy-permitting simulation of the North Paciﬁc at the horizontal resolution of 1/
308 (3 km), Sasaki et al. (2014) qualitatively drew the same conclusion, where the winter (summer) spectrum peaks were found at the wavelength of about 30 km (125 km). Although their peak wavelengths are
somewhat different from those of the present study, mainly because of differences in the 1 km versus 3 km
lateral grid resolution and in the coastal versus offshore setting, the modeling of both works demonstrates
that submesoscale ASCs in the mixed layer develop most energetically in colder months (i.e., April 2011 forthe present analysis). This is because shear instability acts as the primary Ke source to submesoscale variability,with baroclinic instability acting as the secondary Ke source. It should be noticed that the mean KmKe and PeKe
are about the similar magnitudes (not shown), whereas their wave number spectra are quite different because
KmKe varies in space much more signiﬁcantly than PeKe around the Kuroshio path. It is worth noting that theseresults are slightly different from those of Sasaki et al. (2014) and Capet et al. (2008b), who concluded that baroclinic instability is the primary Ke source rather than shear instability, because our study area contains the vigorous background current, viz., the Kuroshio, and the nearshore region much more extensively than those of
Sasaki et al. (2014) and Capet et al. (2008b). Such conditions are, respectively, responsible for enhancing shear
instability associated with intense lateral shear and topographic eddy shedding as suggested in previous
works (e.g., Dong et al., 2007; Kamidaira et al., 2017; Uchiyama et al., 2017a).

5. Vertical 137Cs Flux
ASCs induced by submesoscale eddies are anticipated to promote vertical mixing and resultant verticaltracer transport. The SGF-L2 model manifests that 137Cs mostly remains in the upper ocean (Figure 10d)and serves as the source of the downward 137Cs ﬂux. The vertical diffusive 137Cs ﬂux, 2Ks @z C, is three tofour orders of magnitude smaller than the advective ﬂux; thus, we omit it in the following analysis. Speciﬁcally, on average, the diffusive ﬂux at z 5 2100 m is merely 0.2% of the advective ﬂux. We examined thevertical transport and resultant mixing of 137Cs due to ASCs in relation to the vertical advective 137Cs ﬂux,
which consists of several components. The total vertical advective 137Cs ﬂux (TCF) can be decomposed into
the mean 137Cs ﬂux (MCF) and eddy 137Cs ﬂux (ECF) using a low-pass Butterworth frequency ﬁlter as:
TCF5Cw  Cw1C 0 w 0  MCF1ECF;(5)where C is the concentration of dissolved 137Cs from SGF-L2, and w is the vertical velocity from ROMS-L2.
The covariance term between the eddy and mean components almost vanishes and thus it can be assumed
zero. To quantify the vertical transport by ASCs, the ECF is decomposed further into the mesoscale (m) andsubmesoscale (s) components, using the 2-D Gaussian ﬁlter described in section 4.1, as:
0 w 0 1C 0 w 0  ECF 1ECF :ECF5C 0 w 0  Cm
m
s
m
s s
137

(6)
137

We ﬁrst evaluated the rate of the Cs inventory (viz., the total amount of Cs) attributable to downward
ﬂuxes against that derived from the FNPP1 as the sum of the oceanic direct release and atmospheric deposition over the SGF-L2 domain. In other words, this rate is the amount of 137Cs transported vertically past
the prescribed depth (dc) during a month (or between 13 March and 30 June) divided by the total amountof 137Cs released from 13 March to 30 June. To obtain the subsurface 137Cs inventory in Bq, we integrated
the downward ﬂux in Bq m22 s21 at a certain depth, dc, over the SGF-L2 domain and the 4 month period
between 13 March and 30 June 2011. The result is summarized in Figure 13 for the inventories beneathfour different dc values: (a) at the MLD (ranging between 10 and 80 m on average; Figure 10e), (b)dc 5 100 m, (c) 200 m, and (d) 300 m. Figure 13a shows that approximately 50% of the FNPP1-derived 137Cs

KAMIDAIRA ET AL.

2823

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 13. Ratios of the 137Cs inventory in the SGF-L2 domain against the total 137Cs load from oceanic direct release andatmospheric deposition. Inventories are estimated from the SGF-L2 model as the time-integrated and area-integrated
downward vertical advective ﬂuxes of the total 137Cs ﬂux (TCF, gray), eddy ﬂux (ECF, black), and submesoscale eddy ﬂux(ECFs, red). Area integrals are performed at four different depths, dc, where (a) dc is at the mixed-layer depth (MLD),(b) dc 5 100 m, (c) dc 5 200 m, and (d) dc 5 300 m. Cyan triangles in Figures 13b–13d are the corresponding TCF estimatesfrom the SGF-M model. Note that in Figure 13d, TCF < ECF because upward MCF occurs in deeper areas.

sank below the mixed layer, based on the TCF during the 4 month period. The contribution of the ECF was
42.7%, while that of the ECFs was 26.5%. Hence, the contribution of the submesoscale ECFs is quite substantial, reaching about 62% of the ECF and dominating over the mesoscale ECFm. The maximum release rates
and maximum inventory of 137Cs also appear in April (Figures 2 and 10e), which has the most intensiﬁed
submesoscale variability (section 4.2). Because of the coincidence of the most intensiﬁed submesoscale variability and the highest near-surface 137Cs in April, the ECFs carried 14.5% of the total input below the mixed
layer in April. This amount reduced over time to 3.5% in May and to 0.1% in June.
Figures 13b–13d show that the ECF dominates over MCF in each plot, conﬁrming the importance of eddyinduced vertical transport to the subsurface 137Cs inventory (Figure 10d). The maximum TCF contribution in
April appears for dc 5 100 m (Figure 13b); however, in May, it appears in the deeper inventories fordc 5 200 m (Figure 13c) and 300 m (Figure 13d). This is considered a consequence of a residual effect of theeddy-induced downward transport that introduces a prominent time lag of about 1 month in the verticaltransport behind the increase of near-surface 137Cs. This is because the downward w values associated with
ASCs are decreased considerably below the mixed layer (e.g., Gula et al., 2014; McWilliams et al., 2015). Nevertheless, the contribution of the submesoscale ECFs is still notable because it accounts for 1/3 of the TCF
below the depth of 300 m for the 4 month period (Figure 13d). As the approximate temporal scale range
for SCSs is hours to days (e.g., McWilliams, 2016), the SGF models might underestimate the effects of SCSs
and ASCs on the vertical 137Cs ﬂux because of the under-resolved temporal variability of the daily averaged
input ﬁeld (section 2.2). If SGF-L2 used velocities averaged over shorter periods or as instantaneous values,
the effect of submesoscale ECFs on vertical transport would appear more intense. However, this is beyond
the scope of the present study and sampling frequencies of the input ﬁeld will be examined elsewhere. On

KAMIDAIRA ET AL.

2824

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Figure 14. Spatial distribution of the monthly averaged submesoscale downward 137Cs ﬂux (ECFs) for (a and d) April,(b and e) May, (c and f) June 2011. (top) ECFs at the mixed-layer depth (MLD). (bottom) At z 5 2200 m. Black contours
represent monthly averaged surface velocity magnitude at intervals of 0.25 m s21. Upward ﬂuxes are blanked out.

the other hand, the cyan triangles in Figures 13b–13d represent the corresponding TCF estimated from the
SGF-M. In general, SGF-M underestimates the total downward 137Cs transport by about 1/3–1/2 that of SGFL2 for the 4 month period. This underestimation is gradually recovered for dc 5 100 m (Figure 13b), presumably because of the imposed background vertical eddy diffusivity, Ks, which compensates the eddy-induced
vertical 137Cs transport detectable in SFG-M.
Figure 14 shows the spatial distribution of the monthly averaged submesoscale downward ECFs at the MLD(Figures 14a–14c) and those at z 5 2200 m (Figures 14d–14f) for the 3 months. Note, only the mean downward ECFs are shown here and the upward ECFs are blank. The downward ECFs at the MLD is enhanced
extensively on the continental shelf near the coast off the FNPP1 with meridional elongation to the patternfound between 36.38N and 37.68N in April. A high ECFs tongue develops from the bend of the coastline ataround 378N, 1418E toward the Kuroshio Extension Jet, where the high ECFs is maintained along the Kuroshio front. As ASCs diminish in May and June (n.b., Figures 11d–11f and 12b), ECFs is also suppressed. Consistent with Figure 13c, the ECFs at z 5 2200 m (Figures 14d–14f), which is sufﬁciently below the MLD,
prevails in May (n.b., ratios of the 137Cs inventory are 1.0% in April, 1.6% in May, and 0.6% in June). The ECFs
at z 5 2200 m occurs only offshore of the continental shelf, which is remote from the area of high concentration of 137Cs (Figure 7). Because of the combination of the distance from the source and the retarded
w below the mixed layer, the maximum ECFs at z 5 2200 m lags that at the MLD by about 1 month. The
Kuroshio Extension Jet is found to be at the center of the area in which deep downward transport occurs.

6. Discussion and Summary
We conducted a detailed hindcast of oceanic 3-D dispersal and initial dilution of dissolved 137Cs, released
from the FNPP1 accident, using a double-nested ROMS model embedded within an assimilative oceanic
reanalysis (MOVE-WNP), coupled with an ofﬂine Eulerian 137Cs transport model (the SGF model). The coupled SGF-L2 model had 1 km horizontal grid resolution, which was sufﬁcient to represent SCSs and associated ASCs. A corresponding coarse-resolution model, SGF-M, was conﬁgured with a 10 km horizontal grid
resolution, which enabled us to perform comparative assessments of the effects of SCSs on dispersal anddilution. We carried out extensive model-data comparisons against satellite altimetry data and ﬁeld observations to demonstrate that both the SGF-L2 and SGF-M models could satisfactorily reproduce the Kuroshio,
3-D oceanic structure, and 137Cs concentrations.

KAMIDAIRA ET AL.

2825

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

We exhibited that submesoscale dynamics are essentially important to better representation of the eddyﬁeld and the resultant 137Cs transport. The submesoscale eddy-resolving SGF-L2 outperformed the coarseresolution SGF-M, producing better agreement with the 3-D patterns of 137Cs and an equivalent temporalvariability of the surface 137Cs near the FNPP1. For the initial dilution of the released 137Cs, the SGF-M overestimated the cross-shelf transport, which led to unrealistically prompt entrainment of 137Cs by the Kuroshio, whereas SGF-L2 manifested greater coastal trapping, alongshelf elongation, and temporalpersistence in the near-surface 137Cs distribution. In addition to the improvement in the lateral structure,
SGF-L2 was able to represent the intensiﬁed vertical mixing and transport of 137Cs due to ASCs. According
to our energy conversion analysis, these ASCs occur most intensively on the submesoscale at wavelengthsof 10–20 km, primarily because of shear instability, with baroclinic instability as the secondary mechanism.
The emergence of SCSs has prominent seasonality, in which submesoscale intensiﬁcation occurs in colder
months, which is consistent with previous studies (e.g., Sasaki et al., 2014). The EKE and resultant eddyinduced 137Cs transport were found intensiﬁed most in the coldest month (April 2011) during the study
period, which corresponds to the period immediately after the FNPP1 accident.

Acknowledgments
We are grateful to Norihisa Usui of
MRI, JMA, for his technical support
regarding the data input of MOVEWNP. This study was supported by theJapan Society for the Promotion of
Science grant-in-aid for ScientiﬁcResearch (KAKENHI grants: 15H00977,
15H04049, and 15KK0207). The
supporting information contains
simulation data that enables theﬁgures to be reproduced. The model
output and computing codes used inthis study are available from theauthors upon request (kamidaira.yuki@jaea.go.jp). The oceanic release rate
data of 137Cs used in the oceanic
dispersion simulation are available inthe supporting information. The
atmospheric release rate data of 137Cs
used in the oceanic dispersion
simulation are available at https://www.atmos-chem-phys.net/15/1029/
2015/acp-15-1029-2015-discussion.
html. The MOVE-WNP product may beavailable at the MeteorologicalResearch Institute upon a request. The
SRTM30_Plus data are available at
ftp://topex.ucsd.edu/pub/srtm30_plus/. The J-EGG500 data set is distributedby the JODC at http://jdoss1.jodc.go.
jp/cgi-bin/1997/depth500_ﬁle. The
MSLA (mapped sea-level anomaly) in
delayed-time product was generated
by Ssalto/Duacs and distributed by
AVISO (ftp.aviso.altimetry.fr). The GPVGSM/MSM atmospheric product was
generated by the JMA and it isavailable from the repository at http://
database.rish.kyoto-u.ac.jp/arch/jmadata/data/gpv/original/. The COADS
data are maintained and distributedby NOAA at https://www.esrl.noaa.
gov/psd/data/gridded/data.coads.
1deg.

KAMIDAIRA ET AL.

To investigate the submesoscale inﬂuences on the vertical transport of 137Cs, we introduced frequency andwave number ﬁlters for the decomposition into the mean (seasonal), mesoscale, and submesoscale components. During the ﬁrst 4 month period following the accident, about 50% of the total input of the released137
Cs within the study domain was transported downward below the mixed layer. Surprisingly, 84% of thatwas induced by eddies, with a contribution of 62% from submesoscale eddies. Far below the mixed layer,
downward transport of 137Cs occurred predominantly because of eddies. For instance, 8.8% (2.4%) of thetotal input was transported downward below z 5 2200 m (300 m) by eddies, with a submesoscale contribution of 3.2% (0.7%), which is still about 1/3–1/2 of the total eddy-induced 137Cs ﬂux. Although submesoscaleASCs barely reach such depths below the shallow mixed layer (e.g., Gula et al., 2014), residual effects might
gradually transport 137Cs vertically to depths >200 m. This is inferred from the increase of 137Cs in deeper
areas that is retarded by about 1–2 months in relation to the increase of surface 137Cs.
Suppose if the FNPP1 accident were to occur in a much colder (warmer) season. Based on the outcomes ofthe present study, submesoscale eddies and SCSs would have been more (less) prominent. This could havecaused immediate (retarded) subduction of high levels of 137Cs into the deeper ocean, potentially provoking greater coastal contamination (immediate offshore dilution with possible entrainment by the Kuroshio)
and thus, more signiﬁcant regional (widespread) consequences. Furthermore, another important aspect canbe deduced from the eddy-driven deep transport of 137Cs that occurs at depth along the Kuroshio Extension. This process could lead to enhancement of the 137Cs supply to subsurface waters (viz., Central Mode
Water and Subtropical Mode Water), which subsequently transport 137Cs via associated deep clockwisecirculating gyres on the southern side of the Kuroshio, returning the FNPP1-derived radionuclides to thewestern North Paciﬁc margins in the long term. Indeed, such multiannual subsurface transport processes
have already been detected (e.g., Kaeriyama et al., 2014; Kumamoto et al., 2017). We should also note that
although the submesoscale eddy-resolving SGF-L2 produced the most satisfactory results without any substantial crude tuning of the model parameters, the coarse-resolution SGF-M has a central role in predicting
an immediate ﬁrst guess of the dispersal of radionuclides for much wider areas. However, this requires thatan appropriate vertical eddy diffusivity Ks be determined a priori to account for seasonally varying submesoscale effects.

References
Aoyama, M., Hult, M., Hamajima, Y., Lutter, G., Marissens, G., Stroh, H., et al. (2016b). Tracing radioactivity from Fukushima in the Nnorthern
Paciﬁc Ocean. Applied Radiation and Isotopes, 109, 435–440.
Aoyama, M., Kajino, M., Tanaka, T. Y., Sekiyama, T. T., Tsumune, D., Tsubono, T., et al. (2016a). 134Cs and 137Cs in the North Paciﬁc Ocean
derived from the March 2011 TEPCO Fukushima Dai-ichi Nuclear Power Plant accident, Japan. Part two: Estimation of 134Cs and 137Cs
inventories in the North Paciﬁc Ocean. Journal of Oceanography, 72, 53–65.
Aoyama, M., Tsumune, D., Uematsu, M., Kondo, F., & Hamajima, Y. (2012). Temporal variation of 134Cs and 137Cs activities in surface waterat stations along the coastline near the Fukushima Dai-ichi Nuclear Power Plant accident site. Geochemical Journal, 46, 321–325.
Badin, G., Tandon, A., & Mahadevan, A. (2011). Lateral mixing in the pycnocline by baroclinic mixed-layer eddies. Journal of Physical Oceanography, 41, 2080–2101.
Becker, J. J., Sandwell, D. T., Smith, W. H. F., Braud, J., Binder, B., Depner, J., et al. (2009). Global bathymetry and elevation data at 30 arcseconds resolution: SRTM30_PLUS. Marine Geodesy, 32(4), 355–371.
Boccaletti, G., Ferrari, R., & Fox-Kemper, B. (2007). Mixed layer instabilities and restratiﬁcation. Journal of Physical Oceanography, 37, 2228–
2250.

2826

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Buesseler, K. O., Dai, M., Aoyama, M., Benitez-Nelson, C., Charmasson, S., Higley, K., et al. (2017). Fukushima Daiichi-derived radionuclides in
the ocean: Transport, fate, and impacts. Annual Review of Marine Science, 9, 173–203.
Buesseler, K. O., Jayne, S. R., Fisher, N. S., Rypina, I. I., Baumann, H., Baumann, Z., et al. (2012). Fukushima-derived radionuclides in the ocean
and biota off Japan. Proceedings of the National Academy of Sciences of the United States of America, 109, 5984–5988.
Callies, J., Ferrari, R., Klymak, J. M., & Gula, J. (2015). Seasonality in submesoscale turbulence. Nature Communications, 6, 6862. https://doi.
org/10.1038/ncomms7862
Capet, X., McWilliams, J. C., Molemaker, J. M., & Shchepetkin, A. F. (2008a). Mesoscale to submesoscale transition in the California Current
System. Part I: Flow structure, eddy ﬂux, and observational tests. Journal of Physical Oceanography, 38, 29–43.
Capet, X., McWilliams, J. C., Molemaker, J. M., & Shchepetkin, A. F. (2008b). Mesoscale to submesoscale transition in the California Current
System. Part III: Energy Balance and Flux. Journal of Physical Oceanography, 38, 2256–2269.
Dong, C., McWilliams, J. C., & Shchepetkin, A. F. (2007). Island wakes in deep water. Journal of Physical Oceanography, 37, 962–981.
Estournel, C., Bosc, E., Bocquet, M., Ulses, C., Marsaleix, P., Winiarek, V., et al. (2012). Assessment of the amount of Cesium-137 released into
the Paciﬁc Ocean after the Fukushima accident and analysis of its dispersion in Japanese coastal waters. Journal of Geophysical Research,
117, C11014. https://doi.org/10.1029/2012JC007933
Gula, J., Molemaker, M. J., & McWilliams, J. C. (2014). Submesoscale cold ﬁlaments in the Gulf Stream. Journal of Physical Oceanography, 44,
2617–2643.
Inomata, Y., Aoyama, M., Tsubono, T., Tsumune, D., & Hirose, K. (2016). Spatial and temporal distributions of 134Cs and 137Cs derived fromthe TEPCO Fukushima Daiichi nuclear power plant accident in the North Paciﬁc Ocean by using optimal interpolation analysis. Environmental Science: Processes & Impacts, 18, 126–136.
Isoguchi, O., Shimada, M., & Kawamura, H. (2010). Characteristics of ocean surface winds in the lee of an isolated island observed by synthetic aperture radar. Monthly Weather Review, 139, 1744–1761.
Kaeriyama, H., Shimizu, Y., Ambe, D., Masujima, M., Shigenobu, Y., Fujimoto, K., et al. (2014). Southwest intrusion of 134Cs and 137Cs derived
from the Fukushima Dai-ichi Nuclear Power Plant accident in the Western North Paciﬁc. Environmental Science & Technology, 48, 3120–3127.
Kaeriyama, H., Shimizu, Y., Setou, T., Kumamoto, Y., Okazaki, M., Ambe, D., et al. (2016). Intrusion of Fukushima-derived radiocaesium into
subsurface water due to formation of mode waters in the North Paciﬁc. Sci. Rep., 6, 22010. https://doi.org/10.1038/srep22010
Kamidaira, Y., Uchiyama, Y., & Mitarai, S. (2017). Eddy-induced transport of the Kuroshio warm water around the Ryukyu Islands in the East
China Sea. Continental Shelf Research, 143, 206–218.
Katata, G., Chino, M., Kobayashi, T., Terada, H., Ota, M., Nagai, H., et al. (2015). Detailed source term estimation of the atmospheric release
for the Fukushima Daiichi Nuclear Power Station accident by coupling simulations of atmospheric dispersion model with improved
deposition scheme and oceanic dispersion model. Atmospheric Chemistry and Physics, 15, 1029–1070.
Kawamura, H., Kobayashi, T., Furuno, A., In, T., Ishikawa, Y., Nakayama, T., et al. (2011). Preliminary numerical experiments on oceanic dispersion of 131I and 137Cs discharged into the ocean because of the Fukushima Daiichi Nuclear Power Plant Disaster. Journal of Nuclear Science and Technology, 48, 1349–1356.
Kawamura, H., Kobayashi, T., Furuno, A., Usui, N., & Kamachi, M. (2014). Numerical simulation on the long-term variation of radioactive
cesium concentration in the North Paciﬁc due to the Fukushima disaster. Journal of Environmental Radioactivity, 136, 64–75.
Klein, P., Hua, B.-L., Lapeyre, G., Capet, X., Le Gentil, S., & Sasaki, H. (2008). Upper ocean turbulence from high-resolution 3D simulations.
Journal of Physical Oceanography, 38, 1748–1763.
Kubota, M., Nakata, K., & Nakamura, Y. (1981). Continental shelf waves off the Fukushima coast. Part I: Observations. Journal of the Oceanographic Society of Japan, 37, 267–278.
Kumamoto, Y., Aoyama, M., Hamajima, Y., Nagai, H., Yamagata, T., Kawai, Y., et al. (2017). Fukushima-derived radiocesium in the western
North Paciﬁc in 2014. Journal of Radioanalytical and Nuclear Chemistry, 311, 1209–1217.
Kunze, E., Klymak, J. M., Lien, R. C., Ferrari, R., Lee, C. M., Sundermeyer, M. A., et al. (2015). Submesoscale water mass spectra in the Sargasso
Sea. Journal of Physical Oceanography, 45, 1325–1338.
Large, W. G., McWilliams, J. C., & Doney, S. C. (1994). Oceanic vertical mixing: A review and model with a nonlocal boundary layer parameterization. Reviews of Geophysics, 32, 363–403.
Le Traon, P. Y., Nadal, F., & Ducet, N. (1998). An improved mapping method of multisatellite altimeter data. Journal of Atmospheric and Oceanic Technology, 15, 522–534.
Levy, M., Klein, P., & Treguier, A.-M. (2001). Impact of submesoscale physics on production and subduction of phytoplankton in an oligotrophic regime. Journal of Marine Research, 59, 535–565.
Mahadevan, A., & Archer, D. (2000). Modeling the impact of fronts and mesoscale circulation on the nutrient supply and biogeochemistryof the upper ocean. Journal of Geophysical Research, 105, 1209–1225.
Mahadevan, A., & Tandon, A. (2006). An analysis of mechanisms for submesoscale vertical motion at ocean fronts. Ocean Modelling, 14,
241–256.
Marchesiello, P., McWilliams, J. C., & Shchepetkin, A. (2003). Equilibrium structure and dynamics of the California Current System. Journal of
Physical Oceanography, 33, 753–783.
Marsaleix, P., Auclair, F., & Estournel, C. (2009). Low-order pressure gradient schemes in sigma coordinate models: The seamount test revisited. Ocean Modelling, 30, 169–177.
Marsaleix, P., Auclair, F., Estournel, C., Nguyen, C., & Ulses, C. (2012). Alternatives to the Robert-Asselin ﬁlter. Ocean Modelling, 41, 53–66.
Mason, E., Molemaker, J., Shchepetkin, A. F., Colas, F., McWilliams, J. C., & Sangra, P. (2010). Procedures for ofﬂine grid nesting in regional
ocean models. Ocean Modelling, 35, 1–15.
Masumoto, Y., Miyazawa, Y., Tsumune, D., Tsubono, T., Kobayashi, T., Kawamura, H., et al. (2012). Oceanic dispersion simulations of 137Cs
released from the Fukushima Daiichi Nuclear Power Plant. Elements, 8, 207–212.
McWilliams, J. C. (2016). Submeoscale currents in the ocean. Proceedings of the Royal Society of London A, 472, 20160117. https://doi.org/10.
1098/rspa.2016.0117
McWilliams, J. C., Gula, J., Molemaker, J., Renault, L., & Shchepetkin, A. F. (2015). Filament frontogenesis by boundary layer turbulence. Journal of Physical Oceanography, 45, 1988–2005.
Miyazawa, Y., Masumoto, Y., Sergey, M., Varlamov, S. M., & Miyama, T. (2012). Transport simulation of the radionuclide from the shelf toopen ocean around Fukushima. Continental Shelf Research, 50–51, 16–29.
Miyazawa, Y., Zhang, R., Guo, X., Tamura, H., Ambe, D., Lee, J., et al. (2009). Water mass variability in the western North Paciﬁc detected in
15-year eddy-resolving ocean reanalysis. Journal of Oceanography, 65, 737–756.
Mori, N., Takahashi, T., Yasuda, T., & Yanagisawa, H. (2011). Survey of 2011 Tohoku earthquake tsunami inundation and run-up. GeophysicalResearch Letters, 38, L00G14. https://doi.org/10.1029/2011GL049210

KAMIDAIRA ET AL.

2827

Journal of Geophysical Research: Oceans

10.1002/2017JC013359

Nagai, T., Tandon, A., & Rudnick, D. L. (2006). Two-dimensional ageostrophic secondary circulation at ocean fronts due to vertical mixing
and large-scale deformation. Journal of Geophysical Research, 111, C09038. https://doi.org/10.1029/2005JC002964
Nagai, T., Tandon, A., Yamazaki, H., Doubell, M. J., & Gallager, S. (2012). Direct observations of microscale turbulence and thermohaline
structure in the Kuroshio Front. Journal of Geophysical Research, 117, C08013. https://doi.org/10.1029/2011JC007228
Nuclear Emergency Response Headquarters, Government of Japan (2011). Report of Japanese government to the IAEA ministerial conference
on nuclear safety—The accident at TEPCO’s Fukushima Nuclear Power Stations—VI, Discharge of radioactive materials to the environment.
Retrieved from http://japan.kantei.go.jp/kan/topics/201106/iaea_houkokusho_e.html
P
erez, F. F., Gilcoto, M., & Rıos, A. F. (2003). Large and mesoscale variability of the water masses and the deep chlorophyll maximum in the
Azores Front. Journal of Geophysical Research, 108(C7), 3215. https://doi.org/10.1029/2000JC000360
Pollard, R. T., & Regier, L. A. (1992). Vorticity and vertical circulation at an ocean front. Journal of Physical Oceanography, 22, 609–625.
Qiu, B. (2001). Kuroshio and Oyashio currents. In J. Steele, S. Thorpe, & K. Turekian (Eds.), Encyclopedia of ocean sciences (pp. 1413–1425).
New York, NY: Academic Press.
Roads, J. (2004). Experimental weekly to seasonal U.S. forecasts with the regional spectral model. Bulletin of the American MeteorologicalSociety, 85, 1887–1902.
Rodriguez, E., Morris, C. S., Belz, J. E., Chapin, E. C., Martin, J. M., Daffer, W., et al. (2005). An assessment of the SRTM topographic products(Tech. Rep. JPL DO31639, 143 p). Pasadena, CA: Jet Propulsion Laboratory.
Romero, L., Uchiyama, Y., Ohlman, J. C., McWilliams, J. C., & Siegel, D. A. (2013). Simulations of nearshore particle-pair dispersion in Southern
California. Journal of Physical Oceanography, 443, 1862–1879.
Sasaki, H., Klein, P., Qiu, B., & Sasai, Y. (2014). Impact of oceanic-scale interactions on the seasonal modulation of ocean dynamics by theatmosphere. Nature Communications, 5, 5636. https://doi.org/10.1038/ncomms6636
Shchepetkin, A. F., & McWilliams, J. C. (2005). The regional ocean modeling system (ROMS): A split-explicit, free-surface, topographyfollowing-coordinate oceanic model. Ocean Modelling, 9, 347–404.
Shchepetkin, A. F., & McWilliams, J. C. (2008). Computational kernel algorithms for ﬁne-scale, multiprocess, longtime oceanic simulations.
In R. Temam & J. Tribbia (Eds.), Handbook of numerical analysis: Computational methods for the ocean and the atmosphere (pp. 119–181).
Amsterdam, the Netherlands: Elsevier.
Terada, H., & Chino, M. (2008). Development of an atmospheric dispersion model for accidental discharge of radionuclides with the function of simultaneous prediction for multiple domains and its evaluation by application to the Chernobyl nuclear accident. Journal of
Nuclear Science and Technology, 45, 920–931.
Thomas, L. N., & Taylor, J. R. (2010). Reduction of the usable wind-work on the general circulation by forced symmetric instability. Geophysical Research Letters, 37, L18606. https://doi.org/10.1029/2010GL044680
Tsubono, T., Misumi, K., Tsumune, D., Bryan, F. O., Hirose, K., & Aoyama, M. (2016). Evaluation of radioactive cesium impact from atmospheric deposition and direct release ﬂuxes into the North Paciﬁc from the Fukushima Daiichi nuclear power plant. Deep Sea Research,
Part I, 115, 10–21.
Tsujino, H., Motoi, T., Ishikawa, I., Hirabata, M., Nakano, H., Yamanaka, G., et al. (2010). Reference manual for the Meteorological Research Institute Community Ocean Model (MRI.COM) Version 3 (Tech. Rep. MRI 59). Tsukuba, Japan: Meteorological Research Institute.
Tsumune, D., Tsubono, T., Aoyama, M., & Hirose, K. (2012). Distribution of oceanic 137Cs from the Fukushima Dai-ichi Nuclear Power Plant
simulated numerically by a regional ocean model. Journal of Environmental Radioactivity, 111, 100–108.
Tsumune, D., Tsubono, T., Aoyama, M., Uematsu, M., Misumi, K., Maeda, Y., et al. (2013). One-year, regional-scale simulation of 137Cs radioactivity in the ocean following the Fukushima Dai-ichi Nuclear Power Plant accident. Biogeosciences, 10, 5601–5617.
Uchiyama, Y., Idica, E. Y., McWilliams, J. C., & Stolzenbach, K. D. (2014). Wastewater efﬂuent dispersal in Southern California Bays. Continental Shelf Research, 76, 36–52.
Uchiyama, Y., Kanki, R., Takano, A., Yamazaki, H., & Miyazawa, Y. (2017b). Mesoscale reproducibility in regional ocean modeling with a 3-D
stratiﬁcation estimate based on Aviso-Argo data. Atmosphere-Ocean, 55, 1–18.
Uchiyama, Y., Suzue, Y., & Yamazaki, H. (2017a). Eddy-driven nutrient transport and associated upper-ocean primary production along theKuroshio. Journal of Geophysical Research: Oceans, 122, 5046–5062. https://doi.org/10.1002/2017JC012847
Usui, N., Ishizaki, S., Fujii, Tsujino, Y. H., Yasuda, T., & Kamachi, M. (2006). Meteorological Research Institute multivariate ocean variational
estimation (MOVE) system: Some early results. Advances in Space Research, 37, 806–822.
Wilmott, C. J. (1981). On the validation of models. Physical Geography, 2, 184–194.
Woodruff, S. D., Slutz, R. J., Jenne, R. L., & Steurer, P. M. (1987). A comprehensive ocean-atmosphere data set. Bulletin of the American Meteorological Society, 68, 1239–1250.
Wright, D. G., Thompson, K. R., & Lu, Y. (2006). Assimilating long-term hydrographic information into an eddy-permitting model of theNorth Atlantic. Journal of Geophysical Research, 111, C09022. https://doi.org/10.1029/2005JC003200

KAMIDAIRA ET AL.

2828