Bounding biomass in the Fisher equation

The Fisher-Kolmogorov-Petrovskii-Piskunov equation with a variable growth rate and advection by an incompressible velocity field is considered as a model for plankton dispersed by ocean currents. If the average growth rate is negative then the model has a survival-extinction transition; the location of this transition in the parameter space is constrained using variational arguments and delimited by simulations. The statistical steady state reached when the system is in the survival region of parameter space is characterized by integral constraints and upper and lower bounds on the biomass and productivity that follow from variational arguments and direct inequalities. In the limit of zero-decorrelation time the velocity field is shown to act as Fickian diffusion with an eddy diffusivity much larger than the molecular diffusivity: this allows a one-dimensional model to predict the biomass, productivity, and extinction transitions. All results are illustrated with a simple growth and stirring model.

Figure 1

Center: A snapshot of $\ln (\eta P∕{\gamma }_{\max })$ from a simulation of (1). The plankton concentration along the line $x=0$ is shown on the left-hand side and the growth rate $\gamma \left(y\right)∕{\gamma }_{\max }$ is on the right-hand side. This simulation uses the model in Sec. 2 with $\Gamma =0.1$, $U_{\ast }=10$, $m=1$, ${\tau }_{\ast }=0.25$, and ${\kappa }_{\ast }=1\times {10}^{-4}$.

Figure 2

The model growth rate (2) for various $\Gamma$.

Figure 3

(a) The solid curve shows a Monte-Carlo estimate of the non-dimensional Lyapunov function $\Lambda \left({\tau }_u\right)$ defined in (8); the dashed curve is the fit in (9). (b) Spectra of the inert scalar $C$ in (12) and of the biological tracer $P$ obtained from (1). The spectra are “compensated” by multiplying the variance spectrum by the wave number so that a $k^{-1}$ Batchelor spectrum appears flat. The plankton are growing according to the sinusoidal growth rate in (2) with $\Gamma =0$. The $P$-spectra are line spectra at $y∕\ell =0$, where $\gamma ={\gamma }_{\max }$, and at $y∕\ell =\pi ∕2$, where $\gamma =0$. $C$ is well described by the Batchelor spectrum, while $P$ has a steeper spectral slope. Panels (c) and (d) show snapshots of the $P$ and $C$, respectively; note in panel (c) the contour interval is logarithmic [i.e., panel (c) shows $Z\equiv \ln (\eta P∕{\gamma }_{\max })$]. Both simulations use the parameters ${\kappa }_{\ast }={10}^{-7}$, $U_{\ast }=1$, ${\tau }_{\ast }=2.2214$ and $m=1$.

Figure 4

(a) The instantaneous biomass $\overline{P}\left(t\right)$ shows the exponential growth of the strange eigenmode in (14) followed by the dynamic equilibrium when the nonlinear saturation term becomes important. The initial condition is $P(\mathit{x},0)={10}^{-7}$ and the nondimensional parameters are noted at the bottom of the figure. The dashed line is the theoretical maximum $\overline{P}\left(t\right)=\overline{P}\left(0\right)e^{{\gamma }_{\max }t}$ [11]. (b) Increasing $U_{\ast }$ in (3) decreases $s_{\ast }$ and results in $s_{\ast }<0$ (extinction) if $U_{\ast }$ is greater than about 6.

Figure 5

(a) The first three eigenvalues of (23). Extinction occurs when $s_1<0$, which for ${\kappa }_{\ast }=0.1$ means $\Gamma \lesssim -20$. (b) The solid curve is the $U_{\ast }=0$ extinction transition computed by numerically solving the linearized stability problem (23) with $s_{\ast }=0$; the dashed curve is the approximation in (24).

Figure 6

(a) The eigenvalue $s_{\ast }$ in (23) as a function of ${\kappa }_{\ast }$. The approximation (25) improves as ${\kappa }_{\ast }\to 0$. (b) The eigenmode $\widehat{P}\left(y\right)$ in (23) and the Gaussian approximation from Appendix . The WKB approximation (27) has a singularity near the origin and is not shown. (c) Three solutions for ${\widehat{P}}_y∕\widehat{P}$. Both the Gaussian and WKB approximations are good near the origin, but the Gaussian approximation fails in the region where $\gamma \left(y\right)-s<0$ and is therefore not useful for evaluating the average in (21).

Figure 7

