International Journal of Differential Equations

Volume 2015, Article ID 485860, 8 pages

http://dx.doi.org/10.1155/2015/485860

## Existence and Permanence in a Diffusive KiSS Model with Robust Numerical Simulations

Department of Mathematics and Applied Mathematics, University of the Western Cape, Private Bag X17, Bellville 7535, South Africa

Received 21 July 2015; Revised 27 November 2015; Accepted 1 December 2015

Academic Editor: Timothy R. Marchant

Copyright © 2015 Kolade M. Owolabi and Kailash C. Patidar. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

We have given an extension to the study of Kierstead, Slobodkin, and Skellam (KiSS) model. We present the theoretical results based on the survival and permanence of the species. To guarantee the long-term existence and permanence, the patch size denoted as must be greater than the critical patch size . It was also observed that the reaction-diffusion problem can be split into two parts: the linear and nonlinear terms. Hence, the use of two classical methods in space and time is permitted. We use spectral method in the area of mathematical community to remove the stiffness associated with the linear or diffusive terms. The resulting system is advanced with a modified exponential time-differencing method whose formulation was based on the fourth-order Runge-Kutta scheme. With high-order method, this extends the one-dimensional work and presents experiments for two-dimensional problem. The complexity of the dynamical model is discussed theoretically and graphically simulated to demonstrate and compare the behavior of the time-dependent density function.

#### 1. Introduction

The study of reaction-diffusion problems has gained much attention over the years; such models are largely encountered in various fields and become increasingly useful tools for field engineers, mathematicians, conservation biologists, and population ecologists. In this paper, we describe the current state of the linear Kierstead-Slobodkin [1] and Skellam [2] problem (known as KiSS model) and give possible extension to such reaction-transport problem in a nonlinear form. The first theoretical approach to the KiSS model dates back to the early 1950s, studying the dynamic of reaction-diffusion problem for a population growing on a patch of finite size (say, ) that is surrounded by a deadly environment with infinite mortality [3].

Understanding the conditions that could guarantee the extinction and persistence of species populations in large but finite domains is another focus area of this paper. Further, determination of the critical patch size ensuring the sustenance of population is another vital area to consider, and the critical patch size depends on some factor, namely, the species population in the patch, geometrical patch, boundaries type, and the reproduction kinetics of species population. In what follows, we will present the extended nonlinear KiSS model and determine its critical patch size and reproduction processes with hostile boundaries.

#### 2. Mathematical Analysis of the Main Equation

To start with, we let the critical patch size be and consider the reaction-diffusion problemon the closed interval subject to homogeneous Dirichlet boundary and initial conditionswhere is the species density, is the diffusion coefficient, is the growth rate, and is the critical exponent parameter that determines whether model (1) is linear () or nonlinear (). From the given boundary conditions, we observed that the critical patch size corresponds to the borderline between species population extinction and existence; we thus linearize (1) and (2) about , the washout state. If is stable, we have a total extinction of the species population. But if is unstable (nontrivial case), we have a state that corresponds to the persistence or survival of the species. For simplicity, we relax the nonlinear condition and let . We linearize (1) and (2) about to obtain

For linear dynamics, it is permitted to seek for solutions of the form By using this ansatz in (3), we obtain a simple diffusion equation:and the solution of (3) is given by [4, 5] where are the coefficients to be determined by the Fourier series expansion of . The unstable trivial solution could result in the growth and persistence of the species population if for some values of . The growth rates usually decreased monotonically with , and . By setting the dominant mode , we can determine the fate of species population which guarantee the total washout of the species. The critical patch size known as the KiSS size is defined by the expressionIf , as . The population is wiped out from its initial condition; no nontrivial steady state develops. Instability (bifurcation) could only arise if the patch size is increased. The state loses its stability at the point , and due to nonlinear effect begins to grow and saturates with time . So as , a steady state (1) is determined byIt is obvious that, for any value of , the trivial solution holds for (8). The interest of the present work is in the solutions for which . As a result of the spatial symmetry of (1) and (8), one expects the solution in to be symmetric about the midpoint , and it is reasonable to assume that the species density of nontrivial steady states could reach its maximum point, say at the midpoint , since at the boundaries. Multiplying (8) with and integrating the result with respect to on yield where . If we simplify (9) by separating the variables and integrate it from to with respect to , we obtain We can rescale with to obtainEquation (11) is usually called a time map, from where can be determined implicitly as a function of . We can also change the origin to ; that is, by setting , set ; we have where is the eigenvalue. In the spirit of [4], we have

By taking the limit as , we obtain the critical patch size which is similar to (7). In conclusion, we can say that the nonuniform stationary solution of reaction-diffusion equation (1) exists if and that the trivial solution coexists.

