A Mixed Optimal Control Approach for Upstream Fish Migration

This paper proposes a simple mathematical model forupstream fish migration along rivers. The model describes the fish migration along a river based on a mixed optimal control approach having swimming velocity, school size, and stopping time of migration as control variables. The optimization problem reduces to a variational inequality. Its explicit “viscosity” solution is presented with the dependence of the fish migration on river environment. To prove uniqueness of the solution to the variational inequality requires a constructive argument not based on the conventional theorems. A novel finite difference scheme for solving the variational inequality is also proposed with its convergence results. An application example of the model discusses the upstream migration of Plecoglossus altivelis (Ayu) in Japan, which evaluates the dependence of the fish migration on the habitat quality and provides recommendations for managing river environment. This is an interdisciplinary research between environmental and mathematical fields.


INTRODUCTION
Comprehension and assessment of fish migration are necessary for sustainable development and management of water environmental systems worldwide, examples include those in Asia [1], Europe [2], and North America [3].Anthropogenic degradations of water environment, such as dam and weir constructions and water quality changes, have significantly affected the fish migration through changes in hydrological and hydraulic characteristics of their habitats and migration routes.Fish migration in rivers is highly important from ecological, fishery, and cultural viewpoints.Impacts of environmental changes on life histories of migratory fishes have been assessed based on the willingness to pay [4], ecological services [5], and ecosystems stability [6].The fish migration has been severely affected by physical barriers such as cascading dams [7] and huge hydroelectric dams [8], and the associated environmental changes [9].Evaluation of impacts of environmental changes on the fish migration is the most urgent issue in the research of above-mentioned areas.Mathematical models to evaluate fish migration would help establish a management way of river environment and ecosystems involving migratory fish, however, such research is still rare.
Conventional models for fish migration are mechanistic multi-agent models to simulate detailed interactions among individuals.Such models consider the multi-scale nature of fish school dynamics [10], food-web dynamics [11], and body growth [12].They are necessary for evaluation of microscopic dynamics such as swimming behavior around hydraulic structures.However, they would be inefficient and often too complicated for analyzing macroscopic dynamics such as fish migration along rivers with complex connectivity [13], those in rivers with many hydraulic structures [14], and those in wide brackish waters [15].For approaching the macroscopic fish migration, conceptual mathematical models that track the longitudinal migration dynamics can potentially serve as simpler and more efficient alternatives to the multi-agent models.
If possible, such a model should be biologically plausible as well as mathematically rigorous, however, in most cases, mathematical analysis and practical application have not been carried out simultaneously.This is the motivation of this paper where both the analysis and application are addressed.
The objectives of this paper are to propose and analyze a simple, conceptual mathematical model for upstream fish migration along a 1-D river.The upstream migration is a key biological process for assessing the ecological dynamics of many migratory fishes like salmonids [12,16].The present model can be applied to diadromous and anadromous fish species that have upstream migrations along rivers in their life histories.Our model handles a fish school as a whole: the decision-maker as a synchronized group of individuals.The proposed model describes the fish migration based on a mixed optimal control approach having swimming velocity, school size, and stopping time of migration as control variables.The fish school migrates upstream along the river and stops its migration at some point, namely, a habitat.A performance index containing the cost of migration and the benefit of finding a high-quality habitat is presented, which is maximized by the fish school.This optimization problem ultimately reduces to a Variational Inequality (VI) that governs the optimal migration strategy, which is the equation to be solved.Our approach thus provides a new mathematical description of the fish migration based on the mixed optimal control theory, however, its basic structure is quite natural where the migration arises from decision-making processes by fish.
The present VI admits a "viscosity" solution [17], which is an appropriate weak solution to a wide class of VI's both with [17] and without the diffusion terms [18].Our VI is of the latter type.Its viscosity solution under a simplified case is derived and biological and ecological implications of the solution are discussed with an emphasis on the dependence of the fish migration on river environment.To prove the unique existence of a viscosity solution to the VI requires a constructive argument not relying on the conventional comparison theorems because of its discontinuity.Proving the uniqueness amounts to saying that the exact solution is mathematically rigorous, which is a requirement for biologically plausible solutions.A finite difference scheme for solving the VI under realistic conditions is also proposed with its convergence results.An application example of the model is presented to evaluate upstream migration of the major inland fishery resource Plecoglossus altivelis (P.altivelis, Ayu) in Japan [19].Yoshioka et al. [20] numerically simulated their horizontally 2-D swimming behavior around an existing weir.This paper is based on a mathematical approach for dealing with the upstream fish migration problem, however, it turns out to be an effective candidate for describing a wider range of biological phenomenon from the new standpoint.A mathematical analysis would be inevitable for comprehension and assessment of properties of newly developed mathematical models.The results should then be discussed from both theoretical and practical point of views, which requires an interdisciplinary research framework between mathematics and applied research area.This paper serves as such an example focusing on the problem of fish migration, where the researchers in life and environmental science and mathematics co-work.The analysis results of this paper demonstrate that the concept of viscosity solution is a powerful mathematical tool for analysing fish migration.The present model can be used for identifying the potential area of migration in a river, which would be useful at least from the viewpoint of fishery resources management.This paper is highly mathematical and presents several mathematical analysis results as the main results, however, being different from pure mathematical papers, this paper provides information about how the present model can be used in practice and what kind of result can be obtained from it.The mathematical and numerical analysis results here extend that of the previous research [21].In this paper, more detailed discussion and analysis of the mathematical model and the numerical scheme are presented.
The rest of this paper is organized as follows.The second section presents the mixed optimal control model proposed in this paper and derives its associated VI: the main problem.The exact solution to the VI and its proof of uniqueness, which is a non-trivial and biologically important issue, is presented in this section as well.The third section is devoted to development and verification of a numerical scheme to discretize the VI.The fourth section gives a brief application result of the present mathematical model.The last section gives a summary of this paper and presents future perspectives of our research.

