京都大学防災研究所年報 第 56 号 B 平成 25 年 6 月
Annuals of Disas. Prev. Res. Inst., Kyoto Univ., No. 56 B, 2013

Ocean Circulation Response to Inflow from Abukuma River Outlet
Josko TROSELJ (1), Yosuke YAMASHIKI(2) and Kaoru TAKARA(1) Graduate School of Engineering, Kyoto University(2) Graduate School of Advanced Integrated Studies in Human Survivability, Kyoto University
Synopsis
The focus of this study is to promote relevant modeling survey on Land -Ocean
coupling modeling approach applicable for the bay and estuary zone affected by river
inflow and associated pollution from the Abukuma river basin in Japan. The modelingapproach has been studied by combining hydrological model and ocean circulation
model (MSSG model, JAMSTEC) which runs within supercomputer (ES2), andassociated data manipulation techniques for boundary conditions of hydrological andoceanographic modeling. There were simulated 3 different cases, at first circulation ofocean itself induced by its own temperature and salinity data differences, at second
response of the ocean circulation to inflow from the river outlet, and finally response of
the ocean to river inflow with included rotation of the Earth as Coriolis effect.
Keywords: ocean model, river inflow, coastal zone, integrated model, Coriolis effect1.

Introduction

Modeling of a contact zone between a river
coming from potentially contaminated basin and an
ocean is especially sensitive case for Land-Ocean
mutual interaction due to significant risk of major
environmental disaster which can occur in the caseof contamination of the coastal zone. Therefore, it
is of great importance to study and develop
integrated modeling approach to comprehend thecomplex interaction processes in the contact zone inorder to minimize disaster risk potential, which canconsequently cause undesirable social andeconomical costs. Moreover, all land-induced
activities influence the downstream catchment andrelated fishery production.
In order to gain appropriate computationalmodeling skills to comprehend the complex
interaction processes in contact zone between ariver and ocean, GCOE-HSE Advanced Capstone

Project as a long term internship was undertaken at
Japan Agency for Marine-Earth Science andTechnology (JAMSTEC) in Yokohama for 10
weeks. It was studied a hydro 3-D integrated
modeling approach of MSSG model, which
combines coastal and river zones, in order to extend
my research field into larger scale phenomena in
Adriatic Sea affected by the river water intrusion,
which will be a site study field for later part of my
doctoral course.
Abukuma river outlet was considered as asample site for the contact Land-Ocean zone, andthe ocean circulation response to induced inflow
from the river side was considered. I am planning toapply my knowledge into the similar catchment
area in Croatia, particularly into Rijeka
Metropolitan Region (RMR) and many bays along
the Adriatic Sea coastline, in order to evaluate
potential risk and security level for the downstream
area.

― 87 ―

2.

Geographical position

The Abukuma river watershed is located in thenortheast of Japan mainland, and has an area of
5,400 km2 of which forest and agricultural land use
accounts for 79% and 18% respectively. Most of
the mainstream runs in the middle of the watershedfrom the south to the north. Many tributaries run
eastward to westward. The headwater starts from
Asahidake (elevation of 1,835 m) and reaches down
the Pacific ocean. The elevation of the watershed
ranges from the sea level to more than 2,000 m on
the west side, where a volcanic mountain range islocated. On the contrary, a peneplain is dominanton the east side of the mainstream, forming
mountains of low elevation and gentle river slope(Kinouchi and Musiake, 2008). Figure 1 represents
Abukuma watershed and two lakes located outside
the watershed, and Photo 1 shows satellite image ofAbukuma river mouth.

Photo 1 Abukuma river mouth satellite image(Adityawan et al., 2012)3.

Fig 1 Abukuma watershed and two lakes locatedoutside the watershed (Kinouchi and Musiake,

2008)

Model Description

