UDK/UDC: 519.876.5:551.465.11(520) Pregledni znanstveni članek - Review scientific paper
Prejeto/'Received: 2. 2. 2017 Sprejeto/Accepted: 16. 6. 2017
Calibration of Three-dimensional Model PCFLOW3D - A Comparison WITH MEASUREMENTS FROM THE YATSUSHIRO SEA AND THE SEA OF JAPAN/
Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z meritvami v morju Yatsushiro in Japonskem morju
1 Fakulteta za gradbeništvo in geodezijo, Univerza v Ljubljani, Jamova 2, 1000 Ljubljana
Abstract
PCFLOW3D is a three-dimensional numerical model that simulates the hydrodynamic circulation and transport/dispersion of various pollutants. As the first hydrodynamic (HD) part is crucial for the reliable simulation of further transport/dispersion processes, we have been diligently developing it for several years, calibrating and verifying it at every stage of development. This paper presents two cases where it was possible to draw a comparison with current measurements or observations, and describes several improvements to the model that proved necessary in order to simulate the desired range of complex hydrodynamic situations. In the first case (Case A) we made good use of the available, quite detailed measurements of currents in the Yatsushiro Sea (Japan). Due to the very high and variable current velocities we had to include a better turbulence model (Smagorinsky) to get acceptable agreement with the measurements. In the second case (Case B) we simulated HD circulation in the Sea of Japan, where thermohaline forcing is very important; this called for the model to be completed with a special term, namely so-called topographic stress. In both cases the final agreement is good; in Case B we also describe the simulation of the transport/dispersion of radioactive pollutants in the Sea of Japan.
Keywords: three-dimensional numerical model, current measurements, radioactive pollution, Yatsushiro Sea, Sea of Japan.
PCFLOW3D je tridimenzionalni numerični model za simulacijo hidrodinamičnih tokov in za prenos in širjenje onesnažil. Ker je prvi - hidrodinamični - del zelo pomemben za zanesljivo simulacijo transportno/disperzijskih procesov, smo ga razvijali več let in ga v vsaki fazi umerjali in preverjali. V tem prispevku opisujemo dva primera, kjer je bila možna primerjava z rezultati meritev, in opisujemo nekaj dopolnitev modela, potrebnih za zanesljivo simulacijo nekaterih hidrodinamičnih pojavov. V prvem primeru (primer A) smo izkoristili razpoložljive, zelo detajlne rezultate meritev v morju Yatsushiro (Japonska). Zaradi izredno velikih in spremenljivih hitrosti tokov smo morali vključiti turbulentni model Smagorinsky, da smo dobili sprejemljivo ujemanje izračunov z rezultati meritev. V drugem primeru (primer B) smo simulirali cirkulacijo v Japonskem morju, kjer spremembe gostote vode zaradi spremenljive temperature in
* Stik/Correspondence: rudi.rajar@fgg.uni-lj .si
© Rajar R., Četina M.; This is an open-access article distributed under the terms of the Creative Commons Attribution -NonCommercial - ShareAlike 4.0 Licence.
© Rajar R., Četina M.; Vsebina tega članka se sme uporabljati v skladu s pogoji licence Creative Commons Priznanje avtorstva -Nekomercialno - Deljenje pod enakimi pogoji 4.0.
1 *	v
Rudi Rajar ' , Matjaž Četina
1
Izvleček
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
slanosti pomembno vplivajo na hitrosti toka. Model smo morali dopolniti s t. i. topografsko napetostjo (angl. topographic stress). V obeh primerih smo dosegli dobro ujemanje z rezultati meritev. V drugem primeru opisujemo tudi simulacijo transporta in disperzije radioaktivnih onesnažil v Japonskem morju.
Ključne besede: tridimenzionalni numerični hidrodinamični model, meritve tokov, radioaktivna onesnažila, morje Yatsushiro, Japonsko morje.
1. Short description of the PCFLOW3D model
The three-dimensional model PCFLOW3D was developed over the past three decades at the University of Ljubljana, Faculty of Civil and Geodetic Engineering, Chair of Fluid Mechanics. In fact this is a guasi-three-dimensional (3D) model, which means that it can reliably simulate flows, where the vertical velocity components are of an order of magnitude smaller than the horizontal components. However this condition is almost always fulfilled for flows in lakes and seas, often also with flows in rivers, and the model has been successfully used to solve many practical problems.
The model consists of four modules:
•	Hydrodynamic (HD) module - We use the hydrodynamic equations ("momentum equations") to calculate the velocity field. The model also includes the equations for simulating temperature and salinity distribution. As the variable temperature and salinity influence the density and further on the velocity field, these equations are coupled with their HD equations.
•	Transport-dispersion (TD) module - TD equations describe the transport and dispersion of a "pollutant". This term we use for any parameter that is subjected to TD processes in the water, such as: salinity, heat, oil, mercury, or any other dissolved matter.
•	Bio-chemical (BC) module - This includes the equations describing the bio-chemical reactions of the pollutants. Most often these equations are added to the TD equations as "source terms", since they describe an increase (source) or decrease (sink) in the amount of the pollutant due to biochemical transformations. This module also includes equations to determine heat processes. This is calculated using the same subroutine as any
pollutant, where the temperature is the variable that directly expresses thermal processes.
•	Sediment transport (ST) module - Only
suspended sediment is simulated in the PCFLOW3D model (the bed load is not included). Simulating the transport/dispersion of suspended sediment is important in water-quality models, as usually an important part of the pollutant is bound to the suspended sediment and is transported with it.
The basic equations of model PCFLOW3D are best described in Rajar and Četina (1997).
The expression "forcing factors" is usually used with the same meaning as boundary conditions. These are physical factors that act at the boundaries of the computational field and cause the movement of water. Like most complete 3D HD models, the PCFLOW3D model includes all four boundary conditions:
a)	Wind forcing - At the water surface the shear stress caused by the wind is one of the forcing factors that cause water movement.
b)	Inflow of rivers - The water mass flowing into the computational region will cause water movement through its momentum and viscosity.
c)	Tidal forcing - At the open boundaries of the computational region the tidal water level oscillations cause the flux of water inflowing/outflowing into the region.
d)	Thermohaline forcing - Space and time-varying water density will cause water movement. In the model density is a function of temperature and salinity. Variations of these two parameters can be caused by several phenomena, e.g.: thermal exchange with the atmosphere, inflow of river with variable salinity, inflow of heat from a thermal or nuclear power plant, and more.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
It is important to know the importance of every forcing factor. In several cases we can drop out forcing factors that have a very small or negligible influence on the phenomena.
Several specific adaptations of the model were made during its applications. Some of these adaptations are described with the two cases in the following chapters.
The PCFLOW3D model has already been applied to many studies of surface water quality problems. An adapted version of the model for simulating oil spills is operational in Slovenia. The bio-geo-chemical cycling of mercury, originating from a former mercury mine in Slovenia, has been studied using both 2D and 3D models (Rajar et al., 2000). The model has also been used in the international study on the impact of French underground nuclear tests on the environment, where the outflow of sediments and plutonium from the Mururoa lagoon was simulated for various weather conditions (Rajar, 1998).
2. CASE A: Simulation of currents in the Yatsushiro Sea: detailed calibration of the model
2.1. Introduction
This project begun with research into mercury pollution in the Yatsushiro Sea and Minamata Bay (Rajar et al., 2002; Cetina et al., 2003). Hydrodynamic circulation is the basic phenomenon at work; further on we included a simulation of the transport and dispersion of total mercury, which is a useful tool in analyzing the degree of mercury contamination and proposing measures for remediation.
However, the most important point of the case presented here is a thorough verification of the hydrodynamic (HD) module of the PCFLOW3D model. In most cases it is extremely difficult to obtain all the data necessary for the verification and calibration of 3D models. There is an enormous amount of input data that are necessary for reliable simulations, and rarely are all of these data available.
This case was a rare opportunity where all important data were measured. In the frame of Slovene-Japan collaboration our Japanese colleagues provided us with the data containing very detailed measurements of the current velocities in the Yatsushiro Sea (Minamata Bay being its integral part) and we were able to make a relatively reliable simulation of hydrodynamic circulation and compare the simulated and measured results. There was a kind of lucky coincidence in this case in that only one forcing factor was essential - tidal forcing - which was so extremely strong that it was possible to assume that the influence of the river's inflow and thermohaline forcing were negligible. The wind speed was measured; however its influence was also negligible.
Only the (3D) hydrodynamic module of the PCFLOW3D model was applied and verified in this case.
2.2. Adaptions of the PCFLOW3D model
As in this example the flow velocities varied, with extremely broad limits (from U = 0.0 to circa 3.5 m/s) in space and in time, it was paramount to include a turbulence model that can simulate such a flow. The Smagorinsky equation was applied for the horizontal turbulent viscosity. For the vertical turbulent viscosity constant values for each layer were applied, and parabolic distribution over the depth was assumed. The final results justify the appropriate application of this method. In the basic Smagorinsky equation (Bombač, 2014) the Smagorinsky coefficient is not precisely defined and its value can vary in relatively wide limits (Bombač, 2014). In the manual of 3D numerical model POM (POM, 1991), the value of 0.2 is proposed, while in the program package MIKE 21 (DHI, 2012) the recommended value is between 0.25 and 1.0. Some trial simulations have shown that with the value of 0.2 the simulations are partly unstable, probably due to the extreme velocity gradients. In the final simulations the value of 0.35 showed the best results. Final simulations with the described parameters showed that this turbulence model has very good properties (see Section 2.5).
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
Constant values for each layer were applied for the vertical turbulent viscosity. Parabolic distribution over the depth was assumed with the value of 0.08 in the bottom layer, 0.10 at the surface, and 0.12 as maximum value at middle depth.
As the horizontal velocities were very high and it was possible to assume that the density stratification is negligible, this approximation yielded good results.
2.3. Description of the Yatsushiro Sea and
measurements The Yatsushiro Sea is an almost enclosed sea, along the eastern side of one of the major Japanese islands, namely Kyushu Island (Figure 1). Its length, in the SSE-NNW direction, is 70 km; its width is up to 25 km in the southern part and about 10 km in the northern part. Its northern part is shallow, with the depth being around 10 to 20 m, while in the southern part max. depths are 35 to 60 m.
There are five straits connecting the Yatsushiro Sea with the neighboring seas. Table 1 gives basic dimensions and information about tidal constants in the straits, as they were used in our simulations. Strait C (Hondo Seto) is very narrow and shallow, and the discharge through it is negligible in comparison with the flow through other straits, and it was neglected as boundary opening in all our simulations.
There are very strong tidal level oscillations in the Sea. Table 1 shows e.g. that in the Misumi strait the maximum tidal span is 4.41 m. Measurements of currents (Hydrographic Department, 1978) have shown maximum velocities in Strait A to be about 4.9 knots (2.5 m/s), and in B even as much as 6.9 knots (3.55 m/s).
Measurements were taken with the established method of using flow meters on ships. The time of measurements was registered to synchronize the measured velocities with the tidal cycle.
1
A
tí
Figure 1: Yatsushiro Sea with the marked boundary openings: A: Ushibuka (Kurono Seto); B: Nagashima Kaikyo; C: Hondo Seto; D: Matsushima; E: Misumi Ko.
Slika 1: Morje Yatsushiro z označenimi lokacijami odprtih robov: A: Ushibuka (Kurono Seto); B: Nagashima Kaikyo; C: Hondo Seto; D: Matsushima; E: Misumi Ko.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
Table 1: Main characteristics of the 5 straits of the Yatsushiro Sea, with tidal range (2*ampl.) given for spring tide on 19-20 August 2001.
Preglednica 1: Glavne značilnosti 5 ožin morja Yatsushiro z višino plimovanja (2*amplituda) za visoko plimo dne 19. - 20. avgusta 2001.
Strait	Width (m)	Depth (m)	Tidal range (m)	Time lag (deg.)
A - Ushibuka	400	10	3.30	224.204
B - Nagashima Kaikyo	1800	61	2.72	207.834
C - Hondo Seto	<100	3		
D - Matsushima	1500	7.5	3.60	178.44
E - Misumi	300	26	4.41	178.81
Figure 2: Boundary conditions at open boundaries for Case A2. Z(t) - tidal water level oscillation; U(t) and V(t) - measured velocity components.
Slika 2: Robni pogoji na odprtih robovih za primer A2. Z(t) - oscilacije gladine plimovanja; U(t) in V(t) -merjene komponente hitrosti.
2.4. Input data for simulations
Topography: Relatively accurate topographic data were available from measurements. In the numerical model the computational region, with horizontal dimensions of 46,600 x 56,650 m, was divided into 105 x 122 x 17 control volumes (CVs). The maximum depth of the model was 61 m, divided into 15 + 2 = 17 layers, with the width of 5.0 m in the deepest part and 2.5 m near the surface.
Boundary conditions: As already mentioned, the hydrodynamic circulation was mainly governed by the tidal oscillations. Here we must describe two simulated cases; two different series of boundary conditions (BC) were used.
Case A1: The most logical boundary conditions were to use tidal level oscillations z(t) at the four
straits. However, initial simulations did not agree with the measurements to a great extent. After detailed examination of the phenomenon, we understood that at the openings D and E the topography of several small islands outside of the Yatsushiro Sea is preventing free flow, caused by tidal oscillations. Thus the use of tidal level oscillations as boundary conditions is not physically realistic.
Case A2: Further simulations were made, named Case A2, where we used BC z(t) at two open boundaries (straits); however, at Straits D and E we used velocity components U and V as functions of time as they were determined from the measurements. The functions are shown in Figure 2.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
2.5. Results and comparison with measurements
As was mentioned, the simulations of Case A1 did not yield satisfactory results. Only the results of case A2 are shown here: Figure 3 shows comparison of measured and simulated results at the time of approx. highest inflow velocities, and Figure 4 at the time of highest outflow velocities.
The two figures show good agreement of the simulated and measured current vectors. The minor differences can be assigned to the inevitable inaccuracies of the measurements and to the impossibility of presenting the topography precisely by using the numerical model.
2.6. Conclusions - Case A
This case shows a confirmation of the relative good accuracy of the hydrodynamic part of the PCFLOW3D model. A denser numerical grid, especially near the open boundaries, could still increase the accuracy of the numerical results. Some findings about the application of the turbulence model are also important for further applications. The Smagorinsky principle was found to be relatively simple but effective.
The model applied here was later used as a basis for simulating several phenomena, where it was necessary either to determine the current velocities or transport/dispersion of pollutants. The numerical grid described herein was used to determine the transport/dispersion of mercury from the Minamata Bay to the Yatsushiro Sea (Rajar et al., 2002).
YatsushiroSea CaseA2
Figure 3: Case A2: Comparison of calculated and measured (red) velocity vectors at the time of approx. maximum inflow velocities.
Slika 3: Primer A2: Primerjava izračunanih in merjenih (rdeče) vektorjev hitrosti približno v času največjih vtočnih hitrosti.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
Yatsushiro Sea Case A2
Figure 4: Case A2: Comparison of calculated and measured (red) velocity vectors at the time of approx. maximum outflow velocities.
Slika 4: Primer A2: Primerjava izračunanih in merjenih (rdeče) vektorjev hitrosti približno v času največjih iztočnih hitrosti.
3. CASE B: Modelling of circulation and dispersion of radioactive pollutants in the Sea of Japan
3.1. Introduction
The Sea of Japan is a large, almost enclosed sea between the Japanese islands, Sakhalin Island, the Siberian coast, and Korea. Its length is about 1800 km, width 900 km (see Figure 5a), and the maximum depth is about 3900 m. The surface area is 1.01 x 106 km2 and the volume 1.36 x 106 km3. The relatively weak exchange of water with the Pacific Ocean is enabled through four straits and is described in section 3.3.
A large amount of solid radioactive wastes were deposited in several locations on the sea floor by the former Soviet Union, at a depth of about 3000 m, mainly at stations No. 6, 9, and 10 off Vladivostok (see Figure 5a). The dumped solid wastes consisted of 4 reactor vessels of nuclear
submarines deposited without spent nuclear fuel (178 TBq) as well as waste packed in containers (129 TBq). The total activity dumped in the period 1968-1992 was about 307 TBq. With the aging (corrosion) of containers and reactor vessels, radionuclides from the wastes could be released into the marine environment. The estimated release rates (for all radionuclides together) are about 3 TBq/year. The most important radionuclides that could be released from the wastes are 137Cs and 90Sr. For modelling purposes a hypothetical release was assumed in which 1 TBq of 137Cs would be released in a 90-day period at location No. 9, which should be a very conservative approach. Although several measurements of radionuclide concentration in this sea have shown no increase of radioactivity over the background values, concern exists that there may occur some more significant leakage of radionuclides from the dumped waste. The question is where they might be transported by currents and dispersion, and especially after what
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
time and in what concentrations they would reach the surface layers, where they might contaminate fish and other marine organisms. Measurements and modelling of the transport and dispersion of radionuclides was accordingly carried out at the International Atomic Energy Agency's Marine Environment Laboratory in Monaco (Cetina et al., 2000).
The three-dimensional mathematical model PCFLOW3D has been developed and specifically adapted to simulate the hydrodynamic circulation and the transport and dispersion of radioactive pollutants in lakes, coastal seas, and oceans.
3.2. Specific adaptations of the PCFLOW3D model
The circulation in the very deep Sea of Japan is clearly three-dimensional. For the simulation of the transport and dispersion of radioactive pollutants we adapted and applied an existing three-dimensional (3D) numerical model PCFLOW3D that was developed at the University of Ljubljana. As is shown later on, the initial results of the hydrodynamic simulations have shown specific disagreement with the observed surface and subsurface circulation (Cetina et al., 1996). To improve their accuracy, we had to include simulation of the topographic stress, also named topostress or Neptune effect in the model. The phenomenon is described by Holloway (1992) and Eby and Holloway (1994) as being due to the interaction of eddies with bottom topography, which can exert large systematic forces on the mean circulation. A simple representation of topostress after Eby and Holloway (1994) is used. The model's calibration and determination of the value of the topostress length parameter L for the Sea of Japan is described in Section 3.5.
Mass-transport module: The advection-diffusion equation is solved in the mass transport (submodel. As the simulations had to be run over several decades, the finite volume numerical scheme was considered not to be accurate enough regarding numerical diffusion. The Lagrangian-based particle tracking method was used (Sirca,
1992), which has no numerical diffusion in the traditional sense, as it directly describes the particles' trajectories. As is shown in Section 3.6 applying the method showed its very good characteristics.
Simulation methodology: Due to the necessary very long simulation time (up to 40 years), it was not feasible to run simulations in "real time", i.e. to let unsteady simulations continue over years, since all the forcing factors were time-dependent. The usual procedure is to do time averaging. The databases (NOAA, 1994) provide the temperature and salinity fields for each annual season. Thus we have made hydrodynamic simulations for the four seasons of the year, with all the input data (T-S fields, winds, inflow/outflow currents) being averaged over each season. For each season simulations were run over a period of 60 to 90 days, when a final, nearly steady-state velocity field was obtained.
The four hydrodynamic velocity fields were computed in the following manner. To simulate the thermohaline forcing, which is the most important forcing factor (see section 3.3), there are two possible methodologies: (1) to directly simulate the atmospheric heating and cooling throughout the year and their influence on water circulation (as e.g. in Zavatarelli and Mellor (1995) for the Mediterranean Sea); this is called the prognostic mode, or (2) to use measured data of temperature and salinity distribution in the sea and to force the hydrodynamic circulation with it in the so-called diagnostic mode.
For several reasons we decided to adopt the latter one, the diagnostic mode, with one of the reasons being much shorter computational time.
During the simulations it was proven that the basic features of pollutant dispersion are simulated accurately enough (within the possible order of accuracy of the whole study) only with the winter and summer velocity fields. Therefore the computed hydrodynamic velocity fields for these two seasons were saved on files to be used in the mass-transport model for simulations of the transport and dispersion of radioactive pollutants.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana
Longitude East
Figure 5: a) Topography of the Sea of Japan with the observed pattern of surface currents (Adapted from Yoon, 1991). EKWC: East Korean Warm Current; NKCC: North Korean Cold Current; NSB: Nearshore Branch. Dumping sites No. 6, 9, and 10 are marked with squares. Ml to M3 are locations of current measurements at the bottom of Takematsu et al. (1995); b) Simulated velocities in the surface layer (0 to 10 m) for summer conditions (no wind) with main trajectory lines - with topostress taken into account.
Slika 5: a) Topografija Japonskega morja z opazovanimi površinskimi tokovi (prirejeno po Yoon, 1991). EKWC: Vzhodnokorejski topli tok; NKCC: Severnokorejski hladni tok; NSB: Priobalni tok. Lokacije odlaganja radioaktivnih onesnažil št. 6, 9 in 10 so označene s kvadrati. Ml do M3 so lokacije meritve tokov Takematsu et al. (1995); b) Simulirane hitrosti tokov v površinskem sloju (od 0 do 10 m) za poletne razmere (brez vetra) z glavnimi trajektorijami - topografska napetost je upoštevana.
3.3. Forcing factors in the Sea of Japan's circulation
As the radioactive waste was deposited on the seabed, the bottom currents are the most important for their transport and dispersion. However, due to the great difficulty of taking measurements at a depth of 3000 to 4000 m, there exist only a few measurements of the bottom circulation that could be used for the model calibration and verification. On the contrary, there are sufficient measurements and observations of the surface circulation. To enable comparison of the model results with the available surface circulation pattern the forcing
factors influencing the surface circulation were also taken into account. The following forcing factors were analyzed.
(1) Thermohaline forcing is caused mainly by the intensive cooling of the northern part of the Sea of Japan/East Sea in winter. The maximum temperature difference between the northern and the southern parts of the Sea in winter can reach 14 degrees Celsius. The water density in the NE part of the sea is also increased by ice formation, which causes salt release. Some measurements with chemical and radioactive tracers have shown that this forcing is of primary importance (Tsunogai et
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
al., 1993). The T-S (temperature-salinity) data were taken from the World Atlas (NOAA, 1994), and the values were annotated to all the control volumes for winter and summer conditions.
(2)	The effect of wind on the surface circulation is important, and its influence on bottom currents should also be tested. The winds were therefore accounted for in the simulations.
(3)	Inflow surface currents: The Sea of Japan/East Sea is relatively enclosed, and is connected to the adjacent seas through four straits (Figure 5), of which the Tsushima-Korea Strait between Japan and Korea is the most important. It is about 130 m deep and 170 km wide. Figure 5a shows a schematic surface current chart as obtained by several measurements and observations (from Yoon, 1991). The Tsushima inflow current, which enters through the Tsushima-Korea Strait, brings an important discharge into the Sea of Japan/East Sea. Most measurements from various authors (Yoon, 1991) show that the outflow is about 75% through the Tsugaru Strait and 25% through the Soya Strait; this was also taken into account in the simulations. Mamiya Strait is very shallow (5 m) and has a negligible influence on the global circulation.
Tidal Currents are very weak in this sea, and river inflow and precipitation are low; they all have a negligible influence on global circulation in the Sea and were not taken into account in the simulations.
3.4. Basic data for the simulations
Topography and numerical grid. The computational region was divided into 60 * 72 control volumes. The space step Dx = 15 minutes longitude (from approx. 17.0 km at the southern boundary, to 24.2 km at the northern) and Dy = 15 min latitude (27.8 km approx.). As the computational region extends over 17° 30' in the N-S direction, the Cartesian co-ordinate system was replaced by the polar co-ordinate system, which can relatively well approximate the actual spherical geographical system. The model also takes into account the variable Coriolis parameter.
In the Z-direction (along the depth) the field was divided into 15 active layers; their thickness is smaller at the surface (Dz = 10 m) to more accurately take into account the thermal and salinity stratification and the influence of the wind stress, and the grid is again denser near the bottom (Dz = 200 m) to better simulate the dispersion of radioactive pollutants near the bottom. In between the max. thickness of some layers is 450 m. A maximum depth of 3500 m was supposed in the simulations.
Wind. The data on the wind field over the Sea of Japan/East Sea were obtained from Levitus (1982). Approximate mean winds for winter and summer seasons were evaluated, together with corresponding temperature and salinity fields. Spatially uniform wind was taken into account in all the cases.
Time step. In the hydrodynamic part of the simulation, the time step was limited to Dt = 180 seconds due to numerical stability. In the masstransport sub-model (simulations of dispersion of radioactive pollutants by the particle tracking method), the value of Dt = 0.5 day was used.
Eddy viscosity and diffusivity. For the horizontal eddy viscosity the value of Nh = 1500 m2/s was used in the final computations, and the same value for horizontal diffusivity Dh. The vertical eddy viscosity Nv = 10-4 m2/s was taken into account. The influence of stratification on Nv and Dv, expressed by Richardson number, was accounted for in all the simulations.
For long-term simulations of the transportdispersion of radioactive pollutants (by the particle tracking method) the proper value of the vertical diffusivity Dv is very essential. As there was no possibility of calibration, Dv was determined by a thorough study of available literature. The value of Dv = 4*10-5 m2/s was finally used in the simulations.
As already stated, the final simulations of the transport of radioactive pollutants with the masstransport module were executed with only winter and summer velocity fields due to not very detailed input data, and due to the fact that the circulation in all the four seasons has common basic features.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
3.5.	Calibration and results
Only measurements of surface currents were available for calibrating the numerical model (Figure 5a). A series of simulations were carried out to get acceptable agreement. As the first results did not show acceptable agreement, later on we included the topographic stress into the model. We have tested in detail the influence of this parameter, taking into account the various values of the topostress length parameter L. For this Sea Holloway et al. (1992) determined this value to be L = 6 km, but in our simulations this value seemed to be too high, and the value of L = 4 km was used in the end. Figure 5b shows the finally calculated surface currents (winter conditions) in comparison with the measurements (Figure 5a). The agreement is at least acceptable. The only available direct measurements of bottom currents were published by Takematsu et al. (1995). The measurements were taken by current meters, attached to three moorings, at depths of 1000, 2000, 2500, and 3000 m.
Comparison of the simulated and measured velocities can only be qualitative, as the simulations were run with forcing parameters, averaged over each season, while the measured velocities were obtained from one specific measurement (August 1993 to July 1994). Maximum measured horizontal velocities in the bottom layers, even at the depth of 3000 m, are surprisingly high, up to 10 cm/s, mainly in the spring. It is interesting that the simulations have also shown maximum horizontal velocity components up to 9 cm/s in winter and spring conditions, while in other seasons they are up to 4 cm/s.
3.6.	Transport-dispersion of radioactive pollutants
On the basis of the computed winter and summer hydrodynamic circulation, the transport and dispersion of radioactive pollutants was simulated using the mass-transport sub-model and the particle tracking method. As previously mentioned, the simulation scenario was the continuous release of 1 TBq over 90 days at dumping location No. 9 (see Figure 5a), with no decay taken into account.
Fifty particles were released at each time interval of 0.5 day, which results in a total number of 9000 particles in 90 days.
Figures 6a and 6b show the location of particles in Layer 12, which is most interesting from the point of view of the possible contamination of fish, as it is located at a depth of 100 to 180 m, the deepest layer where there is fishing. Results are shown after 5 years and then after 30 years of simulation, when approximately the maximum concentrations are attained. It can be seen that the activity is dispersed almost all over the Sea of Japan/East Sea in this near-surface layer, but the strongest contamination would be in the region north of the central part. The maximum values of calculated concentration, attained at some locations in this layer, are of the order 10-3 Bq/m3.
Figure 7 is an attempt to present the path of one specific particle in 3 dimensions. The first particles (causing very small concentrations) reach the surface layers already after about 3 years. For the assumed input data, the highest average concentration in the upper 180 m is about 6*10-5 Bq/m3 and it is attained 30 years after the spill of the radioactive waste at the bottom.
3.7. Conclusions - Case B
The following conclusions can be drawn from the study:
•	Besides the study on the possible radionuclide contamination, the present research is also an attempt to contribute some new information on the hydrodynamic and pollutant dispersion processes in the Sea of Japan/East Sea, and on the effect of topographic stress in modelling ocean circulation.
•	Comparison of the simulated velocity fields with observed surface and subsurface currents and with some measurements of bottom currents in the Sea of Japan has shown relatively good verification of the 3D hydrodynamic model PCFLOW3D. The simulated hydrodynamic field in the Sea of Japan is in visibly closer agreement with the observed surface circulation when the topographic stress is taken into account in the modelling.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
Figure 6: a) Simulated transport - dispersion of particles, with topostress taken into account (L = 4 km). 9000 particles released in 90 days from dumping site No. 9. Layer 12 (depth 100-180 m), after 5 years; b) Simulated transport - dispersion of particles, with topostress taken into account (L = 4 km). 9000 particles released in 90 days from dumping site No. 9. Layer 12 (depth 100-180 m), after 30 years.
Slika 6: a) Simulirani razpored delcev z upoštevanjem topografske napetosti (L = 4 km). Iz odlagališča št. 9 je izpuščeno 9000 delcev v 90 dneh. Sloj 12 (na globini 100-180 m), po petih letih simulacije; b) Simulirani razpored delcev z upoštevanjem topografske napetosti (L = 4 km). Iz odlagališča št. 9 je izpuščeno 9000 delcev v 90 dneh. Sloj 12 (na globini 100-180 m), po 30 letih simulacije.
Figure 7: 3D view of 12 years' trajectory of a randomly chosen particle. Start of the simulation is at the beginning of the winter season.
Slika 7: 3D-pogled trajektorije naključno izbranega delca v 12 letih. Začetek simulacije je na začetku zimske sezone.
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
• The final results, i.e. the transport/dispersion of radionuclides that could possibly leak from the containers at the dumping sites on the sea floor, show the following features: the first, very small amount of contamination would reach the surface layers (above 180 m depth, where fishing is possible) after about 3 years. The maximum concentrations, which would be about 10-2 Bq/m3 for the supposed leaking scenario, would be attained after 30 years, but after 20 years already 96% of the maximum concentration would be attained. The radionuclide contamination, however, would everywhere be about two orders of magnitude smaller than the present background values. Therefore, the estimated radiation doses to the local population will be negligible.
4. Final remarks
The two described cases show that the PCFLOW3D model can simulate hydrodynamic circulation (currents) in very different conditions: In Case A (Yatsushiro Sea) the main forcing factor is extremely strong tidal currents, while in Case B (Sea of Japan) the main forcing factor is thermohaline forcing. As the model has been verified for the conditions of strong winds and river inflows in some other cases (Gulf of Trieste, Rajar et al., 2000), the model can be considered to be a good quantitative tool for the further simulation of the transport/dispersion/fate of various pollutants.
References
Bombač, M. (2014). Matematično in fizično modeliranje toka v rekah s poplavnimi področji. (Mathematical and Physical Modeling of Flow in Rivers with Flood Plains). Doctoral Thesis, University of Ljubljana, Faculty of Civil and Geodetic Engineering (in Slovene).
Četina, M., Rajar, R., Širca, A. (1996). Development and Implementation of a Computer Code for Modelling the Dispersion of Radioactive Pollutants around Dumping Sites, Research project for the Intern. Atomic Energy Agency, M.E.L., Monaco, Contract No. 8847/RB/TC, 48 pp.
Četina, M., Rajar, R., Povinec, P. (2000). Modelling of Circulation and Dispersion of Radioactive Pollutants in the Japan Sea, Oceanologica Acta 23(7): 819-836.
Četina, M., Rajar, R, Yano, S., Tada, A., Nakamura, T., Akagi, H. (2003). Hydrodynamic Simulations of The Yatsushiro Sea, Japan. Congress of the International Association for Hydraulic Research, Thesaloniki, Greece.
DHI (2012). MIKE 21 FLOW MODEL, Hydrodynamic Module, User Guide. Denmark, Danish Hydraulic Institute.
Eby, M., Holloway, G. (1994). Sensitivity of a Large-Scale Ocean Model to a Parameterization of Topographic Stress, J. Phys. Oceanogr. 24: 2577-2588.
Holloway, G. (1992). Representing Topographic Stress for Large-Scale Ocean Models, J. Phys. Oceanogr. 22: 1033-1046.
Hydrographic Department (1978). Maritime Safety Agency: Charts of Tidal Currents in Simabara Wan and Yatsushiro Kai. No. 6217, Tokyo, Japan(in Japanese).
Levitus. S. (1982). Climatological Atlas of the World Ocean. NOAA Profess. Paper 13, U.S. Dept. of Commerce, 344 p.
NOAA, U.S. Department of Commerce (1994). World Ocean Atlas, CD No. 1: Objectively Analysed Temperature Fields, CD No. 2: Objectively Analysed Salinity Fields.
POM (1991). A User's Manual for the Princeton Numerical Ocean Model, O'Connor, W., P., Institute for Naval Oceanography, Stennis Space Center, MS, USA.
Rajar, R., Četina, M. (1997). Hydrodynamic and Water Quality Modelling, Ecological Modeling 101: 195-207.
Rajar, R. (1998). Transport of Sediments from Mururoa Lagoon into Ocean. Chapter 3.2 in: The radiological Situation at the Atolls of Mururoa and Fangataufa : Main report, Intern. Atomic Energy Agency, Vienna, 177-179.
Rajar, R., Žagar, D., Širca, A., Horvat, M. (2000). Three-Dimensional Modelling of Mercury Cycling in the Gulf of Trieste, Sci. TotalEnvir. 260: 109-123.
Rajar, R., Žagar, D., Širca, A., Akagi, H., Yano, S., Tomiyasu, T., Horvat, M. (2002). Modelling of Toxic Mercury in Coastal Seas - Case Study: Yatsushiro Sea, Japan, Intern. Envirosoft, Bergen, Norway.
Širca, A. (1992). Modelling of Pollutant Transport with the Particle Tracking Method. Thesis, M.Sc. University
Rajar R., Četina M.: Calibration of Three-dimensional Model PCFLOW3D - A Comparison with Measurements from the Yatsushiro Sea and the Sea of Japan - Umerjanje tridimenzionalnega modela PCFLOW3D - Primerjava z
meritvami v morju Yatsushiro in Japonskem morju _Acta hydrotechnica 30/52 (2017), 37-50, Ljubljana_
of Ljubljana, Faculty of Civil and Geodetic Engineering, 140 p (in Slovene).
Takematsu, M., Ostrovskii, A. G., Kitamura T. (1995). Current features in the Japan Sea Proper Water, personal communication.
Tsunogai, S., Watanabe, Y. W., Harada, K., Watanabe, S., Saito, S., Nakajima, M. (1993). Dynamics of the Japan Sea Deep Water Studied with Chemical and Radiochemical Tracers, Deep Ocean Circul., Phys. and Chem. Aspects, Elsevier, 105-119.
Yoon, J. H. (1991). The Seasonal Variation of the East Korean Warm Current. Reports Instit. for Appl. Mechanics, 38, 23-36.
Zavatarelli, M., Mellor, G. (1995). A Numerical Study of the Mediterranean Sea Circulation, J. Phys. Oceanogr. 25: 1384-1414.