MATHEMATICAL MODEL
This section presents our model and derives and analyses the associated VI.The mathematical notions like function spaces are found in the textbook [22].

Migration dynamics
An Ordinary Differential Equation (ODE) that governs the longitudinal movement of a fish school is presented.The 1-D domain D = (0, L) with the length L > 0 is considered as a river reach bounded by physical barriers placed at the boundary D ∂ that contains the upstream-and downstream-ends x = 0, L. The flow velocity : V D → R is positive and Lipschitz continuous in D .The fish school is considered as a synchronized group, which is seen as a moving point in D. The position of the school at the time t ≥ 0 is denoted as → where the time that Xt firstly hits D ∂ is denoted as T. The swimming of the fish school is considered as a control variable.
Here, umax is the maximum sustained swimming speed, which is assumed to be sufficiently large such that ( ) . The ODE that governs the longitudinal movement of the fish school until it first hits the boundary D ∂ is set as: ( )

Performance index
The performance indexed to be maximized by the fish school, the decision-maker, is presented.It is assumed that all the individuals in the school share a common A Mixed Optimal Control Approach for Upstream ... represents the profit by terminating the migration, which is assumed to be non-negative.The performance index J is set as: with the discount rate δ ≥ 0. Conceptually, δ represents impatience of the fish school: larger δ results in earlier stopping of migration.It would be natural to consider that the existence of a performance index to be maximized during a fish migration process is a result of evolutionary processes.Also, δ can be larger for the migration subject to higher predation pressure: namely higher mortality.In eq. ( 2), the function f ≥ 0 with f (0) = 0 is the hydraulic cost of solo swimming per unit time as a smooth, non-negative, convex, and even function.This f can be identified from laboratory or field experiments on the swimming behavior of individual fishes along a current as demonstrated later.The coefficient a > 0 represents the discount of the hydrodynamic cost by forming a school, which is motivated by the theoretical and experimental results that increasing the size N effectively reduces the hydrodynamic cost per individual [23].The probabilistic model theoretically predicts [24]: On the other hand, the coefficient b > 0 conceptually represents the non-hydraulic cost by forming a school per unit time.Schooling would have many aspects of benefits, such as improvement of navigational performance, hearing perception, and foraging efficiency.While at the same time, schooling would negatively affect passage efficiency, information transmission among individuals, and competitions among the individuals [25].Being different from the sophisticated multi-agent models [26], a simple strategy is employed in this paper.This is because these effects are possibly interacting and difficult to mathematically describe in detail.Therefore, the lumped approach [27] is employed: with k > 0 where its sophisticated parameterization remains as a future topic.This term is necessary so that a fish school is created.For example, removing this term from the performance index J leads to an unreasonable result that forming a fish school with the infinitely large fish school is optimal.The value function : D Φ → R is the maximized performance index defined as: From a biological viewpoint, the value function is an index to be maximized by the fish school through choosing a migration strategy.The optimizers (u, N, τ) to give the supremum of eq. ( 5) are expressed as (uopt, Nopt, τopt), which are functions of x.Finding the optimizers is the goal of the present model.