ETOPO1, 1 Arc Minute Global Relief Modelwas chosen as initial database for bathymetry data,
while World Ocean Atlas 2005 dataset was chosen
as initial 3D database for temperature, salinity,pressure, and velocity field data. The ocean model
was set up rectangular with 200*200 cells withresolution of 300 meters, 60 km range in longitude
axis with 1 day output data for ocean only, and 1
hour output data increment for river input, startingfrom 1st of June 2005. River outlet was positionedon the central west point of the model. In order tosimplify the model physics, the ocean depth was set
flat at 52 meters depth with layers of 1 meterincrements, which is sufficiently deep so that river
plumes never touch the seafloor (Isobe, 2005). The
river side boundary conditions were simulated byextending the model bathymetry throughout theriver mouth at 38.05 N latitude, with constant 1200
meters width and 7 meters depth, and the river flow
input was set to 1300 m3/s with uniform velocity
distribution over the cross section at the entrance
point on the central west side of the domain. The
upper, lower and right borders of the domain haveremained closed. Finally, it was simulated
additional boundary condition into the model,
which included existence of the Coriolis effect.

― 88 ―

Fig 2 Barometry of Japan adjacent area [m]Barometry data of Japan adjacent area fromETOPO1 (Amante, Eakins, 2009) dataset arepresented in Figure 2. It was chosen mainly because
its surveys were collected by several various ocean
research Institutes which revisited limited numberof localities during surveys, and surveys were donewith the most modern technology by using GPS
navigation Acoustic Doppler Current Profiler(ADCP) and databases were integrated together andsmoothed in order to prevent improper steps, so
working with that type of dataset was very simple
and reliable.
World Ocean Atlas 2005 (NODC, 2007)
dataset was chosen as 3D database for temperature,salinity, pressure, and velocity field data mainly
because surveys for the database were done in situ
and it contains many various statistical analyses, aswell as available basic sources of all calculations ofdatasets. The database was made for 33 standard
depths from ocean surface to depth of 5,500 meters.

Fig 4a Position of the model
Temperature data on ocean surface layer of Japan
adjacent area for 1st of June 2005 are presented in
Figure 3. Figures 4a and 4b represent position anddomain for the model, with red line on Figure 4b
representing direction of river inflow perpendicular
to direction of coastline in the river mouth.

Fig 3 Temperature data of ocean surface level ofJapan adjacent area (01/06/05) [K]― 89 ―

Fig 4b Domain of the model (Google Earth, 11 thof February 2013)

3.1 Model equations and methodology
The model is operating by using Navier Stokes
equation (1).

(1)
With simplifications gotten from using finite
elements method, the equation (1) can betransformed into equation (2).

(2)
This equation represents temporal changing of
velocity field data simulated by the model. Aftersetting up initial conditions for the model, it havebeen continuously improved by adding moreboundary conditions, in order to increase precisionof simulated results and to get more reliable output
data. It is done by changing the "f" variable of theequation, which represents external influences tothe hydrodynamic physical processes. Adding more
complex boundary conditions are next steps in mydoctoral study.
Initial conditions of velocity field were static,
and ocean was in non equilibrium quasi stationary
state. Fluid circulation was induced for 3 different
initial conditions, one when ocean circulation was
induced by its own temperature, salinity andpressure differences among adjacent cells, second
when initial side boundary conditions from riverside inflow was simulated together with self
induced ocean, and the most advanced initialconditions when external factor of the Earthrotation was put into the simulation. Equation (3)proves that the ocean with given initial conditions
are in non equilibrium quasi stationary state, while
Figure 5 schematically presents how thecirculation of ocean is transformed from that state
into dynamic state.
Ep + Ek = const;
t=0: Ep1 + Ek1 = Ep2 + Ek2 ; Ek1=Ek2=0
ρ1gh = ρ2gh; ρ1 = ρ2 → Equilibrium state→ T1 ≠ T2, sl1 ≠ sl2 → ρ1 ≠ρ2→ Quasi stationary non equilibrium stateFig 5 Scheme how the ocean circulation process
is induced
The principle of available potential energy
describes mechanism of physical processes for
inducing the ocean circulation process, by using
energy conservation principle, as it is described inequation (4), where Ep(ρ ≠const) is available potential
energy, ρ0 steady potential density (baroclinic
stability) and ρ’ available potential density whichcauses baroclinic instability.
Ep(ρ ≠const)= ρ0gh + ρ’gh