(a) The biomass $⟨P⟩$ as a function of $U_{\ast }$ as determined by Monte Carlo simulation. The extinction transition is at around $U_{\ast }=7.2$. (b) An expanded view of the data in panel (a). (c) Summary of a suite of simulations in which $D_{\ast }=U_{\ast }^2{\tau }_{\ast }∕8$ is fixed as ${\tau }_{\ast }$ is varied. The survival-extinction transition is located via repeated simulations with varying $U_{\ast }$, as in panels (a) and (b). As ${\tau }_{\ast }$ is decreased the survival-extinction transition approaches the solid curve previously shown in Fig. 5b with $D_{\ast }$ playing the role of ${\kappa }_{\ast }$. The velocity field is (3) with $m=8$ to ensure scale separation.

Figure 8

(a) A schematic illustration; below the dashed curve the sufficient condition (38) guarantees survival. However, the actual survival region may be much larger, as indicated by the region below the solid curve. (b) Below each solid curve the sufficient condition (38), obtained by numerical solution of (37), guarantees survival for the indicated value of $\Gamma$ in the sinusoidal growth model (2).

Figure 9

The dashed curves are the basic lower (43) and upper (56) bounds: $\max \{0,⟨\gamma ⟩\}⩽\eta \mathcal{B}⩽(\sqrt{⟨{\gamma }^2⟩}+⟨\gamma ⟩)∕2$. The solid curve is the optimized upper bound (55) and the dashed-dotted curve is the lower bound (48). The circles are the biomass obtained from simulations. The lower bound (48) is tighter than the basic bound (43) only for $\Gamma \lesssim 0.075$.

Figure 10

(a) A time series of $\overline{P}\left(t\right)$ compared with two upper bounds; the dashed line is the simple upper bound in (58) and the straight solid line is the optimized upper bound in (59) using the trial function $f=\exp \left[0.3960\cos \left(k_{1}y\right)\right]$. The dotted line shows the average of the biomass from $t_{\ast }=100$ to 400. (b) and (c) The open circles indicate the biomass obtained from numerical simulation for the indicated parameter values. The dashed lines are the upper and lower bounds in (58) and the solid lines are the bounds in (59) with the upper bound $p$ optimized as in (57). The lower bound in (59) only appears in (c) where it is not very tight but does forbid extinction for small $U_{\ast }$.

Figure 11

(a) Bounding ellipses obtained from (64) and the sinusoidal model in (2); units are nondimensional with ${\gamma }_{\max }=\eta =1$. As $\Gamma$ decreases the growth rate becomes more spatially inhomogeneous and the ellipse expands and the simultaneous bounds become less tight. However, because $\mathcal{P}⩽{\gamma }_{\max }\mathcal{B}$ much of the expanded ellipse is inaccessible. (b) The shaded region indicates the allowed region if $⟨\gamma ⟩>0$. The point $(\mathcal{B},\mathcal{P})={\eta }^{-1}(⟨\gamma ⟩,{⟨\gamma ⟩}^2)$ is indicated by $A$. (c) The shaded region indicates the allowed region if $⟨\gamma ⟩<0$. In this case the point $A$ lies outside the first quadrant and so the constraint (68) is toothless.

Figure 12

(a) The joint bounds illustrated with variable $U_{\ast }$ and fixed ${\tau }_u$. $U_{\ast }$ increases from 0.1 to 50. (b) The joint bound illustrated with variable ${\kappa }_{\ast }$; the maximum value of ${\kappa }_{\ast }$ is 1 for both $U_{\ast }=0$ and 2.

Figure 13

(a) Transects along the line $x=0$ for two solutions of (1) with ${\kappa }_{\ast }=0$ using the same sequence of phases ${\varphi }_x$ and ${\varphi }_y$ in (3) but different values of $\lambda$ in (7). $\lambda$ is varied by holding ${\tau }_u$ and $m$ fixed and varying $U_{\ast }$. The transect with $\lambda =0.16$ is smooth, but the transect with $\lambda =1.03$ is not. (b) The Hölder exponent (69). The vertical dashed line is the line $\lambda =\Gamma$ and coincides with the inflection point of $\alpha \left(\lambda \right)$. (c) Transects along the line $x=0$ for two solutions of (1) with ${\kappa }_{\ast }={10}^{-4}$. The $\lambda =0.16$ transect is very similar to the corresponding transect in (a), but the $\lambda =1.03$ transect is noticeably smoother. (d) The productivity $\mathcal{P}$ as a function of $\lambda$.