Variational Inequality
The dynamic programming principle (Chapter III.4 of [18]) leads to the VI as a governing differential equation of the value function Φ: where : H D × × → R R R is the function called Hamiltonian, and the function : , inf The derivative of Φ with respect to x is denoted as Φ'.Hereafter, the function f is assumed to lead to the inequality: where c1, c2 > 0 are some constants depending on neither x nor p. Also, the condition k > m is assumed to make convex be (fa −1 + b) (u, N).The latter assumption leads to a concave F(x, p) for p ∈ R , and the equation F(x, p) has the two solutions p = p − (x), 0 with p − (x) < 0 for x D ∈ .At least two examples of f that comply with the conditions mentioned above are found, both of which are used in this paper: f1(u) = |u| n+1 with kn > m (condition of forming a school) and ( ) 8) is a key for the unique solvability of eq. ( 6).

Definition of viscosity solutions
Viscosity solutions are appropriate candidates of weak (non-classical) solutions to degenerate elliptic differential equations like the eq.( 6).This VI does not have classical, continuously differentiable solutions in general but has viscosity solutions that are not differentiable and possibly discontinuous.Handling the eq.( 6) requires the concept of viscosity solutions even under simplified conditions as shown in the next sub-section.For the sake of brevity of analysis, the boundary operator : moreover, the operators * * , : and: for the sake of brevity of descriptions.The upper-and lower semi-continuous envelops of : D Φ → R are denoted as Φ * and Φ*, respectively [28].They are defined as: and: respectively, where Bε(x) Following the literature for viscosity solutions [27], the definition of viscosity solutions to the eq.( 6) is stated as follows.
is a viscosity sub-solution (super-solution) to the eq.( 6) at x D ∈ if: for any test functions 1 ( ) is a viscosity solution if it is a viscosity sub-solution as well as a viscosity super-solution.
Note that a classical, sufficiently smooth solution is a viscosity solution.
By definition, a viscosity solution maybe discontinuous at D ∂ .Any smooth, classical solution is a viscosity solution, meaning that the concept of viscosity solution is indeed a weaker than that of the classical solution.An immediate consequence of Definition 1 is that any viscosity solution Φ satisfies Φ ≥ α in D. In fact, if Φ < α at some x D ∈ , then Then, Φ cannot be a viscosity super-solution.This is a contradiction.

Exact viscosity solution
An exact viscosity solution to the eq.( 6) is derived for specific functional forms of the coefficients, which can give biological and ecological implications of the proposed mathematical model.Assume the uniform flow condition V = const < umax and: with the constants A, B, C > 0. This α is motivated by the fact that the quality of habitat measured by the amount of food, such as diatoms, would sharply increase along a river [29].A formal exact solution to the eq.( 6) on D is then found as Figures 1a and 1b: ( ) ( ) ( ) In eq. ( 16), p − is a negative constant determined as a solution of F = 0.The two graphs y = γ(x) and y = α(x) intersect at x = y1 such that y0 < y1 ≤ L if: ( ) Therefore, when eq. ( 18) is satisfied, the found Φ in eq. ( 16) is rewritten for x > y0 as: Similarly, when eq. ( 18) is not satisfied, Φ in eq. ( 16) is rewritten for x > y0 as: The optimal controls can be determined in the sub-domain of D where α(x) < Φ(x).Assume that the minimum of the second term of eq. ( 7) is achieved by the internal solutions ( ) Then, they are positive constants and satisfy: and: To prove that the exact solution in eq. ( 16) is a viscosity solution is not a trivial issue because this Φ is not smooth at x = y1.Since the exact solution satisfies the eq.( 6) where the solution is continuously differentiable, it has to be examined with the condition of viscosity solution at which it is not differentiable or discontinuous.After elementary calculations, it turns out that eq. ( 16) is a viscosity solution to the eq.( 6) by examining it against the condition of viscosity super-and sub-solutions on D .