(4)

The circulating system inside closed borders
will always tend to have uniform distribution oftemperature, salinity and pressure data, in whichcase it would stop circulating. However, even withsmall domain like this and relatively small
differences of initial data among cells, it was maderough estimation that the ocean induced by its own
data differences would circulate for more than 10
years before spending all available potential energy
gotten from the density differences.
3.2 Model set up and results, ocean itself
In addition to already mentioned initial data,
initial bathymetry data were also provided as initialconditions for the model, as shown in Figure 6.

(3)
Fig 6 Initial bathymetry conditions [m]― 90 ―

Fig 7a Longitudinal section of temperature data
at 38.05 latitude, 5 days after starting itself
induced ocean simulation [˚C]Fig 7b Longitudinal section of temperature data
at 38.05 latitude, 30 days after starting itself
induced ocean simulation [˚C]

The itself induced ocean simulation were run by
30 days of simulated time. It was optimal
simulation time for chosen domain size. Results
showed that it is enough time to induce ocean tocirculate in relatively steady dynamic phase, which
was one of my concerns before running the model.
In this work I will present results of simulations
which I found the most interesting among all
available output data. Figure 7a shows longitudinal
section of temperature data at 38.05 latitude 5 days,
while Figure 7b shows it 30 days after the
initialization of ocean circulation process.
3.3 Model results, ocean + river
After checking that itself induced ocean
simulation is producing expected and reasonable
results, river side boundary conditions were put into

the model, as described earlier in this chapter. My
primary focus during the research was gaining
computational skills to be able to navigate inside
the MSSG code and modify my input, so I was not
spending a lot of time in analyzing exact conditions
from the Abukuma river and therefore theestimation of river conditions was done roughly andsimplified according to the information given fromreferences (Matsumoto, 1981 and Yustiani, Mano,
2004). In upcoming part of my research I willchange side river boundary conditions to be more
realistic like on the site study. Figure 8a shows
longitudinal section of salinity data for river
domain at 38.05 latitude 4 hours, while Figure 8b
shows it 11 hours after the initialization of thesimulation. I wanted to get output data about lengthof ocean's water intrusion into river mouth and

Fig 8a Longitudinal section (river domain) of
salinity data at 38.05 latitude, 4 hours after
starting the river simulation [‰]Fig 8b Longitudinal section (river domain) of
salinity data at 38.05 latitude, 11 hours after
starting the river simulation [‰]― 91 ―

Fig 9a Longitudinal section (global domain) of
salinity data at 38.05 latitude, 22 hours after
starting the river simulation [‰]Fig 9b Longitudinal section (global domain) of
salinity data at 38.05 latitude, 44 hours after
starting the river simulation [‰]therefore had chosen longer distance of river flowfrom coastline side than in real conditions. At the
beginning, it was filled by ocean water with salinityof around 34‰, while river water salinity was set
up as 0‰, and I was considering salinity gradientduring time throughout the river channel until it
became more or less uniform, which was sign that it
has reached steady stage of mixing waters for the
given conditions. The data after 4 hours the best
represent period when river water overcame theriver mouth, and effect of bottom friction which
occurs by the way. The data after 11 hours represent
the beginning of more or less steady salinity
gradient at the position, so it represents real contact
zone between river and ocean water which will beoccurring during the whole process of thesimulation. It is also interesting to notice position of
the river water which repressed brackish water very
close to the river mouth, which was logically to
expect due to extremely high river flow used as aboundary condition. It would be interesting to see
level of intrusion of brackish water for different
border conditions, which will be one of concerns of
my future study. Figure 9a shows longitudinal
section of salinity data for global model domain at
38.05 latitude 22 hours, while Figure 8b shows it 44
hours after the initialization of the simulation. It canbe seen that relatively steady phase of the brackish
water body will be held during the whole time ofthe simulation. The ocean + river model
successfully simulated 76 hours of data.