In other physical contexts, (1) is referred to as a blow-up problem [6–9]. It is a known theory from ODEs that solutions may not exist for (global existence) but could tend to for a finite time (local existence). So, similar things happen to reaction-diffusion problems. On a large but finite interval, say , problem (1) can be written in integral form:where and are subject to homogeneous Dirichlet boundary conditions on . The nonlinear term is locally Lipschitz on , and continuous functions vanished at . By argument, if is nonnegative function in , then there exists a nonnegative local solution of (15). A reasonable question to ask is, does as for all or singularity is restricted to a smaller interval? To answer this, we give the following theorem [10].

Theorem 1. *Let and be the solution of the stationary problem Let , where is chosen so that the maximum solution of (15) poses finite existence time, say . If and is sufficiently large, then exist and are finite for all . In other words, singularity (discontinuity) occurs at the point .*

*Proof. *The proof of Theorem 1 can be found in [10].

#### 3. Formulation of Adaptive Methods in Space and Time

Considerations have been given to the use of spectral methods as alternative to conventional finite differences [11], because they are capable of removing the stiffness property in the linear term of (1). Knowing that the reaction-diffusion problem can be split into the linear and nonlinear terms, the use of time and space methods is permitted [12–14]. Also, based on the known integrating factor technique, we are going to formulate the theory of spectral method in two spatial dimensions. By applying the integrating factor technique to Fourier transform systems (1), we havewhere is the double Fourier transforms of species density . In other words, To explicitly remove the adherent stiffness in the second partial derivatives, we let and set , so thatNext, we need to discretize the square domain by considering the equispaced points and in the spatial directions of and . We employ the discrete fast Fourier transform (DFFT) [12, 13, 15] to transform (19) to a system of ordinary differential equations (ODEs):where , , , and . Boundary conditions are now set at extremes of the domain. Also in the experiment, we use a square Fourier node of for . At this stage, we have converted the system to ODEs; also, the stiffness has been removed. It should be mentioned that we can now use any explicit scheme to advance in time.

By adopting the notation of order Runge-Kutta (RK) methods, with time step , we can advance from to for the ODE: where In the interest of the present paper, we apply the general explicit Runge-kutta scheme to (20) for . For brevity, we denote to be in and set the replacement variable as Finally, the* s*-stage RK scheme becomeswith modified term and values at intermediate step The indication here is that one experiments entirely in the spectral domain and inverts transform to get . It should be noted that, by removing the associated stiffness, one can accurately implement the resulting system of ODEs with any higher-order time stepping schemes; see [16–19] for details.

Next, we formulate the fourth-order exponential time-differencing Runge-Kutta (ETDRK4) method by following closely most notations used by Cox and Matthews [16]. By multiplying through by the term , known as the integrating factor, we obtain (17) in the formOn integrating (27) over a single time step in the interval of length , that is, , we havewhere is the linear operator and is the Fourier transform of the nonlinear reaction functions. A direct application of the standard fourth-order Runge-Kutta method leads to a scheme by Cox and Matthew [16], known as the ETDRK4. Numerical experiments have shown that the cancellation errors in their scheme were too pronounced, which caused the method to suffer serious numerical instability and order reduction in its computation. We utilize in this work the modified Krogstad [18, 20] version of the ETDRK4 scheme that has been formulated to overcome the inherent numerical stability with smaller local truncation error,with the stages given as and with functions defined aswhich precisely coincide with the terms in the Lie group methods by Munthe-Kaas [21].

#### 4. Numerical Experiments

In this section, we intend to present the numerical results in one and two dimensions. We expect our results to reflect the mathematical results. We simulate the KiSS model in one and two dimensions with the ETDRK4 method presented in Section 3. Its stability analysis and comparison with other existing higher-order schemes when applied to some reaction-diffusion problems can be found in [18–20, 22] and references therein.

##### 4.1. One-Dimensional Example

In one dimension, consider the nonlinear KiSS modelsubject to the boundary conditions and initial population which has exponential decay as . Mathematically, one expects traveling waves to arise from such initial condition on an infinite domain which we have truncated in the interest of the work here at some large but finite value, . For the extended nonlinear example, we choose . In Figure 1, we illustrate the permanence effect and stability of the existence state via numerical simulations of (32). Ecological parameters , , and are chosen to ensure the long-term survival of the species population. It is noticeable that the ETDRK4 scheme converged even at the expense of variability of the growth rates , , , and for panels (a–d), respectively. It should be noted that the numerical simulation depicts the convergence of the solution to a positive survival state with population density of 1, for all and .