Implications of the exact viscosity solutions
Practical implications of the exact solution of eq. ( 16) are presented focusing on its dependence on the parameters and coefficients.The condition ω ≥ 1 implies that the cost of migration is sufficiently small and/or the benefit of finding a high-quality habitat is sufficiently high.The fish school approaches the upstream-end x = 0 when A is sufficiently high.The school does not migrate upstream (τopt = 0) for all x D ∈ when ω < 1.Therefore, the exact solution in eq. ( 16) implies that the fish school would not migrate upstream along the river if its environment is severely degraded.In particular, y0 is increasing with respect to A, showing that degradation of the river environment results in more downstream termination of the fish migration.Therefore, the shape and the magnitude of the habitat quality α are important for qualitative understanding of the fish migration.This point is investigated in the fourth section as well through a model application.
For the cost function f1 presented in the previous section, assuming the condition kn > m yields that uopt and Nopt are increasing with respect to the flow speed V, which is consistent with the experimental results [30].Similar statement holds for f2.The second of eq. ( 21) states that the swimming speed is increasing with respect to the school size, which is also consistent with the experimental results [30].

Unique solvability
In general, there is no guarantee that the found exact solution is the unique viscosity solution to the eq.( 6).This is not only a mathematical issue but also a practical issue since the viscosity solution should be the plausible solution to the eq.( 6).The issue of unique solvability is discussed below for the cases δ > 0 and δ = 0 separately.In this section, α of eq. ( 15) is employed.
For δ > 0, the eq.( 6) is a unique continuous viscosity solution by the comparison theorem (Theorem 3.17 of [28]) since the eq.( 6) satisfies the continuity and ellipticity conditions {eqs.(3.1), (3.24), and (3.35) of [28]}.The uniqueness of discontinuous viscosity solutions for δ > 0 are unclear at this stage, but it is expected that this subject can be dealt with based on the constructive argument presented below.
For δ = 0, the above comparison argument cannot be used since it requires δ > 0. This is the major difficulty for dealing with the case with δ = 0.The remaining part of this sub-section presents the idea and plan of a constructive proof of unique solvability of the eq.( 6) for δ = 0 in the viscosity sense for the exactly-solvable problem in the previous sub-section.The constructive proof follows a step-by-step manner presented in Propositions 1 through 4. The proof is a bit long, but is presented in this paper for the sake of its self-containedness.
For the sake of brevity of descriptions, the sets D1 and D2 are defined as: It is assumed that viscosity solutions do not have dense non-differentiable points in D. This assumption has a technical aspect, but viscosity solutions having dense non-differentiable points in D can be very irregular and not reasonable from biological viewpoints.
Proposition 1: Φ'(x) = p − in D1.The Rademacher's Theorem (Theorem B.12 [28]) shows that Φ is differentiable almost everywhere in D1 since its viscosity sub-differential is locally bounded in D1 (Proposition 4 of [31]) because of the property in eq. ( 8).This leads to that Φ is decreasing in D since F(x, p) = 0 has the two solutions p = p − (x), 0 with p − < 0 for x D ∈ and α is strictly decreasing.It is shown that a viscosity solution Φ has at most a finite number of non-differentiable points in D1.Non-differentiability of Φ at The case (a) violates the condition for viscosity sub-solutions, and therefore not allowed.On the other hand, the case (b) complies with the condition for viscosity solutions, however, it leads to Φ'(L − 0) = 0 and Φ(L − 0) > 0 due to the strictly decreasing property of α (α' < 0 in D), which violates the condition for viscosity sub-solutions at the boundary x = L. Therefore, a viscosity solution Φ does not have a non-differentiable point in D1.The above result implies Φ' = p − or Φ' = 0 in D1, but the latter again leads to Φ'(L − 0) = 0 and Φ(L − 0) > 0, which violates the condition for viscosity sub-solutions at x = L. Therefore, Φ' = p − in D1.In addition, then Φ is strictly decreasing in D.