3.4 Model results, ocean + river + Coriolis effect
So far this simulation is the most advanced partthat I reached during my doctoral study. As salinity
is factor which has the biggest value's difference
between ocean and river waters, I often used it as
the most valuable indicator about directions oftrajectories of the river water movement. Therefore,
the salinity output data as the best indicator of
influence of the Coriolis effect can be seen at
Figures 10, presented on ocean surface layer 24hours (Fig. 10a) and 48 hours (Fig. 10b) after
inducing the simulation.
4.

Conclusions and follow up

My main concern was about time of putting sideboundary river inflow into the ocean model,
because some time is needed until an ocean reaches
its steady stage of circulation. While analyzing theocean circulation output data, the direction which I
was tending to search for was putting river inflow
into the model on the time when ocean was alreadyrunning for some time period by itself. Therefore,
my goal was trying to simulate much closer natural
conditions where an ocean response is occurring toan river inflow during initial dynamic stage forocean and static stage for river inflow, and notduring initial static stage for both which wassimulated in this work. Finding computational
solution in the model which will solve this problem
is one of my major goals for upcoming research.
I was searching for conditions when it can be― 92 ―

Fig 10a Salinity data on ocean surface layer, 24hours after inducing the simulation with bothriver inflow and Coriolis effect [‰]Fig 10b Salinity data on ocean surface layer, 48
hours after inducing the simulation with bothriver inflow and Coriolis effect [‰]concluded that the ocean reached its steady stage ofcirculation. By my assumptions, the range between
minimal and maximal elevation can be good
indicator of that. When the range reaches more orless constant value, then it can be concluded thatthe ocean circulation has reached its steady stage.
By following that assumption, I made analysis of
maximal and minimal ocean surface elevation, as
well as range between them in every time increment,
as it is shown in Figure 11. It can be concluded thatfrom 9th day after initialization the ocean reaches
its steady stage for the current initial conditions andthis kind of test. As ocean had initial velocity field
zero, its movement was induced only because of
temperature and salinity differences among adjacentcells. Density among adjacent cells was not thesame because of the temperature and salinity
differences, and therefore nor energy field was thesame, which caused instability and initially induced
movement among cells. A cell which has higher
density will have higher mass than the one withlower density because volume of every cell is thesame at the beginning. The cell with higher mass
will be induced by higher gravity force as well as
higher pressure force, which will produce rotational
moment of the force around mid plane between two
adjacent cells induced by gravity force, as well astranslational movement of the cells in direction ofhigher pressure force, the process which will be
repeated in every time increment. Figures 12
represent free ocean surface elevation at the first(Fig. 12a), fifth (Fig. 12b) and twenty-fifth (Fig.

Fig 11 Maximal and minimal ocean surface
elevation, and their mutual range for itself
induced ocean model [m]{m}

Fig 12a Ocean surface elevation, 1 day after
inducing the ocean itself circulation [m]― 93 ―

Fig 12b Ocean surface elevation, 5 days after
inducing the ocean itself circulation [m]Fig 12c Ocean surface elevation, 25 days after
inducing the ocean itself circulation [m]12c) day after inducing the ocean movement.
The other main concern of special interest forme was momentum conservation principle.
Following that principle, every force which isinduced due to temperature and salinity differencesamong adjacent cells will also produce a counter
force on the same line but in the opposite direction.
The important conclusion following this is that totalsum of rotational moments of every force in thesystem around every point will be zero, because
initial conditions were related to static ocean andthere was not any external force entering into themodel, and the induced forces and counter forces
will always produce the same value but the opposite
angle of rotational moment around every point in
the system. The principle becomes much more
complicated when side boundary river input isadded to the system, because it induces continuous
force of water body uniformly distributed at theplace of its appearance. Because of the fact that it isuniformly distributed, it can be concluded that there
exists one line in the system, in the middle point of
the river input, where sum of rotational moments of
the whole system in every time increment will be
zero. I am not sure if this conclusion is correct, butthat is idea which came to my mind based on
previous assumptions, and is for sure something
where my future research will be directed. In order
to find some evidence proving this principle, I madeanalysis of total sum of velocity directions for theavailable output data, and my expectations werethat total sum of both x and y direction velocities
will be around zero, with a little deviations becauseof vertical mixing of the ocean water which was notconsidered in this observation as well as ocean
surface elevations, so not every time the same
velocity gradient produces the same force becauseof the existing mass gradient between adjacent cells.
Figure 13 represents results of the total sum of x
and y velocities during time inside the wholesystem for the ocean model itself. It mostly showed
results of less than 3 millimeters per second, which
is negligible value in comparing with total size ofthe system. The mass conservation principle test
was also done, and it showed average elevation ofocean surface layer of 1,31*10 -4 m over 77,53% of
cells which had water level data, which isequivalent to 3,65*10 5 m3 of water. It is only
0,13‰ of the total volume, which is added to thesurface layer just at the time from first to second

Fig 13 Total sum of x and y velocities inside thewhole system for ocean model itself― 94 ―

time increment of the modeling. On every following
time increments, elevation remains constant which
proves that mass conservation principle isconsistent. The starting lifting of the water level
data can be explained with numerical model error,
because the model itself also makes an calculating
error when moving from static to dynamic phase of
the ocean circulation. One of my concerns forupcoming research will for sure be including
integrated 3D tests of the conservation principles
instead of doing it in 2D level.
For ocean itself model, average surface
temperature was reduced from 17.05 to 16.65 ˚Cwhile the average bottom temperature wasincreased from 13.02 to 13.40 ˚C in the range of 30days. Average ocean surface salinity in the same
timeframe was increased from 33.88 to 33.92 ‰,while average bottom layer salinity was decreased
from 34.23 to 34.21 ‰.
Acknowledgements
This research is part of Ph.D. study of the firstauthor, which is supported by Ministry of
Education,
Culture,Sports,ScienceandTechnology - Japan (MEXT). It is continuous
research started as long term internship in Japan
Agency for Marine-Earth Science and Technology(JAMSTEC), which was supported by Global
Center for Education and Research on Human
Security Engineering for Asian Megacities(GCOE-HSE) program.
I would like to express my honest gratitude toProfessor Dr. Yosuke Yamashiki for initializing
process of my application for the internship and forproviding me with all technical support and advices
which I needed during the time, GCOE-HSE
program leaders for providing logistical support for

the internship, Dr. Keiko Takahashi, program
director of the MESSAGE model for accepting myapplication to complete the internship in JAMSTEC,
as well as Dr. Shinichiro Kida to be my executive
supervisor in behalf of JAMSTEC. They all,together with my Kyoto University supervisor Dr.
Kaoru Takara, made significant contribution to myeducation progress.
References
Adityawan, M.B., Roh, M., Tanaka, H., Farid, M.
(2012): The Effect of River Mouth Morphological
Features on Tsunami Intrusion
Amante, C., Eakins, B.W. (2009): ETOPO1 1
Arc-Minute Global Relief Model: Procedures,
Data Sources and Analysis, NOAA Technical
Memorandum NESDIS NGDC-247, 19 pp
Isobe, A. (2005): Balooning of River-Plume Bulge
and Its Stabilization by Tidal Currents, American
Meteorological Society
Kinouchi, T., Musiake, K. (2008): Simulatinghydrological impact of environmental change in
the Abukuma Watershed, Japan, 4th APHW
Conference
Matsumoto, H. (1981): Developmental Process of
Alluvial Coastal Plain related to the Holocene
Sea-level Change, The Science Reports of Tohoku
University, 7th series (Geography), Vol 31, No. 2
National Oceanographic Data Center (NODC,
2007)
Yustiani, Y.M., Mano, A. (2004): Flood Effect on
Nutrient Distribution in Open Coastal Waters off
the Abukuma River, Asian and Pacific Coasts
2003: Proceedings of the 2nd International
Conference, Makuhari, Japan(Received June 11, revised version July 26, 2013)― 95 ―