NUMERICAL SCHEME
A numerical scheme for solving eq. ( 6) is proposed and verified.

Discretization
Solutions to the eq.( 6) may not be explicitly found in applications, which motivates us to develop a numerical scheme that can approximate its solutions.For this purpose, a finite difference scheme based on exact viscosity solutions to local two-point boundary value problems is presented for solving the eq.( 6) under general conditions.This kind of discretization technique utilizing local exact solutions is called fitting technique and has been applied to elliptic and parabolic problems [32].The domain D is discretized into The discretization procedure below focuses on a vertex xi (1 ≤ i ≤ I − 1) without the loss of generality.
The finite difference scheme for a linear problem is firstly explained and is extended so that the present VI is handled.Consider the linear advection-decay equation: with a piece-wise continuous g.W (≠ 0) and g are approximated in each cell and those in Ci are denoted with the subscript i as Wi and gi, respectively.A finite difference approximation of eq. ( 24) is proposed as: ,0 1 ,1 where Fi,j is expressed with Φi+j−1, Φi+j+1, Wi+j, and gi+j.Consider the local linear two-point boundary value problem in Ci to find : ( ) Its unique viscosity solution in Ci is found as: The present scheme evaluates Fi,j with Ψi as: ,0 1 0 0 1 The boundary conditions are directly specified at x0 = 0 and xI = L. Assembling eq. ( 25) for 1 ≤ i ≤ I -1 with the boundary conditions leads to a system of linear equation with a coefficient matrix whose inverse is positive-definite, showing that the present scheme is unconditionally stable.The scheme reduces to the classical first-order upwind finite difference scheme as 0 R → + .For solving the eq.( 6), a penalty method [33] is employed to reduce a VI to a differential equation that is easier to numerically handle.The penalized counterpart of the eq.( 6) is given as: with the penalty parameter λ > 0, so that it is expressed like eq. ( 24).The original eq.( 6) is formally recoveredwhen λ → +∞ , which is implemented numerically in the penalty method.In the present scheme, the first and second terms of eq. ( 29) are advection-decay terms having solution-dependent coefficients.These coefficients are determined by optimizing the "min" term of eq. ( 7) in each cell.The third term is discretized at each vertex as λ max {αi -Φi, 0}.A standard fixed-point iteration is used for solving eq. ( 29).
The iteration is terminated when the difference of old and updated numerical solutions at each vertex becomes smaller than 10 −8 , which is a sufficiently small value.

Computational accuracy
The present finite difference scheme is degenerate elliptic (Definition 2), Lipschitz continuous (Definition 3), and proper (Definition 6) in the sense of Oberman [34].Therefore, for the linear problem of eq. ( 24), the scheme admits a unique numerical solution (Theorem 8 [34]).In addition, the scheme is monotone and non-expansive in the l ∞ -norm since it is degenerate elliptic (Theorem 3 [34]).The same statements hold for eq.( 29).
Convergence of the scheme for our VI is numerically examined with both continuous and discontinuous viscosity solutions.The domain is uniformly discretized into I cells.For the case δ = 0, the non-oscillatory numerical solutions presented in Figures 1a and 1b demonstrate its accuracy and stability.The numerical solutions accurately capture the non-differentiable and discontinuous points of the solutions.Tables 1 and 2 summarize the computed l ∞ errors for a series of λ and I −1 .In Table 1, "Osc" means non-convergence of numerical solutions.This also applies to the other tables in this paper.The computational results suggest the efficient scaling relationship λ = O(I −1 ), which may not be a sharp relationship but useful for accurate numerical computation.Rigorous proof for this scaling result should be analyzed in future research.
For the case δ > 0, since our VI does not admit exact solutions, the scheme is examined against the following test problem with constant coefficients: which can also be cast into − a penalized form as in eq. ( 29).The viscosity solution to eq. ( 30) is: which is discontinuous at x = L when Ψ < 1, where:  In the numerical computation, the parameter values employed are V = 1, R = 2, L = 1, and A = 1 for the smooth, classical solution case and A = 10 for the non-smooth, viscosity solution case.Tables 3 and 4 summarize the computed l ∞ errors for a series of λ and I −1 .The computational results presented in this sub-section indicate the satisfactory performance of the present finite difference scheme.The error for the smooth, classical solution case is negligible.This is owing to the specialized discretization employed in the present finite difference scheme that the solutions with the exponential profiles like that in eq. ( 32) can be exactly reproduced.

DEMONSTRATIVE APPLICATION
The model is applied to an example of fish migration in Japan.

Study area
The present mathematical model and the verified finite difference scheme are applied to the evaluation of upstream fish migration in Hii River, San-in area, Japan (Figure 2, top left).This river provides indispensable water resources for human activities in this area.The total length of the mainstream and the catchment area of Hii River are 153 km and 2,070 km 2 , respectively.The river flows into the brackish Lake Shinji, which has a downstream brackish lake named Lake Nakaumi.Hydrology of Hi River are reviewed in Somura et al. [35].Both brackish lakes are key habitats for migrating aquatic species that grow and spawn in Hii River.Hii River Fisheries Cooperatives (HRFC) authorizes inland fishery of the middle to upper reaches of the river.
P. altivelis (Ayu) (Figure 3, top left) is one of the common diadromous fishes in Japan and is a major inland fishery resource in the country.The natural P. altivelis has an annual life cycle [36].During autumn, adults spawn eggs in the downstream reaches of the living river and die soon afterwards.Hatched larvae descend to coastal areas of the downstreamwater body of the river, typically a sea or a brackish lake.The larvae grow up to juveniles with feeding on zoo planktons till the next spring.The grown fishes ascend the river toward its midstream and upstream reaches where diatoms (Figure 3, top right), which are staple foods of P. altivelis, are available on the riverbed.They feed on the algae to mature till the coming autumn when they descend the river.
Recently, fish catches of P. altivelis in Japan have been dramatically decreasing due to severe degradations of the river environment and ecosystems.Installing physical barriers, such as huge dams and weirs, into river cross-sections prevents fishes from upstream migration and significantly affects downstream water environment.The population of the fish in some rivers are maintained through intensive artificial hatching of farmed juveniles during the spring season in each year by local fishery cooperatives.The artificially hatched P. altivelis also migrate toward upstream after the hatching and grow till the coming autumn, but has been thought to be unsuccessful in reproduction due to genetic reasons.This would be the main reason why the fish have to be artificially hatched in a river in each spring to maintain their population.
A major issue in releasing the fish is when and where to release farmed P. altivelis along the river and its tributaries so that the released fishes grow well and the fish catch in the river increases.In Hii River, local fishery cooperatives and residents say that the harmful attached algae Cladophora glomerata Kützing (Figure 3, bottom) are significantly growing from the just downstream point toward several km downstream reach of the huge multi-purpose dam named Obara Dam (Figure 2, top right) after its construction in 2011.This dam does not have facilities like fishways that fishes can pass through, serving as a physical barrier.According to residents and the officers of HRFC, the areas at and around which the dam was created, and the harmful algae are growing involved good fishing ground of P. altivelis.They also report that the P. altivelis after the construction of the dam are significantly smaller than those before the construction.In addition, the experimental research indicated that P. altivelis could not digest the harmful algae [37].The Hinobori Weir with the height of 11 m (Figure 2, bottom) for erosion control is the second largest physical barrier in the mainstream.The weir has a vertical-slot fishway longer than 100 m, but its passage efficiency is not clear.In each spring, many P. altivelis are artificially released into the river reach between the Hinobori Weir and Obara Dam.Currently, HRFC is faced with the problem when and where to artificially hatch farmed P. altivelis given the potential shift of the attached algae.This problem is partly approached in this paper.
subject to ( ) which turns out to give f = c3f2(u) with a constant of integration c3 that can be set as c3 = 1 without loss of generality.f2 is twice differentiable and convex in (0, umax), and its derivative diverges at u = umax.A limitation of the above uswim and f is that they cannot simulate upstream fish migration with the flow speed greater than umax.It is conjectured that, above the flow speed umax, the fish is not able to maintain prolonged swimming speed and must perform burst swimming such that the present mathematical model cannot be directly used.What is important here is that the hydraulic cost function f is derived from the swimming behavior of the fish.

Numerical simulation
The river reach from just downstream of Obara Dam (x = 0 km) to just upstream of Hinobori Weir (x = L = 13 km) is set as the domain D, which is discretized into 218 cells with 219 nodes.A steady flow field in the Hii River system involving the river reach is computed with the 1-D shallow water solver [37] specifying the Manning's friction coefficient of 0.05 s/m 1/3 .This numerical solver has already been verified against the benchmark, experimental, and real cases.The flow discharge in the reach is 12.5 to 15.0 m 3 /s, which does not significantly deviate from the observed discharge during April to May.The flow velocity in the reach is utilized as the coefficient V for the present model (Figure 4).The profiles of 0 ≤ α ≤ αmax in D with the maximum value αmax > 0 are specified on the basis of the summation of two hyperbolic-tangential functions in eq. ( 15) whose coefficients were determined from the interviews from the officers and members of HRFC (Figure 4), max α α ≈ at x D ∈ represents that diatoms are reported to be abundant at x, 0 α ≈ represents that diatoms are reported to be sparse at x, and otherwise the diatoms moderately exist.This more accurate identification of α will be carried out in our future research.The other model parameters are set as δ = 10 −6 1/s, m = 1/3, k = 2.0, d = 10 −3 m 2 /s 2 , and λ = 10 6 1/s.
Figure 5 shows the computed area in the river where the fish school terminates its migration, namely the area where Φ = α in D for a range of αmax.The results show that the area expands as αmax decreases.For moderately large αmax, there is a white area around 0.6L < x < 0.9L where the fish school ascends.The results show that hatched fishes in the downstream grey area around 0.9L < x < L would not ascend the river reach.In the other white area around 0 < x < 0.5L, the fish school swims toward the upstream-end of the potential habitat.In addition, a decrease of αmax, namely degradation of the overall A Mixed Optimal Control Approach for Upstream ... Year 2019 Volume 7, Issue 1, pp 101-121 quality of potential habitat of the fish as good fishing grounds, shows that the downstream-end of the potential habitat move toward upstream.Increasing the value of δ leads to the wider black area than that in Figure 5 as shown in Figure 6 for δ = 10 −4 1/s, and δ > 10 −3 1/s results in no upstream migration with Φ = α on D except for near x = 0. Management of river environment is implied to be crucial for fish migration and inland fisheries.The computational results suggest the following management policy for the studied river reach.Especially for the upstream part of the river reach ( 0 0.5 x L < < ), the habitat quality should be improved so that this part serves as a habitat for the fish.For example, cleaning up the riverbed to exterminate the harmful algae can an effective way for that purpose.The results also suggest that the fish should be released in the downstream reach ( 0.5L x L < < ) if the upstream part is not improved.
Figure 5.The computed area where the fish school terminates its migration along the river (lateral axis) with respect to αmax (vertical axis), the area is coloured blue, the parameter value is δ = 10 −5 l/s in this figure Figure 6.The computed area where the fish school terminates its migration along the river (lateral axis) with respect to αmax (vertical axis), the area is coloured blue, the parameter value is δ = 10 −4 1/s in this figure

CONCLUSIONS
This paper presented a mathematical model for upstream fish migration along 1-D rivers from a new, mixed optimal control-based standpoint.The VI that governs the optimal swimming velocity, school size, and stopping time of migration was derived, and its exact solution was derived under a simplified condition.It was shown that unique solvability of the VI is not a trivial issue, and was analyzed from the viewpoint of viscosity solutions.The unique solvability issue was overcome for the simplified case with the help of a constructive argument that does not rely on the conventional comparison arguments.A stable and accurate finite difference scheme based on the fitting technique was developed for numerically solving the VI.The scheme turned out to have satisfactory ability to handle the VI.A numerical application of the present model to upstream migration of P. altivelis in Hii River, Japan suggested potential downstream shift of good fishing ground of the fish due to overgrowth of the harmful attached algae on the riverbed.The results are useful for decision-making for environmental restoration of Hii River, but more useful and reliable results can be obtained if the temperature and water quality influences on the migration are considered in the model.Collecting more hydraulic data of the river and more bioenergetics data of P. altivelis is required to achieve this purpose.
Future research will evaluate migration of P. altivelis in Hii River under more realistic conditions, such as coupling of the present model with the equations of solute transport.Mathematical analysis, unique solvability of the VI, under generalized conditions is also an important research topic.Some drawbacks of the present mathematical model should also be addressed in future.For example, actual river environment can be stochastic and uncertain for the fish, which motivates us to extend the model so that the fish migration is described as a decision-making process under incomplete information.Scaling-up of the microscopic effects, such as the field of pressure around individual fishes, will be a key step toward development of a more biologically and physically plausible model.In addition, further scaling-up of the present model to a lumped 0-D population dynamics model may be useful in practical analysis.Field surveys with the local fisheries cooperatives in and around Hii River will be continued for deeper comprehensions of the river environment and ecology.The presented mathematical modelling framework would potentially serve as a core for integrated assessment of food-water nexus.

ACKNOWLEDGMENT
This research is funded by grants-in-aid for scientific research No. 17K15345, No. 16KT0018, and No. 15H06417 from the Japan Society for the Promotion of Science (JSPS) and Applied Ecology Research Grant No. 2016-02 from Water Resources Environment Center in Japan.A part of this research includes research results from the short-stay researcher projects in 2015, 2016, and 2017 of Institute of Mathematics for Industry, Kyushu University, Japan.The authors thank dr. Yumi Yoshioka of Tottori University for her support in creating a figure panel.The authors also thank all the participants of the Special Session: Math for Water, Soil, and Ecosystems in the international conference SDEWES2017 for their valuable comments and suggestions.

Year 2019 Volume 7 ,
Issue 1, pp 101-121 performance index containing the three terms: benefit of finding a high-quality habitat, hydrodynamic cost, and non-hydrodynamic cost.The population density in the fish school Nt: [0, T] → [0, +∞], referred to as the size in what follows, is a function of the time t.The stopping time of migration is denoted as τ.In what follows,

Figure 1 .
Figure 1.Exact (blue) and numerical solutions with the cell size l = 0.01 L (green), and α (black) for: continuous (A = 1) (a) and discontinuous case (A = 2) with f = f1 and n = 2 (b); the other parameters are V = 1, umax = 5, δ = 0, m = 1/3, k = 0.5, d = 0.5, B = 10, C = 20, λ = 10 6 (the numerical solutions are computed with 100 uniform cells) , the two curves z = Φ(x) and z = α smoothly contact at x = y, showing that α is concave at x = y.Since α is concave for x ≤ y because of its hyperbolic tangential shape, Φ(0) > α(0) follows.This contradicts with the boundary condition Φ(0) = α(0).Therefore, D2 does not involves singletons as disconnected subsets.Proposition 4: Identification of the solution structure.Firstly, note that D1 is a union of open sets and is expressed as 1 1,i i D D = ∪ with (possibly infinitely) a countable number of open sets D1,i.Similarly, note that D2 is a union of closed sets and is expressed as

Figure 4 .
Figure 4. Computed flow velocity V and the modelled coefficient α along the reach Yoshioka, H., et al.A Mixed Optimal Control Approach for Upstream ...

Table 1 .
Relative errors for the continuous viscosity solution with δ = 0

Table 2 .
Relative errors for the discontinuous viscosity solution with δ = 0

Table 3 .
Relative errors for the smooth, classical solution with δ > 0

Table 4 .
Relative errors for the non-smooth, viscosity solution with δ > 0 Yoshioka, H., et al.A Mixed Optimal Control Approach for Upstream ...