% Build: latexmk -pdf main.tex. Bibliography included; no BibTeX step needed. \documentclass[11pt,a4paper]{article} \usepackage[T1]{fontenc} \usepackage{lmodern} \usepackage[margin=27mm]{geometry} \usepackage{amsmath,amssymb,amsthm,mathtools} \usepackage{microtype} \usepackage{graphicx,booktabs,array} \usepackage{enumitem} \usepackage{aliascnt} \usepackage[numbers,sort&compress]{natbib} \renewcommand{\bibfont}{\small} \setlength{\bibsep}{4pt} \usepackage{xcolor} \usepackage[colorlinks=true,linkcolor=blue!45!black,citecolor=blue!45!black,urlcolor=blue!45!black]{hyperref} \usepackage[nameinlink,noabbrev]{cleveref} \hypersetup{pdftitle={Limit Cycles under Concave Per-Capita Growth}, pdfauthor={Zhenyu Sun}, pdfsubject={Exact global cycle counts, polynomial degree, finite jets, and strictly concave predator-prey growth}} \pagestyle{plain} \setlength{\emergencystretch}{2em} \setlist{itemsep=3pt,topsep=5pt} \linespread{1.035} \newtheorem{theorem}{Theorem}[section] \newaliascnt{proposition}{theorem} \newtheorem{proposition}[proposition]{Proposition} \aliascntresetthe{proposition} \crefname{proposition}{proposition}{propositions} \Crefname{proposition}{Proposition}{Propositions} \newaliascnt{lemma}{theorem} \newtheorem{lemma}[lemma]{Lemma} \aliascntresetthe{lemma} \crefname{lemma}{lemma}{lemmas} \Crefname{lemma}{Lemma}{Lemmas} \newaliascnt{corollary}{theorem} \newtheorem{corollary}[corollary]{Corollary} \aliascntresetthe{corollary} \crefname{corollary}{corollary}{corollaries} \Crefname{corollary}{Corollary}{Corollaries} \theoremstyle{definition} \newaliascnt{definition}{theorem} \newtheorem{definition}[definition]{Definition} \aliascntresetthe{definition} \crefname{definition}{definition}{definitions} \Crefname{definition}{Definition}{Definitions} \theoremstyle{remark} \newaliascnt{remark}{theorem} \newtheorem{remark}[remark]{Remark} \aliascntresetthe{remark} \crefname{remark}{remark}{remarks} \Crefname{remark}{Remark}{Remarks} \newcommand{\R}{\mathbb R} \newcommand{\C}{\mathbb C} \newcommand{\N}{\mathbb N} \newcommand{\dd}{\,d} \newcommand{\supp}{\operatorname{supp}} \newcommand{\sgn}{\operatorname{sgn}} \newcommand{\odd}{\mathrm{odd}} \newcommand{\even}{\mathrm{even}} \newcommand{\norm}[1]{\left\|#1\right\|} \title{\bfseries Limit Cycles under Concave Per-Capita Growth} \author{Zhenyu Sun} \date{} \begin{document} \maketitle \vspace{-0.5em} \begin{abstract} We study predator--prey systems with a fixed globally linear functional response and a strictly concave prey per-capita growth rate. For every $N\ge1$ and $K\ge2$, we construct a growth polynomial of degree $2N+2\lceil K/2\rceil+2$ with exactly $N$ positive periodic orbits, all hyperbolic. The construction retains a chosen center $K$-jet and the low-density growth rate and carrying capacity of its limiting profile, subject to an explicit relation between the center height and curvature. A limiting Lyapunov function and return-map estimates at three amplitude scales exclude additional cycles. A dissipative extension gives exact global counts with prescribed limiting center energies. We also prove that profiles with at least $N$ hyperbolic cycles are $C^2$-dense in the admissible class, while preserving the entire germ at the maximum and both ecological endpoint values. The first nonzero canonical odd term determines the sharp differentiability threshold. Further results give parity restrictions for fixed focus jets, degree lower bounds, quantitative polynomial approximation under concavity and interpolation constraints, and amplitude and multiplicity interpolation by Abelian integrals. Finally, a single smooth strictly concave profile has infinitely many nested hyperbolic cycles converging to coexistence. \end{abstract} \smallskip \noindent\textbf{Keywords.} Predator--prey system; exact cycle count; polynomial degree; concavity; finite jet; Abelian integral; differentiability. \noindent\textbf{2020 Mathematics Subject Classification.} 34C07, 34C23, 34C25, 92D25. \section{Introduction and main results} \label{sec:introduction} Consider \begin{equation}\label{eq:original} \dot x=x\bigl(\alpha(x)-\beta y\bigr),\qquad \dot y=y(k\beta x-m),\qquad k,m,\beta>0. \end{equation} Throughout the paper the three constants are fixed. The function $\alpha$ is the prey \emph{per-capita} growth rate; the total prey growth in the absence of predation is $x\alpha(x)$. This distinction matters for shape assumptions. Put \[ a=\frac{m}{k\beta}. \] If $\alpha(a)=c>0$, the unique positive equilibrium is $E_c=(a,c/\beta)$. Svirezhev's ICM report \cite[p.~1682]{Svirezhev1984} asks whether \eqref{eq:original}, with a single-peaked growth rate, can possess a third, a fourth, or further limit cycles. In one regime, a subcritical Hopf bifurcation combined with the global phase portrait yields a two-cycle configuration. The question concerns a globally linear response $\beta x$. In some later ecological literature the name Holling I also denotes a linear response followed by saturation; we use the explicit expression $\beta x$ throughout. Multiple cycles have been constructed under several ecological shape constraints. Kuang \cite{Kuang1988} constructs a Gause-type system with a strictly concave prey isocline and at least three cycles. Wrzosek \cite{Wrzosek1990} obtains arbitrarily many cycles with a unimodal total prey-growth function and a concave functional response. Zegeling and Kooij \cite{ZegelingKooij2021} compare center, Hopf, and singular perturbation mechanisms in broader response families. These results motivate studying the restrictions imposed by fixing the response while retaining strong shape conditions on the actual per-capita growth. Our methods use established planar theory. The energy-coordinate reduction is the generalized Filippov transformation studied by Gasull and Guillamon \cite{GasullGuillamon1995}; their paired-branch criteria include the center and nonoscillation principles used below. The construction of multiple cycles through Abelian integrals has a direct predecessor in Gasull, Guillamon, Li, and Zhang \cite{GasullGuillamonLiZhang1996}. The global zero bound in \cref{sec:moments} is an explicit application of Gavrilov and Iliev \cite{GavrilovIliev2003}. The results developed here concern simultaneous preservation of the fixed predator equation, concavity, local data, ecological endpoint values, and prescribed oscillations. Our explicit degree theorem uses only harmonic moments and a global Lyapunov argument. The Gavrilov--Iliev theorem enters the separate problem of interpolation throughout the regular annulus of an explicit center. \begin{definition}\label{def:admissible} An admissible single-peaked profile at $a$ is a function $\alpha\in C^\infty([0,\infty))$ satisfying \[ \alpha(0)>0,\quad \alpha'(x)>0\ (0a), \] and having a zero $x_1>a$. It is \emph{strictly concave} if $\alpha''(x)<0$ for all $x\ge0$. \end{definition} For a strictly concave profile, the right-hand zero exists automatically: $\alpha'$ is negative and decreasing to the right of $a$, so $\alpha(x)\to-\infty$. Its positive zero is unique. Our first result combines an exact global count with an explicit degree. For $\theta>0$ define the analytic germ at $a$ \[ C_\theta(x)=c-\theta m\left(\frac xa-1-\log\frac xa\right). \] This is a center profile in the canonical coordinates used below. For functions defined near $a$, $j_a^K$ denotes equality of derivatives through order $K$ at $a$. Polynomial degrees in this paper refer to the growth function $\alpha$; a growth polynomial of degree $D$ gives a polynomial vector field of degree $D+1$ in \eqref{eq:original}. \begin{theorem}[Exact global count with a linear degree bound] \label{thm:exact-degree} Fix $k,m,\beta,c,\theta>0$, integers $N\ge1$, $K\ge2$, and distinct numbers $0\theta m(\mathsf H_D-1)$, where $\mathsf H_D=\sum_{n=1}^D1/n$. There exist $\alpha_0>0$, $b>a$, $\rho_0>0$, and an explicitly defined family of real polynomials $\widehat P_\rho$, $0<\rho<\rho_0$, such that \begin{gather*} \deg\widehat P_\rho=D,\qquad \widehat P_\rho''(x)<0\quad(x\ge0),\qquad j_a^K\widehat P_\rho=j_a^KC_\theta,\\ \widehat P_\rho(0)=\alpha_0,\qquad \widehat P_\rho(b)=0. \end{gather*} The system \eqref{eq:original} with growth $\widehat P_\rho$ has exactly $N$ positive periodic orbits, all hyperbolic. On the section $v=0$, $s>0$ of \eqref{eq:energy-coordinate}, their coordinates satisfy \[ s_i(\rho)=\rho\sqrt{2h_i}+O(\rho^2). \] Their stabilities alternate, the outermost is attracting, and the coexistence equilibrium is a stable weak focus for even $N$ and an unstable weak focus for odd $N$. \end{theorem} The numbers $\alpha_0,b$ are determined by the chosen limiting polynomial and then remain fixed throughout the family. For every fixed $c>0$ and each $(N,K)$, one can satisfy the displayed height condition by choosing $\theta>0$ sufficiently small. For a preassigned center curvature, the same degree conclusion is subject to that condition. The proof, including the exact endpoint correction and exclusion of every additional periodic orbit, is given in \cref{sec:explicit-polynomials}. Two restrictions accompany this degree bound. First, when a prescribed jet fixes the stability of a focus, the total number of hyperbolic positive periodic orbits must be even for a stable focus and odd for an unstable focus (\cref{prop:parity}). The center jet in \cref{thm:exact-degree} leaves that stability free at a higher order. Second, degree bounds for arbitrary prescribed jets or for approximation of arbitrary baselines depend on the sizes of those data. \Cref{sec:degree} gives explicit lower bounds and a constructive data-dependent upper bound. This is a sufficient degree bound for the chosen center jet. The proof separates the local creation of cycles from their global exclusion. In energy coordinates, the even part of the growth profile defines a reversible center, and the odd part determines the leading return displacement. Harmonic moments prescribe its simple zeros at a shrinking scale. A Lyapunov function for the limiting polynomial forces every possible cycle toward the equilibrium; return-map estimates then exclude cycles at both smaller and larger relative amplitudes. The second result treats arbitrary local growth profiles. \begin{theorem}[Localized realization with preservation of the germ] \label{thm:main-density} Let $\alpha$ be an admissible single-peaked profile at $a$. For every $N\ge1$ there are admissible profiles $\alpha_\rho$, $\rho\downarrow0$, such that: \begin{enumerate}[label=(\roman*)] \item $\alpha_\rho=\alpha$ on a neighborhood of $a$ and outside an interval of length $O(\rho)$ centered at $a$; \item $\norm{\alpha_\rho-\alpha}_{C^2([0,\infty))}\to0$; \item \eqref{eq:original} with $\alpha_\rho$ has at least $N$ distinct hyperbolic positive limit cycles converging to $E_c$. \end{enumerate} The neighborhood on which the germs agree may depend on $\rho$. If $\alpha$ is globally strictly concave, every $\alpha_\rho$ can be chosen globally strictly concave. The values $\alpha(0)$, $\alpha(a)$, the zero $x_1$, and the entire growth profile away from the shrinking interval are unchanged. \end{theorem} The $C^2$ conclusion can be refined. In the energy coordinate introduced in \cref{sec:normal-form}, let the first nonzero odd term of the growth profile have degree $d$. Then $d\ge3$ is odd, and the construction converges in $C^{d-1}$. A sufficiently small $C^d$ perturbation that fixes the lower jet has no cycle in a fixed sufficiently small neighborhood of the equilibrium. This sharp threshold is proved in \cref{sec:localized,sec:threshold}. Polynomial profiles can retain additional ecological data. \begin{theorem}[Strictly concave polynomial realization] \label{thm:main-polynomial} Let $\alpha$ be an admissible globally strictly concave profile at $a$, with positive zero $x_1$. Fix integers $N\ge1$ and $K\ge2$, a number $R>x_1$, and $\varepsilon>0$. There is a real polynomial $\widehat\alpha$ such that \begin{gather*} \widehat\alpha''(x)<0\quad(x\ge0),\qquad j_a^K\widehat\alpha=j_a^K\alpha,\qquad \widehat\alpha(0)=\alpha(0),\qquad \widehat\alpha(x_1)=0,\\ \norm{\widehat\alpha-\alpha}_{C^2([0,R])}<\varepsilon, \end{gather*} and \eqref{eq:original} has at least $N$ distinct hyperbolic positive limit cycles. Its unique maximum is $a$, and its unique positive zero is $x_1$. \end{theorem} The degree in this approximation theorem can increase with the jet, the approximation tolerance, and the chosen hyperbolicity margins. \Cref{sec:degree} quantifies this dependence. Its conclusion is an at-least-$N$ statement for an arbitrary baseline. Exact global counts at specified regular center energies are obtained by the dissipative extension theorem in \cref{sec:exact-extension}, which also admits polynomial approximation with its count preserved. The local realization theorem has a smooth infinite counterpart. \begin{theorem}[One profile with infinitely many cycles] \label{thm:main-infinite} For every fixed $k,m,\beta>0$ and $c>0$, there is an admissible globally strictly concave $C^\infty$ profile $\alpha$, with maximum at $a$ and $\alpha(a)=c$, for which \eqref{eq:original} has infinitely many distinct nested hyperbolic limit cycles converging to $E_c$. All derivatives of $\alpha$ at $a$ agree with those of \[ \alpha_{\mathrm{cen}}(x) =c-m\left(\frac xa-1-\log\frac xa\right). \] \end{theorem} The explicit center yields an additional finite-dimensional theory. For odd perturbations of degrees $2r+3,\ldots,2r+2N+3$, the cyclicity of every compact regular subannulus is exactly $N$. The limiting energies, and even the actual section coordinates and multiplicities, can be prescribed. Finally, \cref{sec:responses} treats a fixed nonlinear response $p$. Under the simple break-even hypothesis, the condition $ap'(a)<2p(a)$ permits strictly concave per-capita realizations; a response with elasticity at least two everywhere excludes cycles for every strictly concave positive low-density growth profile. \paragraph{Organization.} \Cref{sec:normal-form,sec:global-dynamics} give the coordinate reduction, uniform confinement, and the global parity law. \Cref{sec:explicit-polynomials} proves \cref{thm:exact-degree}. \Cref{sec:moments,sec:exact-extension} treat prescribed regular amplitudes and exact global extension. \Cref{sec:localized,sec:threshold,sec:global,sec:degree} develop local density, interpolation, and quantitative degree estimates. \Cref{sec:infinite,sec:responses,sec:computations} give smooth infinite realization, fixed nonlinear responses, and reproducible computations. % --- normal_form.tex --- \section{Energy coordinates and an explicit center} \label{sec:normal-form} \subsection{Basic geometry and periodic averages} The positive quadrant is invariant. For an admissible profile, every positive solution is bounded forward in time: $x(t)$ is bounded above by $\max\{x(0),x_1\}$, and $z=kx+y$ satisfies \[ \dot z=kx\alpha(x)-my =kx(\alpha(x)+m)-mz\le C-mz \] for a constant $C$ on this bounded prey interval. The Jacobian at $E_c=(a,c/\beta)$ has trace $a\alpha'(a)$ and determinant $mc$. In particular, at a growth maximum located at $a$, its eigenvalues are $\pm i\sqrt{mc}$. \begin{proposition}[Identities for periodic coexistence] \label{prop:periodic-identities} For any positive nonconstant periodic orbit of period $T$, write $\langle u\rangle=T^{-1}\int_0^T u(t)\,dt$. Then \[ \langle x\rangle=a,\qquad \beta\langle y\rangle=\langle\alpha(x)\rangle,\qquad \log\mu=\int_0^T x(t)\alpha'(x(t))\,dt, \] where $\mu$ is its transverse Floquet multiplier. If $\alpha$ is strictly concave, then \[ 0<\beta\langle y\rangle<\alpha(a). \] More precisely, if $0<\kappa\le-\alpha''\le K$ on the prey projection of the orbit, then \[ \frac\kappa2\langle(x-a)^2\rangle \le \alpha(a)-\beta\langle y\rangle \le\frac K2\langle(x-a)^2\rangle. \] \end{proposition} \begin{proof} Integrate the equations for $\log y$ and $\log x$ over one period. In the divergence integral, these two identities cancel all terms except $x\alpha'(x)$, proving the multiplier formula. Strict Jensen inequality gives the strict average-growth bound: $x$ cannot be constant on a nonconstant orbit, since the prey equation would then also force $y$ to be constant. The two quantitative inequalities follow by Taylor's theorem at $a$ and $\langle x-a\rangle=0$. The prey projection contains $a$, by the first identity, so both curvature bounds apply throughout each Taylor segment. \end{proof} These identities constrain the means and stabilities of all the cycles constructed below, even when their number is arbitrarily large. \subsection{Canonical reduction} The following reduction is a specialization of the generalized Filippov transformation \cite[\S3]{GasullGuillamon1995}. We include the calculation to fix the time and energy conventions. \begin{lemma}[Exact energy-coordinate reduction]\label{lem:reduction} Suppose $\alpha(a)=c>0$. Define \begin{equation}\label{eq:energy-coordinate} s(x)=\sgn(x-a) \sqrt{2m\left(\frac xa-1-\log\frac xa\right)}, \qquad v=\log\frac{\beta y}{c}. \end{equation} Then $s$ is a smooth increasing diffeomorphism from $(0,\infty)$ onto $\R$, with $s(a)=0$ and $s'(a)=\sqrt m/a$. The positive time change \[ \frac{d\tau}{dt}=xs'(x) =\frac{m(x/a-1)}{s(x)} \] transforms \eqref{eq:original} into \begin{equation}\label{eq:canonical} s'=f(s)-c(e^v-1),\qquad v'=s, \qquad f(s)=\alpha(x(s))-c. \end{equation} The time factor extends at $x=a$ to $\sqrt m$. \end{lemma} \begin{proof} The expression inside the square root has a double zero at $a$ with positive quadratic coefficient. The signed square root is therefore smooth there. Differentiation gives $sxs'(x)=m(x/a-1)$. Both factors $s$ and $x/a-1$ have the same sign. Finally, \[ \dot s=xs'(x)(\alpha(x)-ce^v),\qquad \dot v=m(x/a-1)=xs'(x)s, \] which proves the transformed equations. Since $s(x)\to-\infty$ as $x\downarrow0$ and $s(x)\to+\infty$ as $x\to\infty$, the transformation is global on the positive quadrant. \end{proof} At a nondegenerate maximum of $\alpha$, \[ f(0)=f'(0)=0,\qquad f''(0)=\frac{a^2}{m}\alpha''(a)<0. \] Write $f=E+O$ for its even and odd parts. \begin{lemma}[Even-part first integral]\label{lem:even-integral} For every smooth even $E$ with $E(0)=0$, the field $s'=E(s)-c(e^v-1)$, $v'=s$ has a nondegenerate center at the origin. There is a smooth local function $H_E$ with a nondegenerate minimum such that, for the full field \eqref{eq:canonical}, \begin{equation}\label{eq:odd-Lyapunov} \frac{dH_E}{d\tau}=K_w(s^2/2,v)\,sO(s),\qquad K_w>0. \end{equation} \end{lemma} \begin{proof} Put $w=s^2/2$ and write $E(s)=\mathcal E(w)$ for $w\ge0$. The smooth function $\mathcal E$ has a smooth extension across zero. Consider the scalar nonautonomous equation \[ \frac{dw}{dv}=\mathcal E(w)-c(e^v-1). \] Let $K(w,v)$ be the value at $v=0$ of the solution through $(w,v)$. Smooth dependence gives $K_w>0$, $K(w,0)=w$, and \[ K_v+[\mathcal E(w)-c(e^v-1)]K_w=0. \] Set $H_E(s,v)=K(s^2/2,v)$. Its quadratic part is $(s^2+cv^2)/2$, so its sufficiently small positive levels are ovals. For the even field it is a first integral; the vector field has no other zero near the origin, hence these ovals are periodic. For the full field, substitution into the displayed equation gives \eqref{eq:odd-Lyapunov}. \end{proof} This is the local form of the paired-energy-branch center and nonoscillation criteria in \cite[Theorems A and B]{GasullGuillamon1995}. In particular, if $O(s)=bs^d+o(s^d)$, where $d$ is odd and $b\ne0$, then $sO(s)$ has one strict sign away from $s=0$ in a sufficiently small neighborhood. There are no periodic orbits in that neighborhood. For analytic $f$, the equilibrium is a center if and only if $f$ is even. For the explicit profile \begin{equation}\label{eq:center-profile} \alpha_{\mathrm{cen}}(x) =c-\frac12s(x)^2 =c-m\left(\frac xa-1-\log\frac xa\right), \end{equation} one has \[ \alpha_{\mathrm{cen}}'(a)=0,\qquad \alpha_{\mathrm{cen}}''(x)=-\frac{m}{x^2}<0. \] Its canonical system has the first integral \begin{equation}\label{eq:center-integral} H_c(s,v)=\frac12e^v s^2+\frac c2(e^v-1)^2. \end{equation} Indeed, differentiating \eqref{eq:center-integral} along $s'=-c(e^v-1)-s^2/2$, $v'=s$, gives zero. The profile \eqref{eq:center-profile} is used locally; global smooth strictly concave profiles agreeing with it near $a$ are constructed in \cref{sec:global}. For the explicit moment calculations we set $c=1$. The general case reduces to this normalization by $s=\sqrt c\,u$ and the time change $d\sigma=\sqrt c\,d\tau$. All amplitude statements in \cref{sec:moments} use the normalized energy. % --- global_dynamics.tex --- \section{Global confinement, return maps, and parity} \label{sec:global-dynamics} Throughout this section, $k,m,\beta>0$ are fixed, $a=m/(k\beta)$, and $\alpha$ is an admissible profile with its maximum at $a$. Write $c=\alpha(a)$, $\alpha_0=\alpha(0)>0$, and let $x_1>a$ be its positive zero. The results below use the single-peak assumptions; global strict concavity is not required. \subsection{A common compact set for all periodic orbits} \begin{lemma}[Uniform periodic-orbit bounds] \label{lem:cycle-compactness} Consider a family of admissible profiles satisfying \[ \alpha(0)\ge d_0>0,\qquad \alpha(a)\le c_1,\qquad x_1\le B. \] Every positive nonconstant periodic orbit of every member of this family is contained in \begin{equation}\label{eq:cycle-rectangle} Q=[\ell,B]\times[y_*,M],\qquad M=\frac{kB(c_1+m)}m, \quad \ell=\frac a2\exp\!\left(-\frac{2\beta M}m\right), \quad y_* =\frac{d_0}{2\beta} \left(\frac{\ell}{a}\right)^{2m/d_0}. \end{equation} In particular, no sequence of such periodic orbits can approach either coordinate axis or escape to infinity. \end{lemma} \begin{proof} At a largest prey coordinate on a periodic orbit, $\alpha(x)=\beta y>0$, whence $xa$. This proves the lower predator bound. \end{proof} The lemma needs only uniform bounds on three scalar quantities. In particular, it applies both to families retaining the ecological endpoint data exactly and to sufficiently small perturbations of those data. For a fixed profile, one may take $d_0=\alpha_0$, $c_1=c$, and $B=x_1$. \begin{lemma}[Confinement by a limiting Lyapunov function] \label{lem:lyapunov-confinement} Let $X_n\to X_*$ in $C^1$ on a neighborhood of a compact set $Q$, and let $\Gamma_n\subset Q$ be periodic orbits of $X_n$. Suppose that a smooth function $H$ satisfies $X_*H\le0$ there. Let $\mathcal I$ denote the largest compact invariant subset of $\{X_*H=0\}\cap Q$. If $H\le h_*$ on $\mathcal I$, every Hausdorff limit of the $\Gamma_n$ is contained in $\{H\le h_*\}\cap Q$. In particular, if $\mathcal I=\{E\}$ and $H$ has its unique minimum on $Q$ at $E$, then $\Gamma_n\to\{E\}$ in Hausdorff distance. \end{lemma} \begin{proof} Take a subsequence converging in Hausdorff distance to a nonempty compact set $L\subset Q$. Continuous dependence of solutions, or equivalently compactness of the solutions through points of $\Gamma_n$, shows that $L$ is invariant under $X_*$ in both time directions. Choose $z\in L$ maximizing $H$. For $t\le0$, monotonicity gives $H(\phi_t(z))\ge H(z)$, while maximality gives the reverse inequality. The complete backward orbit through $z$ therefore lies in $\{X_*H=0\}\cap Q$, and $H$ is constant on it. Its alpha-limit is a nonempty compact invariant subset of this zero-dissipation set, hence is contained in $\mathcal I$. It follows that $H(z)\le h_*$ and consequently $H\le h_*$ on $L$. In the stated special case, $H(z)=H(E)$, so uniqueness of the minimum forces $L=\{E\}$. The same argument applies to every convergent subsequence. \end{proof} We will also use the following neighborhood formulation. Suppose $X_*$ has a closed center disc $\mathcal D$ such that $\mathcal I=\mathcal D$ and $H\le h_*$ on $\mathcal D$. Within an admissible class having a common compact set \eqref{eq:cycle-rectangle}, every open neighborhood of $\{H\le h_*\}\cap Q$ contains all periodic orbits of every sufficiently $C^1$-close field. Otherwise a sequence of excluded periodic orbits would contradict \cref{lem:lyapunov-confinement}. This is the form needed when a local center is extended by a dissipative exterior. \subsection{Persistence and a global transverse section} \begin{lemma}[Interior persistence and a global return map] \label{lem:global-return} Every positive solution is eventually contained in a compact subset of the positive quadrant. The half-line \[ \Sigma=\{(a,y):y>c/\beta\} \] has a smooth first-return map $P:(c/\beta,\infty)\to(c/\beta,\infty)$ with $P'>0$. Every positive periodic orbit meets $\Sigma$ exactly once, and $P(y)a$ so close to $a$ that $\alpha(b)>3\alpha_0/4$. While $x\le b$ and $y0. \] This proves interior persistence without an assumption on the stability of the coexistence equilibrium. Start now at a point of $\Sigma$. Since $\dot x<0$ there, the trajectory enters $xa$. The same argument, now with $y$ increasing, gives a further finite crossing, necessarily on $\Sigma$. All these crossings are transverse. Smooth dependence of the flow and the implicit function theorem give a smooth first-return map. Its derivative is positive; explicitly, the standard variational formula expresses it as a positive exponential multiplied by the ratio of the two negative transverse velocities. Analytic dependence gives the analytic assertion when $\alpha$ is analytic. Every positive periodic orbit surrounds $E_c$, by the planar index formula, since there is no other equilibrium in the positive quadrant. Along $\Sigma$ every crossing points to the left. The intersections of an oriented Jordan curve with a ray transverse to it alternate between entry and exit as one moves along the ray. The ray from $E_c$ to infinity begins inside the curve and ends outside it. Since all crossing orientations here agree, there is exactly one intersection. Thus periodic orbits correspond precisely to fixed points of $P$. Finally, a trajectory starting at $(a,y)$ satisfies $x(t)0$. With $M=kx_1(c+m)/m$, scalar comparison gives \[ z(t)\le M+(z(0)-M)e^{-mt}M. \] At the return to $\Sigma$ the prey coordinate is again $a$, so $P(y)0$, the cycles are nested, the outermost cycle is attracting, and their stabilities alternate. Moreover, \[ E_c\text{ stable}\ \Longrightarrow\ N\text{ even}, \qquad E_c\text{ unstable}\ \Longrightarrow\ N\text{ odd}. \] Here $N=0$ is allowed in the stable case. If the first nonzero canonical odd coefficient is $O(s)=bs^d+o(s^d)$, with $d$ odd, then $b<0$ gives the stable case and $b>0$ the unstable case. A prescribed jet determining this coefficient therefore determines the parity of any such global configuration. \end{proposition} \begin{proof} All cycles surround the unique equilibrium, so disjointness of trajectories makes them nested. By \cref{lem:global-return}, $P(y)-y$ has precisely $N$ simple zeros and is negative for large $y$. Its sign alternates at the zeros. The outermost zero is therefore attracting, and the stabilities alternate inward. Near $c/\beta$ the sign is negative for a stable focus and positive for an unstable focus, proving the parity assertion. The relation between $b$ and stability follows from \cref{lem:even-integral}: $\dot H_E$ has the sign of $b$ away from $s=0$, and the only invariant subset of the local zero-dissipation line is the equilibrium. \end{proof} The hypothesis concerns \emph{all positive periodic orbits}. For a smooth growth rate, finiteness of the isolated limit cycles alone does not exclude a separate period annulus. Such an annulus would contribute an interval of fixed points to $P$ and invalidate the simple-zero sign count. The constructions below exclude all additional periodic orbits, so this distinction causes no ambiguity. \begin{proposition}[No analytic center and finiteness of analytic cycles] \label{prop:analytic-finite} If an admissible profile is real analytic on $(0,\infty)$, its coexistence equilibrium is a weak focus of finite order. It has only finitely many positive periodic orbits, each isolated. In particular, a globally strictly concave admissible polynomial cannot possess a period annulus or infinitely many positive limit cycles. If all its cycles are hyperbolic, their number satisfies \cref{prop:parity}. \end{proposition} \begin{proof} The signed energy coordinate is a real analytic diffeomorphism from $(0,\infty)$ to $\mathbb R$, including across $x=a$. Hence $f(s)=\alpha(x(s))-c$ is real analytic on $\mathbb R$. If the equilibrium were a center, the analytic center criterion following \cref{lem:even-integral} would make $f$ even near zero. The identity theorem would then make it even on all of $\mathbb R$. As $s\to-\infty$, one has $x(s)\to0$, and therefore $f(s)\to\alpha_0-c$. Evenness would imply $\alpha(x)\to\alpha_0>0$ as $x\to\infty$, contradicting strict decrease beyond the positive zero $x_1$. Thus the analytic odd part has a first nonzero term of some finite odd degree $d\ge3$. By \cref{lem:even-integral}, the small oval levels of the even-part first integral are one-sided barriers. A periodic orbit can neither remain inside such an oval nor enter and leave it. Consequently a neighborhood of the equilibrium is met by no positive periodic orbit, and the equilibrium is a stable or unstable weak focus. The global return map in \cref{lem:global-return} is analytic and is not the identity, since $P(y)x_{1,*}$, there exists $\epsilon>0$ such that every admissible profile $\widetilde\alpha$ satisfying \[ j_a^{d-1}\widetilde\alpha=j_a^{d-1}\alpha_*, \qquad \|\widetilde\alpha-\alpha_*\|_{C^d([0,R])}<\epsilon \] also has exactly $N$ positive periodic orbits, all hyperbolic. The endpoint values need not be retained exactly. \end{proposition} \begin{proof} The fixed jet retains $c=\alpha_*(a)$ and hence the same coexistence equilibrium and section $\Sigma$. For sufficiently small $\epsilon$, one has $\widetilde\alpha(0)\ge\alpha_*(0)/2$ and $\widetilde\alpha(R)<0$. Admissibility then gives $\widetilde x_10$, define \begin{equation}\label{eq:theta-center} C_\theta(x)=c-\theta m\bigl(u-\log(1+u)\bigr) =c-\frac\theta2s(x)^2,\qquad x>0. \end{equation} This local growth law is strictly concave, with $C_\theta''(x)=-\theta m/x^2$. The associated canonical field has the first integral \begin{equation}\label{eq:theta-energy} H_\theta(s,v)=\frac12e^{\theta v}s^2+ c\left[\frac{e^{(\theta+1)v}-1}{\theta+1} -\frac{e^{\theta v}-1}{\theta}\right]. \end{equation} The second term has derivative $ce^{\theta v}(e^v-1)$ and its unique minimum is zero at $v=0$. Thus $H_\theta\ge0$, with equality only at $(0,0)$. For a general growth profile, \begin{equation}\label{eq:theta-dissipation} H_\theta'=e^{\theta v}s\bigl[\alpha(x(s))-C_\theta(x(s))\bigr]. \end{equation} These statements follow by direct differentiation. When $\theta=1$, \eqref{eq:theta-energy} agrees with \eqref{eq:center-integral}. Fix $N\ge1$ and $K\ge2$, and put \begin{equation}\label{eq:degree-parameters} q=\left\lceil K/2\right\rceil,\qquad d=2q+1,\qquad D=d+2N+1=2(N+q+1). \end{equation} In particular, $D$ is even, $d>K$, and $D-2>K$. Let $\mathsf H_D=\sum_{n=1}^D1/n$, and assume \begin{equation}\label{eq:height-condition} c>\theta m(\mathsf H_D-1). \end{equation} The degree-$D$ Taylor polynomial of $C_\theta$ at $a$ is \begin{equation}\label{eq:taylor-polynomial} A_D(x)=c+\theta m\sum_{n=2}^D\frac{(-1)^{n+1}}n u^n. \end{equation} It satisfies $A_D(a)=c$, $A_D'(a)=0$, and $A_D(0)=c-\theta m(\mathsf H_D-1)>0$. \begin{lemma}[Global concavity and a signed Taylor remainder] \label{lem:even-taylor} The polynomial $A_D$ is strictly concave on the entire real line. For $u>-1$, \begin{equation}\label{eq:signed-taylor-remainder} A_D(x)-C_\theta(x) =-\theta m\int_0^u\frac{t^D}{1+t}\,dt. \end{equation} The difference in \eqref{eq:signed-taylor-remainder} has sign opposite to $u$ whenever $u\ne0$. \end{lemma} \begin{proof} Differentiation gives \[ -A_D''(x)=\frac{\theta m}{a^2}S_D(u),\qquad S_D(u)=\sum_{j=0}^{D-2}(j+1)(-u)^j. \] For $u<0$ every summand is nonnegative and the constant term is positive. For $u\ge0$, the finite geometric-sum identity gives \begin{equation}\label{eq:taylor-curvature} S_D(u)=\frac{1+D u^{D-1}+(D-1)u^D}{(1+u)^2}>0. \end{equation} The derivative of $A_D-C_\theta$ with respect to $u$ is $-\theta m u^D/(1+u)$. Its value at zero is zero, proving \eqref{eq:signed-taylor-remainder}. Since $D$ is even and $1+t>0$ on the integration segment, the integral has the sign of $u$. \end{proof} \subsection{Moment interpolation and the polynomial formula} For an analytic germ $G$ at $a$, denote its Taylor polynomial of order $L$ by $T_{a,L}G$. Define \begin{equation}\label{eq:Q-j} Q_j=T_{a,D-1}\bigl[s(x)^{d+2j}\bigr],\qquad 0\le j\le N. \end{equation} The polynomial $Q_j$ has degree at most $D-1$ and vanishes at $a$ to order $d+2j$. In particular, \begin{equation}\label{eq:Q-top} Q_N=m^{(D-1)/2}u^{D-1}. \end{equation} Let $00$, set \begin{equation}\label{eq:raw-polynomial} P_\rho=A_D+\delta\sum_{j=0}^N b_j\rho^{2(N-j)}Q_j, \qquad 0\le\rho\le1. \end{equation} Its coefficientwise limit is \begin{equation}\label{eq:limit-polynomial} P_0=A_D-\delta Q_N. \end{equation} To choose $\delta$ it suffices to compute the finite constants \[ L_0=\sum_{j=0}^N|b_j|\,|Q_j(0)|,\qquad L_2=\sum_{j=0}^N|b_j| \sup_{x\ge0}\frac{|Q_j''(x)|}{-A_D''(x)}, \] and impose \begin{equation}\label{eq:delta-choice} \delta L_0<\frac{c-\theta m(\mathsf H_D-1)}2, \qquad \delta L_2<\frac12. \end{equation} The suprema are finite: their denominators are positive, and their numerators have smaller degree. They can also be bounded by coefficient estimates or by the finite critical sets of rational functions. It follows that, uniformly for $0\le\rho\le1$, \[ P_\rho(0)>0,\qquad -P_\rho''\ge-\tfrac12A_D''>0. \] Every $P_\rho$ has degree $D$, maximum $c$ at $a$, and a unique positive zero to the right of $a$; its $K$-jet equals that of $C_\theta$. \subsection{Fixing the endpoint data at the same degree} Let $b>a$ be the positive zero of $P_0$, put $z=(b-a)/a>0$, and write $E_\rho=P_\rho-P_0$. All coefficients of $E_\rho$ are $O(\rho^2)$. Define the corrected family \begin{equation}\label{eq:corrected-polynomial} \widehat P_\rho=P_\rho+t_0(\rho)u^{D-2}+t_1(\rho)u^{D-1}, \end{equation} where \begin{equation}\label{eq:endpoint-linear-system} \begin{pmatrix}1&-1\\z^{D-2}&z^{D-1}\end{pmatrix} \begin{pmatrix}t_0\\t_1\end{pmatrix} =-\begin{pmatrix}E_\rho(0)\\E_\rho(b)\end{pmatrix}. \end{equation} The determinant is $z^{D-2}(1+z)>0$, and the solution is \begin{equation}\label{eq:endpoint-correction} t_1=\frac{E_\rho(0)-E_\rho(b)/z^{D-2}}{1+z},\qquad t_0=-\frac{zE_\rho(0)+E_\rho(b)/z^{D-2}}{1+z}. \end{equation} Thus $t_0,t_1=O(\rho^2)$, and, exactly, \begin{equation}\label{eq:fixed-endpoint-data} \widehat P_\rho(0)=P_0(0),\qquad \widehat P_\rho(b)=0. \end{equation} The correction has degree at most $D-1$ and vanishes to order at least $D-2>K$ at $a$. Hence it preserves both the degree and the required jet. Each change from $P_0$ has coefficients $O(\rho^2)$ and degree at most $D-1$. Its second derivative, divided by the strictly positive polynomial $-P_0''$ of degree $D-2$, is uniformly small on $[0,\infty)$. Consequently $\widehat P_\rho''<0$ globally for all sufficiently small $\rho$. The fixed values in \eqref{eq:fixed-endpoint-data} are the data of the chosen limiting polynomial; they need not be specified independently in advance. \subsection{A uniform return-map lemma} The following elementary form of first-order averaging makes the amplitude estimates precise. The even terms are retained in an exact center before the odd perturbation is varied. \begin{lemma}[Return maps near varying reversible centers] \label{lem:uniform-return} Consider fields \[ S'=B_\nu(S)-g_\nu(V)+\epsilon\Phi_\nu(S),\qquad V'=S, \] where $B_\nu$ is even, $B_\nu(0)=g_\nu(0)=0$, $B_\nu\to0$, $g_\nu\to cV$, and $\Phi_\nu\to\Phi$ in $C^4$ on a fixed neighborhood of the compact annulus in question. For $\nu$ sufficiently small, the even field has a period annulus containing any fixed compact harmonic annulus under consideration. On the section $V=0$, $S>0$, use $h=S^2/2$. Its return displacement obeys \begin{equation}\label{eq:uniform-return} \frac{\Delta(\nu,\epsilon,h)}\epsilon \longrightarrow \oint_{(S^2+cV^2)/2=h}\Phi(S)\,dV \end{equation} in $C^1$ on compact positive energy intervals as $(\nu,\epsilon)\to0$. The assertion also holds along any sequence of the indicated fields and amplitudes. \end{lemma} \begin{proof} The scalar-flow construction in \cref{lem:even-integral}, with $g_\nu(V)$ replacing $c(e^V-1)$, yields an exact first integral $H_\nu$ normalized by $H_\nu(S,0)=S^2/2$. On the compact sets in question it converges to $(S^2+cV^2)/2$. Alternatively, solutions close to the harmonic ovals cross the fixed line $S=0$ twice, and reversibility under $(S,V)\mapsto(-S,V)$ closes each orbit. Continuous dependence and transversality give uniform bounds on return times and their derivatives. The displacement is identically zero at $\epsilon=0$. Its first derivative there is obtained by integrating $K_w S\Phi_\nu(S)$ along an unperturbed oval, with $K_w\to1$. This converges in $C^1$ to the right-hand side of \eqref{eq:uniform-return}. Taylor's theorem in the independent parameter $\epsilon$ gives a uniform $O(\epsilon^2)$ remainder. The same estimates apply to sequences; no smooth joint parametrization of competing amplitude ratios is required. \end{proof} \subsection{Exclusion at all amplitude scales} \begin{proof}[Proof of \cref{thm:exact-degree}] The construction \eqref{eq:corrected-polynomial} already proves the degree, jet, concavity, and endpoint assertions. It remains to count all its positive periodic orbits. First, \eqref{eq:signed-taylor-remainder} and \eqref{eq:Q-top} imply \[ s\bigl[P_0(x(s))-C_\theta(x(s))\bigr]<0\qquad(s\ne0). \] Thus the derivative of \eqref{eq:theta-energy} along the limiting field is strictly negative for $s\ne0$. The only invariant subset of that line is the equilibrium: there $s'=-c(e^v-1)$ and $v'=0$. The endpoint values, maximum, and positive zero of $\widehat P_\rho$ are fixed, so \cref{lem:cycle-compactness} confines every periodic orbit to one compact subset of the positive quadrant. Since $\widehat P_\rho\to P_0$ smoothly there, \cref{lem:lyapunov-confinement} shows that all possible periodic orbits converge to $E_c$ as $\rho\downarrow0$. Let $O_\rho$ be the odd part of $\widehat P_\rho(x(s))-c$. On one fixed neighborhood of zero it has the uniform analytic factorization \begin{equation}\label{eq:odd-full-expansion} \begin{split} O_\rho(s)={}&\delta\sum_{j=0}^N b_j\rho^{2(N-j)}s^{d+2j}\\ &+\rho^2s^{D-1}\mathcal A(s^2,\rho) +s^{D+1}\mathcal B(s^2,\rho), \end{split} \end{equation} where $\mathcal A,\mathcal B$ and all fixed finite derivatives are bounded for small $s,\rho$. Indeed, the difference between $Q_j(x(s))$ and $s^{d+2j}$ starts at degree $D$. Since $D$ is even, its odd part starts at degree $D+1$. The difference $A_D-C_\theta$ also starts at degree $D+1$. Finally, $u(s)=s/\sqrt m+O(s^2)$ and the correction coefficients are $O(\rho^2)$, so the odd part of the endpoint correction starts at $\rho^2s^{D-1}$. Every positive periodic orbit encloses the unique equilibrium and meets the ray $v=0$, $s>0$ exactly once: the Jordan-ray argument of \cref{lem:global-return} applies since $v'=s>0$ on this ray. Let $r$ be that intersection coordinate. The preceding confinement gives $r\to0$. We exclude the two extreme ratios of $r$ and $\rho$. \emph{Large relative amplitudes.} Suppose a sequence of cycles has $r/\rho\to\infty$. Set $s=rS$, $v=rV$; the initial point becomes $(1,0)$. Separate the even part of the scaled field and use the independent odd amplitude \[ \epsilon=\delta r^{D-2}. \] By \eqref{eq:odd-full-expansion}, its normalized odd force converges smoothly on compact sets to $b_NS^{D-1}$. The lower monomials have factors $(\rho/r)^{2(N-j)}$; the two remainder terms are bounded by constants times $\rho^2/\delta$ and $r^2/\delta$, respectively. The even field converges to the harmonic oscillator. Applying \cref{lem:uniform-return} at energy $1/2$, the displacement divided by $\epsilon$ tends to \[ b_NC_{q+N}(1/2)^{q+N+1}\ne0. \] This contradicts periodicity. \emph{Small relative amplitudes.} Suppose instead that $r/\rho\to0$. Use the same spatial scaling but take \[ \epsilon=\delta\rho^{2N}r^{d-1}. \] The normalized odd force now converges to $b_0S^d$. The higher monomials carry factors $(r/\rho)^{2j}$, while the two errors are bounded by \[ C\frac{\rho^2}{\delta}(r/\rho)^{2N},\qquad C\frac{r^2}{\delta}(r/\rho)^{2N}. \] They tend to zero. The normalized displacement tends to $b_0C_q(1/2)^{q+1}\ne0$, again a contradiction. Consequently there are constants $00$ by \eqref{eq:delta-choice}, and let $b=1+z>1$ be the positive zero of $A_6-\delta u^5$. In this case \eqref{eq:endpoint-correction} reduces to \[ t_0=\frac{5\delta\rho^2(2z-1)}{3z},\qquad t_1=-\frac{5\delta\rho^2(11z+12)}{36z}. \] The corrected polynomial has the particularly simple form \[ \widehat P_\rho=A_6-\delta u^5+ \frac{5\delta}{3}\rho^2u^3(1+u)(1-u/z). \] It is globally strictly concave for small $\rho$, has $\widehat P_\rho(0)=31/20+\delta$, positive zero $b$, and the prescribed quadratic center jet. Its system has exactly one positive periodic orbit, an attracting hyperbolic cycle with $s(\rho)=\sqrt2\rho+O(\rho^2)$. \begin{table}[t] \centering \caption{Degrees supplied by the explicit construction. Each row retains the indicated center jet and gives exactly $N$ positive periodic orbits, subject to \eqref{eq:height-condition}.} \label{tab:degree-examples} \begin{tabular}{rrrr} \toprule Number $N$ & Jet order $K$ & First free odd degree $d$ & Degree $D$\\ \midrule 1&2&3&6\\ 2&2&3&8\\ 3&2&3&10\\ 4&2&3&12\\ 3&4&5&12\\ 3&5&7&14\\ \bottomrule \end{tabular} \end{table} For a fixed $c,\theta$, condition \eqref{eq:height-condition} is a sufficient condition for this Taylor construction at the stated degree. For any fixed $c>0$ and each pair $(N,K)$, a sufficiently small $\theta>0$ satisfies it. This latter choice changes the preserved center curvature. The degree estimate is an existence upper bound for this specified jet; \cref{sec:degree} quantifies the additional costs of other interpolation and approximation data. % --- moments.tex --- \section{Prescribed oscillations and constrained Abelian integrals} \label{sec:moments} Consider the local family \begin{equation}\label{cyc:eq-family} s'=-\bigl(e^v-1\bigr)-\tfrac12s^2+F(s),\qquad v'=s, \qquad F(s)=\sum_{j=0}^{N}\lambda_j s^{2r+2j+3}, \end{equation} where $r,N\ge0$ are integers. Put $q=e^v-1$ and $Y=q+s^2/2$. Then the equations are exactly \begin{equation}\label{cyc:eq-duffing} s'=-Y+F(s),\qquad Y'=s-\tfrac12s^3+sF(s). \end{equation} The unperturbed Hamiltonian is \[ K(s,Y)=\frac{Y^2}{2}+\frac{s^2}{2}-\frac{s^4}{8}. \] Its central period annulus is $0-1$, which holds throughout this annulus. Indeed, on the oval $K=h$ one has $q\in[-\sqrt{2h},\sqrt{2h}]$. For $j\ge0$ define \[ J_j(h)=\oint_{K=h}s^{2j+3}\,dY =2\int_{-\sqrt{2h}}^{\sqrt{2h}} \left(\frac{2h-q^2}{1+q}\right)^{j+3/2}dq, \qquad 00$. For completeness, the moments also admit the convergent representation \begin{equation}\label{eq:hypergeom} J_j(h)=\pi c_jh^{j+2} {}_2F_1\!\left(\frac j2+\frac34,\frac j2+\frac54;j+3;2h\right), \qquad c_j=2^{-j-1}\binom{2j+4}{j+2}. \end{equation} Indeed, put $q=\sqrt{2h}\,t$ in the integral. Expanding $(1+\sqrt{2h}\,t)^{-j-3/2}$ and integrating term by term removes the odd terms; the beta integral gives the displayed hypergeometric series. Consequently $J_j(h)=\pi c_jh^{j+2}(1+O(h))$, with analytic remainders. Along the perturbed field, $K'=s(1+q)F(s)$. The first return displacement is therefore $\oint F(s)\,dq=\oint F(s)\,dY$ to first order, since the difference is the integral of the exact form $sF(s)\,ds$. We use the usual Abelian-integral return-map principle, as in \cite{GasullGuillamonLiZhang1996}. \begin{theorem}[A constrained consequence of the Gavrilov--Iliev theorem]\label{cyc:thm-ect} For every pair of integers $r,N\ge0$, the ordered family $(J_r,J_{r+1},\ldots,J_{r+N})$ is an extended complete Chebyshev system on $(0,1/2)$. More strongly, a nonzero real linear combination of these functions has at most $N$ zeros, counted with multiplicity, in $\mathbb C\setminus([1/2,\infty)\cup\{0\})$. \end{theorem} \begin{proof} Set $s=2x$, $Y=2y$, and $e=h/2$. The rescaled Hamiltonian is $H(x,y)=y^2+x^2-x^4$. For its central ovals, define \[ B_k(e)=-\oint_{H=e}x^{2k}y\,dx, \] using the counterclockwise orientation. Thus $B_k>0$ on $(0,1/4)$ and \begin{equation}\label{cyc:eq-rescaling} J_j(h)=(2j+3)2^{2j+4}B_{j+1}(h/2). \end{equation} Integration of $d(x^{2k+1}y^3)$ along an oval yields the recurrence \begin{equation}\label{cyc:eq-recurrence} (2k+7)B_{k+2}=(2k+4)B_{k+1}-(2k+1)eB_k, \qquad k\ge0. \end{equation} Moreover, direct differentiation gives the Picard--Fuchs system \begin{equation}\label{cyc:eq-pf} \begin{pmatrix}B_0\\ B_1\end{pmatrix} =A(e)\begin{pmatrix}B_0'\\ B_1'\end{pmatrix}, \qquad A(e)=\begin{pmatrix} 4e/3&-2/3\\ 2e/15&4e/5-4/15 \end{pmatrix}. \end{equation} Here \[ \det A(e)=\frac{4e(4e-1)}{15},\qquad \operatorname{tr}A(e)=(\det A(e))', \] and the eigenvalues of $A'$ are $4/3$ and $4/5$. We apply Theorem~1 and Proposition~4 of Gavrilov and Iliev \cite{GavrilovIliev2003}. Their parameters are $\lambda=3/4$, $\mu=5/4$, and $\lambda^*=3/4$. For $\sigma\ge3/4$, let \[ \begin{aligned} \mathcal V_\sigma=\bigl\{I(e)&=P(e)B_0(e)+Q(e)B_1(e):\ P,Q\in\mathbb R[e],\\ &I(e)=O(|e|^\sigma)\text{ sectorially at infinity}\bigr\}. \end{aligned} \] Their theorem gives at most $\dim\mathcal V_\sigma$ zeros, counted with multiplicity, in $D=\mathbb C\setminus[1/4,\infty)$ for every nonzero element of this space. Their dimension formula is \[ \dim\mathcal V_\sigma= \lfloor\sigma-3/4\rfloor+\lfloor\sigma-5/4\rfloor+2. \] All hypotheses are displayed in \eqref{cyc:eq-pf}: there are two distinct finite singular values $0,1/4$, and the oval integrals are analytic at $0$. For $L\ge0$, put $\sigma_L=(2L+3)/4$. The Fuchsian exponents at infinity give $B_0=O(|e|^{3/4})$ and $B_1=O(|e|^{5/4})$. Recurrence \eqref{cyc:eq-recurrence} therefore implies $B_k=O(|e|^{(2k+3)/4})$. Consequently \[ \operatorname{span}(B_0,\ldots,B_L)\subset\mathcal V_{\sigma_L}. \] Since $B_k(e)=b_ke^{k+1}(1+O(e))$ with $b_k>0$, these $L+1$ functions are linearly independent. The displayed dimension formula gives $\dim\mathcal V_{\sigma_L}=L+1$, so the inclusion is equality. Now take $L=r+N+1$. A linear combination of $B_{r+1},\ldots,B_{r+N+1}$ is an element of $\mathcal V_{\sigma_L}$ and vanishes to order at least $r+2$ at $e=0$. The complex zero bound is $L+1=r+N+2$, leaving at most $N$ zeros away from the origin. Rescaling by \eqref{cyc:eq-rescaling} proves the stronger assertion. Applying the same argument to every prefix proves the ECT property. \end{proof} \paragraph{Elementary local alternative.} The global theorem uses the cited complex nonoscillation result. A fully elementary local statement follows from \[ J_j(h)=\pi c_jh^{j+2}(1+O(h)),\qquad c_j=2^{-j-1}\binom{2j+4}{j+2}. \] For $0\le\ell\le N$, \[ W(J_r,\ldots,J_{r+\ell})(h)= \left(\prod_{j=0}^{\ell}\pi c_{r+j}\right) \left(\prod_{j=1}^{\ell}j!\right) h^{(\ell+1)(r+2)}(1+O(h)). \] Thus all Wronskians are positive on a sufficiently short positive interval. This alternative is useful when only local analyticity of a more general restoring force is available. \begin{corollary}[Interpolation of oscillation amplitudes]\label{cyc:cor-interpolation} Fix $r,N\ge0$ and distinct energies $00$, all derivatives of the prey growth function through order $2r+2$ at $a$ agree exactly with those of the fixed center profile $1-s(x)^2/2$. In particular its value and negative second derivative are fixed independently of the number of cycles. Taking $r$ large preserves any prescribed finite jet of this center profile while realizing arbitrarily many cycles at arbitrary energies in the period annulus. The prescribed jet in this statement belongs to the fixed center profile. \Cref{sec:localized} treats arbitrary nondegenerate maximum profiles by localized smooth perturbations. \begin{remark} The cyclicity bound concerns regular ovals on a compact subannulus in the specified finite-dimensional family. The boundary two-saddle polycycle and the dynamics outside the globally extended local model are separate bifurcation problems. Prescribed positions throughout the whole annulus are a statement about the local canonical family. A simultaneous globally concave realization with $\alpha(0)>0$ uses sufficiently small ovals. For example, the left tangent intercept of the center profile at $\ell0$. \subsection{A concave extension with a dissipative exterior} Keep $k,m,\beta,c>0$ fixed and put \[ u(x)=\frac{s(x)}{\sqrt c} =\operatorname{sgn}(x-a) \sqrt{\frac{2m}{c}\left(\frac xa-1-\log\frac xa\right)}. \] After the positive time change $d\sigma=\sqrt c\,d\tau$, the center profile \eqref{eq:center-profile} becomes \[ u'=-(e^v-1)-\frac12u^2,\qquad v'=u. \] Its normalized first integral is \begin{equation}\label{eq:extension-energy} H(u,v)=\frac{H_c(\sqrt c\,u,v)}c =\frac12e^vu^2+\frac12(e^v-1)^2. \end{equation} Thus the energies below are exactly the normalized energies used in \cref{sec:moments}. \begin{lemma}[A strictly dissipative concave extension] \label{lem:dissipative-extension} For every sufficiently small $b>0$ there exists a globally strictly concave admissible profile $\alpha_b$ such that \[ \alpha_b(x)=\alpha_{\mathrm{cen}}(x) \quad\text{when }|u(x)|\le b, \] and \begin{equation}\label{eq:extension-order} \alpha_b(x)>\alpha_{\mathrm{cen}}(x)\quad(u(x)<-b), \qquad \alpha_b(x)<\alpha_{\mathrm{cen}}(x)\quad(u(x)>b). \end{equation} It has its maximum $c$ at $a$ and satisfies $\alpha_b(0)>0$. If $0\frac m{x^2}\quad(x>R). \] It can be chosen small and constant near zero, and positive and constant for all sufficiently large $x$. Smooth flat transitions at $L$ and $R$ give the stated strict inequalities off the matching interval. On the left, first take $b$ small enough that $\int_L^a m/x\,dx0$. Its positive right-hand zero exists and is unique because the profile is strictly concave. Twice integrating the curvature comparisons from $a$ proves both exact matching on $[L,R]$ and \eqref{eq:extension-order}. In the normalized coordinates, differentiation of \eqref{eq:extension-energy} gives \begin{equation}\label{eq:extension-Lyapunov} H'=e^vu\,\frac{\alpha_b(x(u))-\alpha_{\mathrm{cen}}(x(u))}{c} \le0, \end{equation} with strict inequality whenever $|u|>b$. Every periodic orbit must therefore remain in the strip $|u|\le b$, where the field agrees exactly with the center. To identify the compact invariant part of that strip, put $q=e^v-1$ and $p=e^{v/2}u$. Then \[ H=\frac{p^2+q^2}{2},\qquad q>-1, \qquad (p',q')=\sqrt{1+q}\,(-q,p). \] The center ovals $00$ such that $2h_00$, with energy coordinate $h=(e^v-1)^2/2$. On $h_0\le h\le h_b$, the unperturbed first return is the identity, and smooth dependence gives \begin{equation}\label{eq:extension-return-expansion} \Delta(h,\varepsilon) :=P_\varepsilon(h)-h =\varepsilon M(h)+O(\varepsilon^2) \quad\text{in }C^1. \end{equation} The variation formula is the one used in \cref{sec:moments}; the cutoff is identically one on a neighborhood of every unperturbed oval under consideration. Around each simple zero $h_i$, the implicit function theorem gives one and only one zero of $\Delta$ for small positive $\varepsilon$. On the compact complement of these neighborhoods, $M$ is bounded away from zero, so \eqref{eq:extension-return-expansion} gives no other zeros. There are exactly $N$ simple fixed points on this annulus. \emph{A collar of the last center oval.} The orbit $H=h_b$ is regular. Its return map extends smoothly to a fixed collar $h_b\le h\le h_b+\zeta$, with $h_b+\zeta<1/2$, and to small positive and negative $\varepsilon$. For the baseline, the Lyapunov inequality \eqref{eq:extension-Lyapunov} gives $\Delta(h,0)\le0$ throughout this collar. At its inner boundary, \[ \partial_\varepsilon\Delta(h_b,0)=M(h_b)<0. \] By continuity, after reducing $\zeta$ and the allowed amplitude, $\partial_\varepsilon\Delta<0$ on the whole collar. Thus $\Delta(h,\varepsilon)<0$ for every sufficiently small positive $\varepsilon$ and every $h$ in the collar. There are no periodic orbits with section coordinate there. \emph{The remaining exterior.} The profiles converge smoothly to $\alpha_b$ and retain $\alpha_b(0)$, $c$, and $x_{1,b}$, so \cref{lem:cycle-compactness} gives a common compact set for all their periodic orbits. By \cref{lem:dissipative-extension,lem:lyapunov-confinement}, every such periodic orbit lies in $\{Hx_{1,b}$ and integers $J\ge d-1$, $L\ge d$, there are globally strictly concave polynomials converging to $\alpha_\varepsilon$ in $C^L([0,R])$ and satisfying exactly \[ j_a^J\alpha_n=j_a^J\alpha_\varepsilon, \qquad \alpha_n(0)=\alpha_b(0), \qquad \alpha_n(x_{1,b})=0, \] such that, for all sufficiently large $n$, their positive periodic orbits consist of exactly $N$ hyperbolic cycles. Consequently, for every finite jet of the specified center profile \eqref{eq:center-profile}, with arbitrary fixed $c>0$, an exact global polynomial realization exists. This statement asserts finiteness of the required degree; the explicit linear degree estimate is furnished separately by \cref{thm:exact-degree} under its stated height condition. \end{corollary} \begin{proof} Apply \cref{prop:polynomial-interpolation} with the stated jet and endpoint constraints. The resulting polynomials are globally strictly concave, have their maximum $c$ at $a$, and have the same positive zero $x_{1,b}$. They converge in $C^d$ because $L\ge d$, and retain the jet through $d-1$ because $J\ge d-1$. The target profile has exactly $N$ hyperbolic periodic orbits and a nonzero canonical odd coefficient of degree $d$ by \cref{thm:exact-extension}. All hypotheses of \cref{prop:exact-count-stability} therefore hold, proving that every sufficiently accurate polynomial approximant has exactly $N$ positive periodic orbits, all hyperbolic. To retain the jet of \eqref{eq:center-profile} through any chosen finite order $K$, take an odd $d>K$ in \cref{thm:exact-extension} and then choose $J\ge d-1$. The construction works for every $c>0$ because the matching interval of the smooth extension can be chosen sufficiently small. \end{proof} The differentiability requirement in this transfer has a specific purpose. Mere $C^1$ closeness preserves the already hyperbolic cycles, but does not exclude new cycles approaching a degenerate focus. The $C^d$ estimate together with the fixed lower jet preserves the leading odd sign and supplies the missing local barrier. The approximation retains the specified finite jet. % --- density.tex --- \section{Localized creation of cycles with preservation of the germ} \label{sec:localized} Consider the canonical system \begin{equation}\label{eq:canonical-density} \dot s=f(s)-c(e^v-1),\qquad \dot v=s,\qquad c>0, \end{equation} where $f$ is smooth near zero and \[ f(0)=f'(0)=0,\qquad f''(0)<0. \] Write $f=E+O$, where $E$ is even and $O$ is odd. Suppose first that \begin{equation}\label{eq:first-odd-term} O(s)=b s^d+O(s^{d+2}),\qquad b\ne0,\qquad d=2r+1\ge3. \end{equation} \begin{theorem}[Localized cycle creation]\label{thm:localized-density} Let $N\ge1$ and let $00$, the function \[ f_\rho(s)=f(s)+\rho^dR(s/\rho) \] has the following properties. \begin{enumerate} \item $f_\rho=f$ both on a neighborhood of zero and outside an interval $|s|0$, their intersection coordinates satisfy \[ s_i(\rho)=\rho\sqrt{2h_i}+O(\rho^2). \] \item Their multipliers have the form \[ \mu_i(\rho)=1+\rho^{d-1}\bigl(M_0'(h_i)+O(\rho)\bigr), \] where $M_0'(h_i)\ne0$. The signs of these derivatives alternate; the first has sign $-\operatorname{sgn}(b)$. \end{enumerate} \end{theorem} \begin{proof} We first construct $R$. For the harmonic oscillator \[ \dot S=-cV,\qquad \dot V=S, \] let $\gamma_h$ denote the oval $(S^2+cV^2)/2=h$. Its odd moments are \begin{equation}\label{eq:harmonic-moments} \oint_{\gamma_h}S^{2j+1}\,dV=K_jh^{j+1},\qquad K_j=\frac{\pi}{\sqrt c\,2^j}\binom{2j+2}{j+1}>0. \end{equation} Choose a compact interval $I\Subset(0,\infty)$ containing all $h_i$ in its interior. There are unique numbers $a_1^0,\ldots,a_N^0$ such that \begin{equation}\label{eq:target-melnikov-polynomial} bK_rh^{r+1}+\sum_{\ell=1}^Na_\ell^0K_{r+\ell}h^{r+\ell+1} =bK_rh^{r+1}\prod_{i=1}^N(1-h/h_i). \end{equation} Take a smooth even cutoff $\chi$ that equals one on an interval containing the $S$-projections of all $\gamma_h$, $h\in I$, and has compact support. Take another smooth even cutoff $\theta_\kappa$ that vanishes for $|S|\le\kappa$ and equals one for $|S|\ge2\kappa$. The matrix \[ A_{i\ell}(\kappa) =\oint_{\gamma_{h_i}} \theta_\kappa(S)S^{2(r+\ell)+1}\,dV \] converges as $\kappa\downarrow0$ to $K_{r+\ell}h_i^{r+\ell+1}$. The latter matrix is invertible: after factoring positive row and column factors, it is a Vandermonde matrix. Consequently, for all sufficiently small fixed $\kappa$, there are unique $a_\ell(\kappa)$, converging to $a_\ell^0$, for which \[ M_0(h_i)=0,\qquad M_0(h):=\oint_{\gamma_h} \left[bS^d+\theta_\kappa(S) \sum_{\ell=1}^Na_\ell(\kappa)S^{2(r+\ell)+1}\right]dV. \] Moreover, $M_0$ and its first derivative converge uniformly on $I$ to the right-hand side of \eqref{eq:target-melnikov-polynomial}. This follows directly from the parametrization \[ S=\sqrt{2h}\cos\vartheta,\qquad V=\sqrt{2h/c}\sin\vartheta. \] The cutoff error is supported where $|S|\le2\kappa$; differentiation in $h$ remains harmless because $h$ is bounded away from zero. In particular, the prescribed zeros remain simple and retain the derivative signs in \eqref{eq:target-melnikov-polynomial}. Fix such a $\kappa$ and put \[ R(S)=\chi(S)\theta_\kappa(S) \sum_{\ell=1}^Na_\ell(\kappa)S^{2(r+\ell)+1}. \] Now rescale $s=\rho S$, $v=\rho V$. The modified system becomes \[ \dot S=B_\rho(S,V)+\rho^{d-1}D_\rho(S),\qquad \dot V=S, \] where \[ B_\rho(S,V)=\frac{E(\rho S)}\rho -c\frac{e^{\rho V}-1}\rho, \qquad D_\rho(S)=\frac{O(\rho S)}{\rho^d}+R(S). \] Both expressions extend smoothly to $\rho=0$, with \[ B_0(S,V)=-cV,\qquad D_0(S)=bS^d+R(S). \] For every sufficiently small $\rho$, the unperturbed field $(B_\rho,S)$ is reversible under $(S,V)\mapsto(-S,V)$ and has an elliptic linearization at the origin. Its local orbits are periodic. One elementary proof uses the two successive intersections of a near-harmonic orbit with the fixed line $S=0$; reversibility closes the orbit. Equivalently, the scalar-flow construction of $K$ in \cref{lem:even-integral} gives a smooth first integral $H_\rho$ normalized by $H_\rho(S,0)=S^2/2$; as $\rho\to0$ it converges smoothly on the required compact set to $(S^2+cV^2)/2$. Consequently, on any fixed compact annulus containing the ovals indexed by $I$, the periodic-orbit construction and the transversality of the section are uniform in $\rho$. Introduce an independent perturbation parameter $\epsilon$ and write \[ \dot S=B_\rho(S,V)+\epsilon D_\rho(S),\qquad \dot V=S. \] On the section $V=0$, $S>0$, use $h=S^2/2$. Its displacement is smooth in $(\rho,\epsilon,h)$ and vanishes identically when $\epsilon=0$. Uniformly in $C^1(I)$, \[ \Delta(\rho,\epsilon,h) =\epsilon\bigl(M_\rho(h)+O(\epsilon)\bigr), \qquad M_\rho=M_0+O(\rho). \] At $\rho=0$, the formula for $M_0$ follows by integrating the derivative of $(S^2+cV^2)/2$ along $\gamma_h$. Substituting $\epsilon=\rho^{d-1}$ gives \[ \rho^{1-d}\Delta(\rho,\rho^{d-1},h)=M_0(h)+O(\rho) \quad\hbox{in }C^1(I). \] The implicit function theorem yields one hyperbolic cycle near each $h_i$, with $h_i(\rho)=h_i+O(\rho)$ and the stated multipliers. Finally, for each $0\le j\le d-1$, \[ \left\|\frac{d^j}{ds^j}(f_\rho-f)\right\|_\infty \le \rho^{d-j}\|R^{(j)}\|_\infty. \] This proves the norm estimate. Because $R$ vanishes near zero, $R'(S)/S$ is bounded, and \[ \frac{f_\rho'(s)-f'(s)}s =\rho^{d-2}\frac{R'(s/\rho)}{s/\rho}\longrightarrow0 \] uniformly on the support. Since $f'(s)/s\to f''(0)<0$, the derivative sign is preserved. All germ and support assertions are immediate. \end{proof} \begin{corollary}[Strictly concave prey growth]\label{cor:concave-density} Fix $k,m,\beta>0$ and set $a=m/(k\beta)$. Let $\alpha$ be any admissible smooth prey growth function with maximum at $a$ and $\alpha(a)>0$. If its canonical odd part has a first nonzero term of degree $d$, then for every $N$ there are growth functions $\alpha_\rho$ with at least $N$ hyperbolic positive limit cycles, equal to $\alpha$ near $a$ and outside a shrinking neighborhood of $a$, and converging to $\alpha$ in $C^{d-1}$. In particular, the carrying capacity, the value at zero, the tail, and the entire germ at the maximum are unchanged. If $\alpha''<0$ globally, the functions $\alpha_\rho$ may also be chosen globally strictly concave. Hence profiles arbitrarily close in $C^2$ to an admissible quadratic profile $\alpha(x)=c-b(x-a)^2$, $c>ba^2$, can have arbitrarily many hyperbolic cycles. \end{corollary} \begin{proof} Use the smooth local coordinate $s(x)=\operatorname{sgn}(x-a)\sqrt{2m(x/a-1-\log(x/a))}$ and put $f(s)=\alpha(x(s))-\alpha(a)$. Add $\rho^dR(s(x)/\rho)$ to $\alpha$. The coordinate change preserves the germ, support, and norm conclusions. Since $d\ge3$, the perturbation is $C^2$-small. Strict concavity is preserved on its shrinking support, where the original second derivative is bounded above by a negative constant; outside the support the growth function is unchanged. For the quadratic profile the first canonical odd term has degree three: \[ f(s)=-\frac{ba^2}{m}s^2 -\frac{2ba^2}{3m^{3/2}}s^3+O(s^4). \] \end{proof} \section{The sharp differentiability threshold} \label{sec:threshold} \begin{proposition}\label{prop:sharp-topology} Under \eqref{eq:first-odd-term}, there are $s_0>0$, a neighborhood $U$ of the origin, and $\epsilon_0>0$ such that every smooth $\widetilde f$ satisfying \[ j_0^{d-1}\widetilde f=j_0^{d-1}f, \qquad \|\widetilde f-f\|_{C^d([-s_0,s_0])}<\epsilon_0 \] has no periodic orbit contained in $U$ for its canonical system. Thus the $C^{d-1}$ approximation assertion cannot generally be upgraded to $C^d$ while preserving the indicated jet. \end{proposition} \begin{proof} Shrink the $s$-interval so that $O(s)/s^d$ has the sign of $b$ and absolute value at least $|b|/2$. Taylor's theorem and the jet condition give \[ |(\widetilde f-f)_{\rm odd}(s)| \le \frac{\|\widetilde f-f\|_{C^d}}{d!}|s|^d. \] For a sufficiently small norm bound, the odd part $\widetilde O$ therefore has the sign of $bs$ for $s\ne0$. Apply \cref{lem:even-integral} to the even part $\widetilde E$. Its local function $H_{\widetilde E}$ satisfies \[ \dot H_{\widetilde E}=K_w\,s\widetilde O(s),\qquad K_w>0. \] This has the fixed sign of $b$ away from $s=0$, and its integral over a nonconstant periodic orbit cannot vanish. The scalar flow defining $K$ exists on a common neighborhood under the stated small $C^d$ bound, by continuous dependence and uniform $C^2$ bounds on the even part. \end{proof} \subsection{The flat odd-part case} If $O$ is flat at zero, then for every finite $k$ and every $N$ one can still create at least $N$ hyperbolic cycles while preserving the entire germ and making $\|f_\rho-f\|_{C^k}\to0$. Choose an odd integer $d=2r+1>\max\{k,2\}$ and use the same outer and inner cutoffs. This time the leading monomial must also be inserted: \[ R_\kappa(S)=\chi(S)\theta_\kappa(S) \sum_{\ell=0}^N a_\ell(\kappa)S^{d+2\ell}, \qquad a_0(\kappa)=1. \] Solve the $N$ equations $\oint_{\gamma_{h_i}}R_\kappa(S)\,dV=0$ for $a_1(\kappa),\ldots,a_N(\kappa)$. At $\kappa=0$ their coefficient matrix is the invertible Vandermonde-type matrix used above, and the target displacement is $K_rh^{r+1}\prod_i(1-h/h_i)$. For small fixed $\kappa$ the coefficients converge, the prescribed zeros remain simple, and their derivative signs are unchanged. Now $O(\rho S)/\rho^d\to0$ smoothly on compact sets, including smooth extension in $\rho$ at zero by flatness. The limiting displacement is $\oint R_\kappa(S)\,dV$. The same two-parameter return-map proof gives the cycles, and the norm estimate is $O(\rho^{d-k})$. Taking $k$ successively larger and choosing the corresponding scales sufficiently small gives a sequence converging in the smooth Fr\'echet topology, with any prescribed sequence of finite lower bounds on cycle counts. \begin{proof}[Proof of \cref{thm:main-density}] If the canonical odd part has a first nonzero Taylor coefficient, \cref{thm:localized-density,cor:concave-density} apply, and $d\ge3$ gives $C^2$ convergence. If all its Taylor coefficients vanish, apply the flat case with $k=2$. The perturbations vanish near the maximum and outside their shrinking support; hence the low-density value, positive zero, and the germ are exactly preserved. On the support, strict concavity follows from $C^2$ convergence and the negative curvature of the original profile near its maximum. Outside the support the profile is unchanged. \end{proof} \begin{remark} Germ preservation means that, for each member of the approximating family, there is a neighborhood on which the functions agree. The size of that neighborhood shrinks with the perturbation scale. In the finite-order odd case, the family converges in $C^{d-1}$ and need not converge in $C^d$. These are the topology and locality features used by the realization theorem. \end{remark} % --- global_realization.tex --- \section{Global concavity and polynomial realization} \label{sec:global} Local oscillation constructions can be extended without changing the vector field on the selected compact cycles. The following elementary extension lemma will also be used for nonlinear responses. \begin{lemma}[Concave extension]\label{lem:concave-extension} Let a smooth local function $A$ be given near $a>0$, with $A(a)>0$, $A''<0$ locally, and \[ L:=A(a)-aA'(a)>0. \] After restricting its local domain if necessary, it extends to a smooth function on $[0,\infty)$ that is strictly concave, positive at zero, has a unique nondegenerate maximum at a positive point, and has a unique positive zero to the right of that maximum. \end{lemma} \begin{proof} Write $A_0=A(a)$, $D=A'(a)$, and prescribe a global smooth weight $w>0$ agreeing with $-A''$ on a small interval around $a$. Define \[ \widetilde A(x)=A_0+D(x-a)-\int_a^x(x-t)w(t)\,dt. \] Then $\widetilde A=A$ locally and $\widetilde A''=-w<0$. The two conditions needed on the left extension are \[ \int_0^a t w(t)\,dt0. \] The first follows by restricting the prescribed interval about $a$ and taking a sufficiently small positive background elsewhere. If $D<0$, the second can be obtained by adding a positive bump of mass larger than $-D$ close to zero. Its contribution to the first, weighted integral can be made arbitrarily small by placing its support sufficiently close to zero. Smooth positive interpolation at the ends of the prescribed interval does not affect these strict inequalities. On the right choose $w$ eventually to be a positive constant. Thus $\widetilde A(0)>0$, $\widetilde A'(0)>0$, $\widetilde A'(x)\to-\infty$, and $\widetilde A(x)\to-\infty$. Strict concavity gives the unique nondegenerate maximum and the claimed zero. \end{proof} Taking $A=\alpha_{\mathrm{cen}}$ from \eqref{eq:center-profile} gives $A(a)-aA'(a)=c>0$. Thus there are admissible globally strictly concave profiles that agree exactly with this integrable model on a sufficiently small neighborhood of $a$. \subsection{Polynomial approximation with exact constraints} \begin{proposition}\label{prop:polynomial-interpolation} Let $\alpha\in C^\infty([0,\infty))$ be strictly concave, with $\alpha(a)=c>0$, $\alpha'(a)=0$, $\alpha(0)>0$, and $\alpha(x_1)=0$ for $x_1>a$. Fix $R>x_1$, a finite jet order $K\ge2$, and a finite approximation order $L$. Then $\alpha$ can be approximated on $[0,R]$ in $C^L$ by polynomials $\alpha_n$ satisfying, exactly, \[ \alpha_n''(x)<0\quad(x\ge0),\qquad j_a^K\alpha_n=j_a^K\alpha,\qquad \alpha_n(0)=\alpha(0),\qquad \alpha_n(x_1)=0. \] In particular, any finite collection of hyperbolic positive limit cycles is preserved by a sufficiently accurate such approximation. \end{proposition} \begin{proof} Put $w=-\alpha''>0$ on $[0,R]$, choose $0K-2$. Introduce \[ p_{n,t}(x)=p_n(x) \bigl[1+t_0(x-a)^{2M}+t_1(x-a)^{2M+1}\bigr], \qquad w_{n,t}=p_{n,t}^2+d. \] These weights are globally positive polynomials and retain the required jet at $a$. The two equations to be imposed are \[ \int_0^a xw_{n,t}(x)\,dx=c-\alpha(0),\qquad \int_a^{x_1}(x_1-x)w_{n,t}(x)\,dx=c. \] At $t=0$, the Jacobian of the two left-hand sides has rows $(A_-,B_-)$ and $(A_+,B_+)$, where \[ A_->0,\quad B_-<0,\qquad A_+>0,\quad B_+>0. \] Indeed, apart from positive weights, the entries are integrals of $2p_n^2(x-a)^{2M}$ and $2p_n^2(x-a)^{2M+1}$ over the left and right intervals, respectively. Thus the determinant is strictly positive. View the two moments as a smooth map of $(p,t)$, with $p$ in a neighborhood of $q$ in $C([0,R])$. At $(q,0)$ the target moments hold and the displayed Jacobian is invertible. The implicit function theorem with $p$ as a Banach-space parameter gives a continuous solution $t(p)$ near $q$, with $t(q)=0$. Hence $t_n=t(p_n)\to0$ solves both moment equations exactly. Define \[ \alpha_n(x)=c-\int_a^x(x-u)w_{n,t_n}(u)\,du. \] Then $\alpha_n''=-w_{n,t_n}<0$ globally, $\alpha_n(a)=c$, and $\alpha_n'(a)=0$. The weight jet gives the full requested growth-rate jet, and the two moment equations give the two endpoint values. The differentiable approximation follows from the construction. Every positive periodic orbit of the target system lies in $00$. Thus $C^1$ approximation on this interval suffices for the persistence of its finitely many hyperbolic cycles. \end{proof} \begin{proof}[Proof of \cref{thm:main-polynomial}] First apply \cref{thm:main-density} to the prescribed profile, taking the perturbation small enough that the resulting strictly concave profile $\alpha_* $ is within $\varepsilon/2$ of $\alpha$ in $C^2([0,R])$. It has $N$ hyperbolic cycles, agrees with $\alpha$ near $a$, and has the same values at $0$ and $x_1$. Apply \cref{prop:polynomial-interpolation} to $\alpha_*$, retaining the jet through order $K$ and the two endpoint values. Choose the approximation within $\varepsilon/2$ in $C^2$ and within the $C^1$ persistence neighborhood of all $N$ cycles. The resulting polynomial has all the asserted properties. \end{proof} \begin{corollary}[Fixed quadratic data]\label{cor:quadratic-polynomial} Fix $b,c>0$ with $c>ba^2$, and let $\alpha_0(x)=c-b(x-a)^2$. For every $N\ge1$ and $K\ge2$ there are globally strictly concave polynomials with at least $N$ hyperbolic cycles such that \[ \alpha(0)=c-ba^2,\quad \alpha(a)=c,\quad \alpha''(a)=-2b,\quad \alpha^{(j)}(a)=0\ (3\le j\le K), \quad x_1=a+\sqrt{c/b}. \] They can be chosen arbitrarily close to $\alpha_0$ in $C^2$ on any prescribed compact interval containing $[0,x_1]$. \end{corollary} \begin{remark} Hyperbolicity also permits a small change in the predator mortality or the linear predation coefficient after the construction. The cycles persist while the break-even density moves away from the growth-rate maximum. Thus the resulting oscillations persist in systems with a hyperbolic coexistence equilibrium as well. \end{remark} % --- quantitative_degree.tex --- \section{Cycle number, retained jets, and polynomial degree} \label{sec:degree} The degree in \cref{thm:exact-degree} is explicit and linear in the number of cycles and the retained jet order. Quantitative approximation of a prescribed growth profile involves additional data. This section separates these two questions, gives necessary degree bounds, and makes the approximation in \cref{prop:polynomial-interpolation} effective. Throughout this section, \[ \norm{g}_{C^j(I)}=\max_{0\le i\le j}\sup_{x\in I}|g^{(i)}(x)|. \] \subsection{The degree forced by a center jet} \begin{proposition}\label{prop:center-jet-degree} Fix $\theta>0$ and $K\ge2$. Any polynomial $P$ satisfying $P''<0$ on $[0,\infty)$ and $j_a^KP=j_a^KC_\theta$, where $C_\theta$ is the center profile in \eqref{eq:theta-center}, has degree at least \begin{equation}\label{eq:center-jet-degree} D_{\mathrm{jet}}=2\left\lceil\frac K2\right\rceil. \end{equation} If $c>\theta m(\mathsf H_{D_{\mathrm{jet}}}-1)$, this lower bound is attained by an admissible globally strictly concave polynomial. Here $\mathsf H_j=\sum_{i=1}^j1/i$. \end{proposition} \begin{proof} For $u=(x-a)/a$, the center germ is \[ C_\theta(x)=c+\theta m \sum_{n=2}^{\infty}\frac{(-1)^{n+1}}n u^n. \] Its derivative of every order $n\ge2$ at $a$ is nonzero. Hence $\deg P\ge K$. If $K$ is odd and $\deg P=K$, the prescribed jet determines $P$ completely, and its leading coefficient is positive. Then $P''$ is positive for all sufficiently large $x$, contrary to the hypothesis. This proves \eqref{eq:center-jet-degree}. For an even integer $d$, the Taylor polynomial \[ T_d(x)=c+\theta m\sum_{n=2}^d\frac{(-1)^{n+1}}n u^n \] satisfies \[ -T_d''(x)=\frac{\theta m}{a^2} \sum_{j=0}^{d-2}(j+1)(-u)^j>0. \] For $u<0$ every summand is nonnegative and the constant term is positive. For $u\ge0$, the sum is \[ \frac{1+d u^{d-1}+(d-1)u^d}{(1+u)^2}>0. \] Moreover, $T_d(a)=c$, $T_d'(a)=0$, and $T_d(0)=c-\theta m(\mathsf H_d-1)$. Taking $d=D_{\mathrm{jet}}$ proves the attainment assertion under the stated height condition. \end{proof} Consequently, the explicit degree \eqref{eq:degree-bound} is $D_{\mathrm{jet}}+2N+2$. This is the degree cost of the construction in \cref{thm:exact-degree}; it does not assert that this degree is optimal among all profiles with $N$ cycles. The height condition also has a precise role. For fixed $c>0$ one may choose $\theta$ sufficiently small for each $N,K$, whereas fixing $\theta$ fixes the curvature of the center jet and retains the condition $c>\theta m(\mathsf H_D-1)$ in that construction. Its two fixed endpoint values are the values determined by the limiting polynomial in \cref{thm:exact-degree}. Approximation with arbitrary prescribed admissible endpoint data is covered below by a data-dependent bound. \subsection{A sharp lower bound from the left endpoint} \begin{proposition}[The curvature budget at the left endpoint] \label{prop:endpoint-degree-budget} Let $P$ have degree $D\ge2$, and suppose \[ P''(x)<0\quad(x\ge0),\qquad P(a)=c,\quad P'(a)=0,\quad P''(a)=-\kappa<0, \quad P(0)=\alpha_0>0. \] Then \begin{equation}\label{eq:endpoint-degree-budget} \frac{c-\alpha_0}{\kappa a^2}\ge b_D, \qquad b_D= \begin{cases} \displaystyle\frac4{D(D+2)},&D\text{ even},\\[4pt] \displaystyle\frac4{(D-1)(D+3)},&D\text{ odd}. \end{cases} \end{equation} In particular, \begin{equation}\label{eq:endpoint-degree-inverse} D\ge\sqrt{1+\frac{4\kappa a^2}{c-\alpha_0}}-1. \end{equation} The constants $b_D$ are sharp infima under the stated shape, value, slope, and curvature conditions, with the endpoint value allowed to vary. For the center curvature $\kappa=\theta m/a^2$, the denominator on the left of \eqref{eq:endpoint-degree-budget} is $\theta m$. \end{proposition} \begin{proof} Set $w=-P''$ and $q(t)=w(at)/\kappa$. Then $q>0$ on $[0,\infty)$, $q(1)=1$, and $\deg q=d=D-2$. Integrating the curvature twice gives \begin{equation}\label{eq:normalized-curvature-budget} \frac{c-\alpha_0}{\kappa a^2}=\int_0^1tq(t)\dd t. \end{equation} We first establish two elementary polynomial facts. Every real polynomial nonnegative on $[0,\infty)$ can be written as \begin{equation}\label{eq:halfline-squares} q(t)=A(t)^2+tB(t)^2, \qquad \deg A\le\lfloor d/2\rfloor, \quad \deg B\le\lfloor(d-1)/2\rfloor. \end{equation} To see this, factor $q$ over the reals. Its positive real roots have even multiplicity. Apart from squares and a positive constant, the factors are powers of $t$, factors $t+\lambda$ with $\lambda>0$, and irreducible quadratics $(t-u)^2+v^2$ with $v\ne0$. The first two types have the required form. For the last, put $\rho=\sqrt{u^2+v^2}$ and use \[ (t-u)^2+v^2=(t-\rho)^2+2(\rho-u)t. \] The class is closed under multiplication because \[ (A^2+tB^2)(C^2+tE^2)=(AC-tBE)^2+t(AE+BC)^2. \] The degree bounds in \eqref{eq:halfline-squares} follow because the leading coefficients of the two nonzero summands are positive. For $r=1,2$ and any polynomial $p$ of degree at most $n$, one has \begin{equation}\label{eq:weighted-endpoint-polynomial} \int_0^1t^rp(t)^2\dd t \ge\frac{p(1)^2}{(n+1)(n+r+1)}. \end{equation} Indeed, the Rodrigues polynomials \[ L_j^{(r)}(t)=\frac1{j!t^r} \frac{d^j}{dt^j}\bigl[t^{j+r}(t-1)^j\bigr],\qquad 0\le j\le n, \] have degree $j$, value one at $t=1$, and are orthogonal for the weight $t^r$ on $[0,1]$. Integration by parts $j$ times, followed by the beta integral, gives \[ \int_0^1t^r\bigl(L_j^{(r)}(t)\bigr)^2\dd t =\frac1{2j+r+1}. \] Expanding $p=\sum_{j=0}^n a_jL_j^{(r)}$ and applying Cauchy--Schwarz therefore yields \[ p(1)^2\le \left(\sum_{j=0}^n\frac{a_j^2}{2j+r+1}\right) \left(\sum_{j=0}^n(2j+r+1)\right) =(n+1)(n+r+1)\int_0^1t^rp(t)^2\dd t. \] If $d=2n$ with $n\ge1$, \eqref{eq:halfline-squares} and \eqref{eq:weighted-endpoint-polynomial} give \begin{align*} \int_0^1tq(t)\dd t &\ge\frac{A(1)^2}{(n+1)(n+2)} +\frac{B(1)^2}{n(n+2)}\\ &\ge\frac{A(1)^2+B(1)^2}{(n+1)(n+2)} =\frac1{(n+1)(n+2)}. \end{align*} For $n=0$, $q=1$ and the same conclusion follows directly. If $d=2n+1$, both $A$ and $B$ have degree at most $n$, so \[ \int_0^1tq(t)\dd t \ge\frac{A(1)^2}{(n+1)(n+2)} +\frac{B(1)^2}{(n+1)(n+3)} \ge\frac1{(n+1)(n+3)}. \] Substitution of $d=D-2$ proves \eqref{eq:endpoint-degree-budget}. Since both parity-dependent bounds are at least $4/[D(D+2)]$, \eqref{eq:endpoint-degree-inverse} follows. For sharpness, define the endpoint reproducing kernel \[ K_n^{(r)}(t,1)=\sum_{j=0}^n(2j+r+1)L_j^{(r)}(t). \] Take, respectively, for even and odd $d$, \[ q(t)=\left[\frac{K_n^{(1)}(t,1)}{K_n^{(1)}(1,1)}\right]^2, \qquad q(t)=t\left[\frac{K_n^{(2)}(t,1)}{K_n^{(2)}(1,1)}\right]^2. \] These polynomials are nonnegative on the entire half-line, have degree exactly $d$, satisfy $q(1)=1$, and attain the corresponding bound for nonnegative curvature. The polynomials \[ q_\varepsilon(t)=\frac{q(t)+\varepsilon}{1+\varepsilon}, \qquad\varepsilon>0, \] are strictly positive on $[0,\infty)$, retain their degree and the normalization at one, and approach the same bound. Integrating $P_\varepsilon''(x)=-\kappa q_\varepsilon(x/a)$ with $P_\varepsilon(a)=c$, $P_\varepsilon'(a)=0$ gives admissible profiles for all sufficiently small $\varepsilon$ whenever $c>\kappa a^2b_D$. This proves the claimed sharpness as an infimum in the strictly concave class. \end{proof} The budget in \eqref{eq:normalized-curvature-budget} explains why degree estimates retaining independently prescribed ecological values must depend on dimensionless ratios such as $(c-\alpha_0)/(\kappa a^2)$. Additional derivatives in a prescribed jet can impose further restrictions. \subsection{The size of an arbitrary prescribed jet} \begin{proposition}\label{prop:jet-size-degree} Let $P$ have degree $D\ge2$, with $P''<0$ on $[0,\infty)$, $P(a)=c>0$, $P'(a)=0$, and $P(0)>0$. Put $n=D-2$. Then, for every $j\ge2$, \begin{equation}\label{eq:jet-size-degree} |P^{(j)}(a)|\le \frac{4^{j-1}c}{a^j}n^{2j-4}(n+1)^2. \end{equation} For $j=2$, the factor $n^0$ is interpreted as one. In particular, \[ -P''(a)\le\frac{4c(n+1)^2}{a^2},\qquad |P'''(a)|\le\frac{16c n^2(n+1)^2}{a^3}. \] \end{proposition} \begin{proof} For $w=-P''>0$, \[ \int_0^a xw(x)\dd x=c-P(0)0$ with $c>\kappa a^2/2$, a smooth compactly supported function $\psi$ with $\psi(0)=0$, $\psi'(0)\ne0$, and $\epsilon>0$ small enough that $\kappa+\epsilon\psi>0$. For large $L$, set \[ w_L(x)=\kappa+\epsilon\psi(L(x-a)),\qquad A_L(x)=c-\int_a^x(x-z)w_L(z)\dd z. \] These profiles are globally strictly concave, have their maximum $c$ at $a$, and have the same curvature $A_L''(a)=-\kappa$. Since \[ \int_0^a xw_L(x)\dd x=\kappa a^2/2+O(L^{-1}), \] their values at zero are positive for all sufficiently large $L$. Their curvature is constant outside a bounded interval, so each has a unique positive zero to the right of $a$. However, $|A_L'''(a)|=\epsilon L|\psi'(0)|\to\infty$. Any polynomial with the same three-jet must therefore have degree tending to infinity by \cref{prop:jet-size-degree}, irrespective of any requested cycles. \end{proof} \subsection{The degree cost of approximation near a fixed focus} \begin{proposition}\label{prop:degree-accuracy-focus} Let $A$ be a fixed admissible polynomial of degree $D_0$, with maximum at $a$, and suppose that the first nonzero odd term of $A(x(s))-A(a)$ has degree $d\ge3$. Let $I$ be a compact prey interval whose interior contains $a$. There are a neighborhood $U$ of the coexistence equilibrium and a constant $C>0$ such that every polynomial $P$ of degree $D$ satisfying \[ j_a^{d-1}P=j_a^{d-1}A,\qquad \norm{P-A}_{C^2(I)}\le\varepsilon \] and having a periodic orbit contained in $U$ obeys \begin{equation}\label{eq:degree-accuracy-focus} \varepsilon\max\{D,D_0,1\}^{d-2}\ge C. \end{equation} Consequently, for sufficiently small $\varepsilon$, $D\ge C_1\varepsilon^{-1/(d-2)}$. For an admissible quadratic $A(x)=c-b(x-a)^2$, this gives $D\ge C_1/\varepsilon$ when its two-jet is retained and a cycle is required in the fixed neighborhood. \end{proposition} \begin{proof} Reduce $I$, if necessary, to $[a-\ell,a+\ell]\Subset(0,\infty)$. Put $r=d-2$ and $Q=P''-A''$, whose degree is at most $n=\max\{D,D_0,2\}-2$. The interior Bernstein inequality \cite[eq.~(1)]{Sroka2022} reads \[ |Q'(x)|\le \frac{n\norm{Q}_{C^0(I)}}{\sqrt{\ell^2-(x-a)^2}}. \] Apply it successively on the intervals $I_j=[a-\ell_j,a+\ell_j]$, where \[ \ell_j=\ell\left(1-\frac j{2r}\right),\qquad 0\le j\le r. \] Since $\ell_{j-1}^2-\ell_j^2\ge\ell^2/(2r)$, each step costs at most $\sqrt{2r}\,n/\ell$. Thus \[ \norm{Q^{(j)}}_{C^0(I_r)} \le\left(\frac{\sqrt{2r}\,n}{\ell}\right)^j \norm{Q}_{C^0(I)},\qquad 0\le j\le r. \] The case $n=0$ is immediate, since all derivatives of $Q$ vanish. Combining these estimates with the assumed bounds for orders zero and one gives \begin{equation}\label{eq:interior-bernstein-iteration} \norm{P-A}_{C^d(I_r)} \le B_r\varepsilon\max\{D,D_0,1\}^{d-2},\qquad B_r=\max\{1,\sqrt{2r}/\ell\}^{r}. \end{equation} Choose a fixed $s$-interval whose image under $x=x(s)$ lies in $I_r$. Composition with this analytic coordinate multiplies the $C^d$ bound by a finite constant depending only on the coordinate and the interval; it also preserves equality of the lower jets. The no-cycle neighborhood in \cref{prop:sharp-topology} therefore applies whenever the right-hand side of \eqref{eq:interior-bernstein-iteration} is sufficiently small. Its contrapositive proves \eqref{eq:degree-accuracy-focus}. For the quadratic profile, \[ A(x(s))-c=-\frac{ba^2}{m}s^2 -\frac{2ba^2}{3m^{3/2}}s^3+O(s^4), \] so $d=3$, giving the final assertion. \end{proof} The polynomial-baseline hypothesis is essential to this proof: $P''-A''$ must be a polynomial before Bernstein's inequality can be applied. For an arbitrary smooth baseline, quantitative approximation also requires control of its derivatives or its modulus of continuity. \begin{proposition}\label{prop:oscillatory-baseline-degree} Even with the entire germ at $a$, the values at $0$ and $a$, and the positive zero fixed, there is a family of admissible smooth strictly concave profiles for which polynomial $C^2$ approximation to a fixed positive tolerance requires unbounded degree. \end{proposition} \begin{proof} Start with $A(x)=c-\kappa(x-a)^2/2$, where $c>\kappa a^2/2$, and write $b=a+\sqrt{2c/\kappa}$. Choose disjoint compact intervals in $(a,b)$, a smooth cutoff $0\le\chi\le1$ supported in the first and equal to one on a nondegenerate interval $J$, and a nonzero nonnegative smooth function $\psi$ supported in the second. Set \[ \gamma_L= \frac{\int_a^b(b-x)\chi(x)\sin(Lx)\dd x} {\int_a^b(b-x)\psi(x)\dd x},\qquad w_L=\kappa+\frac\kappa4\bigl[\chi(x)\sin(Lx)-\gamma_L\psi(x)\bigr]. \] Integration by parts gives $\gamma_L=O(L^{-1})$. Thus $w_L>0$ for large $L$. Define $A_L(x)=c-\int_a^x(x-z)w_L(z)\dd z$. The supports lie away from $a$, so every $A_L$ has the same germ at $a$ and the same value at zero as $A$. The definition of $\gamma_L$ gives $A_L(b)=0$. Strict concavity makes $b$ its unique positive zero. On $J$, $w_L=\kappa+\kappa\sin(Lx)/4$. If $\norm{P-A_L}_{C^2([0,b])}<\kappa/8$, then $-P''-\kappa$ has alternating signs at the successive extrema of $\sin(Lx)$ in $J$. There are $L|J|/\pi+O(1)$ such extrema. The intermediate value theorem forces at least that number minus one of distinct zeros of $-P''-\kappa$. Hence $\deg P\ge L|J|/\pi-O(1)$, proving the assertion. \end{proof} \subsection{A quantitative Hermite--Bernstein construction} We now give explicit sufficient degree bounds for a specified smooth target and a specified approximation tolerance. Both the finite jet and the two endpoint values are imposed exactly, while strict concavity holds on the entire half-line. \begin{theorem}\label{thm:quantitative-polynomial-approximation} Let $\alpha$ be admissible and globally strictly concave, with $\alpha(a)=c$, $\alpha'(a)=0$, and positive zero $x_1>a$. Fix $R>x_1$, integers $K,L\ge2$, and $\varepsilon>0$. The constants defined below give an integer $n\ge1$ and a polynomial $\widehat\alpha$ such that \begin{gather} \widehat\alpha''(x)<0\quad(x\ge0),\qquad j_a^K\widehat\alpha=j_a^K\alpha, \qquad \widehat\alpha(0)=\alpha(0),\qquad \widehat\alpha(x_1)=0,\label{eq:quantitative-constraints}\\ \norm{\widehat\alpha-\alpha}_{C^L([0,R])}\le\varepsilon, \qquad \deg\widehat\alpha\le2n+4\max\{K,L\}. \label{eq:hermite-bernstein-degree} \end{gather} All constants involve only finite derivatives, integrals, and norms of the prescribed target on $[0,R]$. \end{theorem} Here are the constants and the construction. Put \[ w=-\alpha'',\qquad 0<\zeta<\min_{[0,R]}w,\qquad q=\sqrt{w-\zeta},\qquad H=\max\{K,L\}-2. \] Choose the smallest even integer $J>H$, and define $g_0=(x-a)^J$, $g_1=(x-a)^{J+1}$. For a continuous function $p$ and $t=(t_0,t_1)\in\R^2$, set \[ p_t=p(1+t_0g_0+t_1g_1),\qquad w_{p,t}=p_t^2+\zeta, \] \[ \Phi(p,t)= \begin{pmatrix} \displaystyle\int_0^a xw_{p,t}(x)\dd x-[c-\alpha(0)]\\[5pt] \displaystyle\int_a^{x_1}(x_1-x)w_{p,t}(x)\dd x-c \end{pmatrix}. \] The matrix $J_0=D_t\Phi(q,0)$ is invertible: its first row has signs $(+,-)$ and its second has signs $(+,+)$. Let \begin{align*} Q&=\norm{q}_{C^0([0,R])},& P_*&=Q+1,& B&=\sup_{[0,R]}(|g_0|+|g_1|),\\ W&=\tfrac12\max\{a^2,(x_1-a)^2\},& \Lambda&=\norm{J_0^{-1}}_{\infty\to\infty},& A&=W(2Q+1), \end{align*} \begin{equation}\label{eq:bernstein-correction-constants} \begin{gathered} F=2AB,\qquad G=2WP_*^2B^2,\qquad \tau=\min\{1,(4\Lambda G)^{-1}\},\\ \eta_0=\min\{1,(4\Lambda F)^{-1},\tau/(2\Lambda A)\}, \qquad C_t=2\Lambda A. \end{gathered} \end{equation} For the differentiable estimate, define \begin{align*} A_*&=\norm{q^{(H+1)}}_{C^0([0,R])} \max_{0\le j\le H}\frac{R^{H-j+1}}{(H-j)!},\\ Q_H&=\norm{q}_{C^H([0,R])},\\ B_H&=\max_{0\le j\le H}\sum_{\ell=0}^j\binom j\ell \left(\norm{g_0^{(j-\ell)}}_{C^0([0,R])} +\norm{g_1^{(j-\ell)}}_{C^0([0,R])}\right),\\ C_{q,H}&=1+(Q_H+1)C_tB_H,\\ C_H&=\max\{1,R,R^2/2\}\,2^HC_{q,H}(2Q_H+C_{q,H}). \end{align*} Finally choose \begin{equation}\label{eq:bernstein-n-choice} \eta=\min\{\eta_0,\varepsilon/C_H\},\qquad n=\max\left\{1,\left\lceil(A_*/\eta)^2\right\rceil\right\}. \end{equation} The extra lower bound $n\ge1$ includes the case $A_*=0$. \begin{proof}[Proof of \cref{thm:quantitative-polynomial-approximation}] Let $B_ng$ be the degree-$n$ Bernstein polynomial of $g=q^{(H)}$ on $[0,R]$. If $Z$ has the binomial distribution with parameters $n$ and $x/R$, then $B_ng(x)=\mathbb E[g(RZ/n)]$. The Lipschitz bound and Cauchy--Schwarz give \[ |B_ng(x)-g(x)| \le\norm{q^{(H+1)}}_{C^0}\sqrt{\frac{x(R-x)}n} \le\frac{R}{2\sqrt n}\norm{q^{(H+1)}}_{C^0}. \] Subtract the error at $a$, and integrate the resulting polynomial $H$ times, using the exact derivatives of $q$ at $a$ as the integration constants. The resulting polynomial $p_n$ has degree at most $n+H$, retains $j_a^Hq$, and satisfies \begin{equation}\label{eq:bernstein-hermite-errors} \norm{p_n^{(j)}-q^{(j)}}_{C^0([0,R])} \le\frac{R^{H-j+1}}{(H-j)!\sqrt n} \norm{q^{(H+1)}}_{C^0([0,R])},\qquad 0\le j\le H. \end{equation} For $H=0$, this just means the Bernstein polynomial corrected by a constant. Thus \eqref{eq:bernstein-n-choice} ensures $\norm{p_n-q}_{C^H}\le\eta$. We next solve the two moment equations with quantitative control. Write $p=p_n$. Since $\eta\le1$, one has $\norm{p}_{C^0}\le P_*$. The moment residual satisfies \[ \norm{\Phi(p,0)}_\infty\le W(2Q+1)\eta=A\eta. \] For $\norm{t}_\infty\le\tau$, differentiating the moments gives \[ \norm{D_t\Phi(p,t)-J_0}_{\infty\to\infty} \le 2W(2Q+1)B\eta+2WP_*^2B^2\tau =F\eta+G\tau. \] Consequently, the map $\mathcal T(t)=t-J_0^{-1}\Phi(p,t)$ has Lipschitz constant at most $1/2$ on the closed $\tau$-ball. Also $\norm{\mathcal T(0)}_\infty\le\Lambda A\eta\le\tau/2$, so the ball is invariant. The contraction theorem gives a fixed point with \begin{equation}\label{eq:quantitative-moment-solution} \Phi(p,t)=0,\qquad \norm{t}_\infty\le2\Lambda A\eta=C_t\eta. \end{equation} The product rule, the definition of $B_H$, and \eqref{eq:quantitative-moment-solution} yield \[ \norm{p_t-q}_{C^H}\le \eta+(Q_H+1)C_tB_H\eta=C_{q,H}\eta. \] Using $\norm{uv}_{C^H}\le2^H\norm{u}_{C^H}\norm{v}_{C^H}$, we obtain \[ \norm{p_t^2-q^2}_{C^H} \le2^HC_{q,H}(2Q_H+C_{q,H})\eta. \] Define the growth polynomial by \[ \widehat\alpha(x)=c-\int_a^x(x-z)\bigl[p_t(z)^2+\zeta\bigr]\dd z. \] Its second derivative is negative on the whole half-line. Integrating the preceding error estimate twice, with matching value and slope at $a$, gives \[ \norm{\widehat\alpha-\alpha}_{C^{H+2}([0,R])} \le C_H\eta\le\varepsilon. \] Because $J>H$, multiplication by $1+t_0g_0+t_1g_1$ retains the jet of $p_n$ through order $H$. Squaring then retains the curvature jet, and integration retains the growth jet through order $H+2\ge K$. The two zero moments give the values at $0$ and $x_1$ exactly. This proves \eqref{eq:quantitative-constraints} and the approximation assertion. Finally, \[ \deg p_t\le n+H+J+1,\qquad \deg\widehat\alpha\le2n+2H+2J+4. \] Since the smallest even integer $J>H$ satisfies $J\le H+2$, the last expression is at most $2n+4H+8=2n+4\max\{K,L\}$, proving \eqref{eq:hermite-bernstein-degree}. \end{proof} For a specified target with a finite family of hyperbolic cycles, one may take $L=2$ and choose $\varepsilon$ below a persistence tolerance for those cycles. An exact-count argument additionally specifies the required displacement margins outside the cycle neighborhoods and the finite-jet condition near the equilibrium; the higher-$L$ form of \cref{thm:quantitative-polynomial-approximation} supplies the needed approximation without losing the endpoint or concavity constraints. The degree depends explicitly on these tolerances and on the derivatives of the target. \Cref{cor:no-uniform-arbitrary-jet-degree,prop:degree-accuracy-focus,prop:oscillatory-baseline-degree} show why such dependence cannot generally be removed. % --- infinite_cycles.tex --- \section{Infinitely many cycles in one smooth profile} \label{sec:infinite} The finite realization theorem can be iterated. The central issue is persistence: later perturbations supported at smaller values of $|s|$ still affect the outer cycles, because those cycles also cross the strip where $s$ is small. We choose the later perturbations so that their combined $C^1$ norm remains within the persistence neighborhood of every previously selected cycle. \begin{proof}[Proof of \cref{thm:main-infinite}] Choose a globally strictly concave extension $\alpha_b$ of \eqref{eq:center-profile} by \cref{lem:concave-extension}, with $\alpha_b(0)>0$. Its positive zero is denoted by $x_1$. Restrict the coordinate patch to a compact interval $I\Subset(0,x_1)$ containing $a$, on which $\alpha_b=\alpha_{\mathrm{cen}}$. The normalization at the end of \cref{sec:normal-form} permits us to work with $c=1$. For general $c$, a perturbation $R(u)$ in the normalized coordinate $u=s/\sqrt c$ corresponds to the growth perturbation $cR(s(x)/\sqrt c)$; flatness, smallness, and the following persistence argument are unchanged. In the Duffing coordinates of \cref{sec:moments}, write \[ s'=-Y+R(s),\qquad Y'=s-\frac{s^3}{2}+sR(s),\qquad K(s,Y)=\frac{Y^2}{2}+U(s),\quad U(s)=\frac{s^2}{2}-\frac{s^4}{8}. \] Choose $r_1>0$ so small that $|s|\le r_1$ lies inside the coordinate patch and $U$ is strictly increasing on $[0,r_1]$. \emph{One insertion.} Suppose that the current perturbation $R_*$ vanishes on $|s|\le r$. Choose a regular oval $K=h_*$ whose maximal $|s|$ is less than $r/4$, and a compact interval $E$ about $h_*$ whose ovals all have maximal $|s|0$ is defined by $U(A(h))=h$, then for odd $G$, \begin{equation}\label{eq:cutoff-moment} \oint_{K=h}G(s)\,dY =4\int_0^{A(h)} \frac{G(s)U'(s)}{\sqrt{2(h-U(s))}}\,ds. \end{equation} The cutoff changes this integral only on $[0,2\delta]$. For $h\in E$ the denominator stays away from zero there, whereas $P(s)U'(s)=O(s^4)$. The integral difference, and its derivative in $h$, are therefore $O(\delta^5)$. Thus, for small fixed $\delta$, the cutoff moment still has a simple zero inside $E$. For every sufficiently small nonzero $\varepsilon$, the perturbation $R_*+\varepsilon G_\delta$ consequently creates a hyperbolic cycle near that oval. This follows from the $C^1$ return-map expansion \[ \Delta(h,\varepsilon) =\varepsilon\oint_{K=h}G_\delta(s)\,dY+O(\varepsilon^2). \] The new perturbation remains zero on $|s|\le\delta$, and its amplitude can be made arbitrarily small in any prescribed finite collection of derivative norms. \emph{Diagonal construction.} Let $H_n=\varepsilon_nG_n$ and $R_n=\sum_{j=1}^nH_j$. Choose the insertion bands successively so that \[ \supp H_n\subset\{\delta_n\le|s|<3r_n/4\},\qquad r_{n+1}<\min(\delta_n/2,r_n/4). \] The next insertion is then always made in an exact center region. Choose the associated energy intervals, and compact annular neighborhoods $V_n$ of the inserted cycles, pairwise disjoint and converging to the origin. After stage $n$, choose $\eta_n>0$ such that any change of $R_n$ of $C^1([-r_1,r_1])$ norm less than $\eta_n$ preserves a hyperbolic cycle in $V_n$, with the same stability. This is possible by hyperbolic periodic-orbit persistence. At every later insertion choose the nonzero amplitude so small that \begin{equation}\label{eq:diagonal-budget} \|H_n\|_{C^{n+2}([-r_1,r_1])} \le 2^{-n-2}\min\{1,\kappa,\eta_1,\ldots,\eta_{n-1}\}, \end{equation} where the absent persistence radii at the first stage are omitted. Cycle creation places no lower positive bound on the amplitude, so these constraints are compatible with every insertion. For each fixed derivative order, the tail of \eqref{eq:diagonal-budget} is summable. Hence $R=\sum_{n\ge1}H_n$ converges in $C^\infty$. Each summand vanishes near zero, and therefore $R^{(j)}(0)=0$ for every $j$. For every fixed $n$, \[ \|R-R_n\|_{C^1} \le\sum_{j>n}\|H_j\|_{C^1}<\eta_n. \] Thus the limiting system retains a hyperbolic cycle in every $V_n$. They are distinct, surround the same equilibrium, and converge to it. \emph{Preservation of concavity.} Composition with the fixed energy coordinate is bounded in $C^2$ on $I$. Choose $\kappa$ sufficiently small that the total growth-rate perturbation has second derivative of absolute value less than $\tfrac12\min_I(-\alpha_b'')$. Add this perturbation to $\alpha_b$ on the patch, extending it by zero outside. It is supported strictly inside the patch. The resulting global $\alpha$ is smooth and strictly concave, and has exactly the full jet of $\alpha_{\mathrm{cen}}$ at $a$. Its value at zero and its positive zero $x_1$ are unchanged. The coordinate changes and positive changes of time preserve the hyperbolicity and nesting of the cycles. \end{proof} \begin{corollary}[Nonanalyticity and limiting periods] The profile in \cref{thm:main-infinite} is nonanalytic at $a$. The original-time periods $T_n$ and nontrivial Floquet multipliers $\mu_n$ of its shrinking cycles satisfy \[ T_n\longrightarrow\frac{2\pi}{\sqrt{mc}},\qquad \mu_n\longrightarrow1. \] \end{corollary} \begin{proof} An analytic profile with the specified full jet would equal $\alpha_{\mathrm{cen}}$ near $a$ and have a center neighborhood, contradicting the presence of isolated cycles approaching the equilibrium. The linearization in the original time has eigenvalues $\pm i\sqrt{mc}$. On shrinking circles the angular flow converges uniformly to that of the linear oscillator, which gives the period limit. The multiplier limit follows from the divergence formula, since the divergence tends to its equilibrium value zero and the periods stay bounded. \end{proof} Individual cycles are hyperbolic, but the family is not uniformly hyperbolic near its accumulation point. The polynomial approximation of \cref{sec:global} preserves any finite selected subfamily; it does not retain this entire infinite family in one polynomial model. % --- fixed_response.tex --- \section{Fixed nonlinear responses and an elasticity threshold} \label{sec:responses} The localization argument allows the functional response to be fixed in a considerably larger class. Consider \begin{equation}\label{eq:general-response} \dot x=x\alpha(x)-p(x)y=p(x)(h(x)-y),\qquad \dot y=yg(x),\qquad h(x)=\frac{x\alpha(x)}{p(x)}. \end{equation} Here $p,g$ are fixed smooth functions on the positive prey axis, $p(x)>0$ there, and $a>0$ is a simple predator break-even point: $g(a)=0$, $g'(a)>0$. Smooth extension to the boundary may be imposed when required for the ecological model; all the cycles below stay in the interior. The usual conversion relation is $g=kp-m$. The classical generalized Filippov transformation \cite[\S3]{GasullGuillamon1995} is \begin{equation}\label{eq:general-coordinate} s(x)=\sgn(x-a)\sqrt{2\int_a^x\frac{g(u)}{p(u)}\,du}, \qquad v=\log(y/c),\qquad c=h(a)>0. \end{equation} Near $a$, $s$ is smooth and increasing, with $s'(a)=\sqrt{g'(a)/p(a)}$. Since $ss'=g/p$, the time factor $d\tau/dt=p(x)s'(x)=g(x)/s(x)$ is positive and gives \[ s'=h(x(s))-ce^v,\qquad v'=s. \] Thus the previous canonical constructions prescribe the prey isocline $h$. To control the actual per-capita rate $\alpha=(p/x)h$, an additional condition is needed. Define the response elasticity by \[ \mathcal E_p(x)=\frac{xp'(x)}{p(x)}. \] \begin{theorem}[Strictly concave realization below elasticity two] \label{thm:elasticity} If $\mathcal E_p(a)<2$, then for every $N\ge1$ there is a smooth globally strictly concave per-capita growth rate $\alpha$, positive at zero, with a unique positive maximum and a positive carrying capacity, for which \eqref{eq:general-response} has at least $N$ hyperbolic positive limit cycles near the coexistence equilibrium at $x=a$. Conversely, a globally strictly concave $\alpha$ with $\alpha(0)>0$ can have a nonhyperbolic coexistence equilibrium at this simple break-even point only if $\mathcal E_p(a)<2$. \end{theorem} \begin{proof} Put $r(x)=p(x)/x$. For sufficiency, first prescribe locally \[ \alpha_{\mathrm{loc}}(x) =r(x)\left[c-\frac12s(x)^2\right]. \] At $a$ this function satisfies \begin{align*} \alpha_{\mathrm{loc}}''(a) &=c r''(a)-r(a)\frac{g'(a)}{p(a)} =c r''(a)-\frac{g'(a)}a,\\ \alpha_{\mathrm{loc}}(a)-a\alpha_{\mathrm{loc}}'(a) &=c\left[\frac{2p(a)}a-p'(a)\right]>0. \end{align*} Choose $c>0$ sufficiently small that the second derivative is negative, and restrict to an interval on which strict concavity holds. The canonical isocline is the explicit center $h=c-s^2/2$. The normalization $s=\sqrt c\,u$ reduces it to \cref{sec:moments}. Choose $N$ sufficiently small target energies so that all corresponding ovals project strictly inside this interval. The moment interpolation theorem gives $N$ simple Melnikov zeros. A sufficiently small nonzero amplitude creates $N$ hyperbolic cycles, while preserving the two strict inequalities above and local strict concavity. For the extension, choose in advance an interval so short that the weighted curvature integral on the left in the proof of \cref{lem:concave-extension} is less than half the displayed positive quantity. Then choose the ovals inside a smaller interval and the perturbation sufficiently small. The same extension construction retains the entire prey projection of these cycles. It produces a global strictly concave function with $\alpha(0)>0$, a positive maximum, and a positive right zero, without changing those cycles. For necessity, the Jacobian at $(a,h(a))$ has trace $p(a)h'(a)$ and determinant $p(a)h(a)g'(a)>0$. Nonhyperbolicity therefore implies $h'(a)=0$. Strict concavity gives \[ x\alpha'(x)<\alpha(x)-\alpha(0)\qquad(x>0). \] On the other hand, \begin{equation}\label{eq:isocline-derivative} h'(x)=\frac{x\alpha'(x)+(1-\mathcal E_p(x))\alpha(x)}{p(x)}. \end{equation} At $a$, the condition $h'(a)=0$ consequently yields \[ (2-\mathcal E_p(a))\alpha(a)>\alpha(0)>0, \] which proves $\mathcal E_p(a)<2$. \end{proof} The maximum of $\alpha$ need not occur at $a$ for a nonlinear response. In the local center construction, $\alpha'(a)=c(p(x)/x)'|_{x=a}$. Thus the fixed-response statement distinguishes the maximum of per-capita growth from the maximum of the prey isocline. \begin{proposition}[A global exclusion condition] \label{prop:elasticity-exclusion} Suppose $\alpha$ is concave, $\alpha(0)>0$, and has a positive zero $x_1$. If $\mathcal E_p(x)\ge2$ for every $00$ imply $\alpha>0$ on $(0,x_1)$ and $\alpha<0$ to the right of $x_1$. At the maximum prey coordinate of a positive periodic orbit, $x\alpha(x)=p(x)y>0$; thus every such orbit is contained in $00$. Each row admits the realization theorem at any positive simple break-even point. In particular, varying the response function is unnecessary for these local many-cycle constructions: it can be fixed in advance. % --- computations.tex --- \section{Algebraic verification and numerical illustrations} \label{sec:computations} The constructions above are analytic existence and counting proofs. The separate computational supplement contains exact symbolic checks of selected instances of the degree construction and numerical illustrations of the Abelian-integral interpolation. \subsection{Exact algebra for the degree construction} The script \texttt{quantitative\_computations/verify\_quantitative\_algebra.py} constructs $A_D$, $Q_j$, and $b_j$ with exact symbolic arithmetic for the six pairs $(N,K)$ in \cref{tab:degree-examples}. It verifies the degree and jet identities, the lowest nonzero order of every $Q_j$, the harmonic-moment identity \eqref{eq:polynomial-target-M}, the curvature formula \eqref{eq:taylor-curvature}, and the signed Taylor remainder. The same-degree correction is checked against its two endpoint equations and its vanishing order at $a$. The degree-six example is recorded separately with rational coefficients. The outputs contain every polynomial and coefficient, a degree table, and the individual verification results. These checks verify algebraic identities for the six selected pairs $(N,K)$. The existence of a sufficiently small parameter interval, and the exclusion of all additional periodic orbits on that interval, are supplied by the analytic proof of \cref{thm:exact-degree}. \subsection{Moment interpolation for three and four cycles} The following computations illustrate the regular-annulus moment interpolation and locate selected cycles by direct integration. They use floating-point arithmetic and do not supply a certified global cycle count at the displayed parameter value. Set $k=m=\beta=c=1$ and use the local profile \begin{equation}\label{eq:numerical-profile} \alpha(x)=1-\frac{s(x)^2}{2}+\delta s(x)^3P_N(s(x)^2), \qquad s(x)=\sgn(x-1)\sqrt{2(x-1-\log x)}. \end{equation} Write $P_N(z)=\sum_{j=0}^Nb_jz^j$ and normalize $b_N=-1$. For $N=3$, prescribe Melnikov zeros at $0.05,0.10,0.15$; for $N=4$, prescribe them at $0.04,0.08,0.12,0.16$. Solving $\sum b_jJ_j(h_i)=0$ with 80 decimal digits gives approximately \begin{align*} P_3(z)={}&0.00462584243209378 -0.0975516389901764z+0.575793634905719z^2-z^3,\\ P_4(z)={}&-0.000734820723917012 +0.0220639437260988z-0.209685998923226z^2\\ &\quad+0.783964939740499z^3-z^4. \end{align*} The complete coefficients are included in the accompanying CSV and JSON files. Independent quadrature of the defining moment integral agrees with \eqref{eq:hypergeom} with a computed relative discrepancy below $3.2\times10^{-81}$ at the check energy $h=0.117$. \Cref{fig:melnikov} shows the prescribed zeros. \begin{figure}[htbp] \centering \includegraphics[width=\linewidth]{figures/melnikov_prescribed_zeros.pdf} \caption{The interpolating Melnikov functions for three and four selected cycles. Each is divided by $h^2$ and then by its maximum absolute value on the displayed interval. These positive normalizations preserve zeros and their crossing signs. Filled blue markers indicate the attracting branches for positive $\delta$; open orange markers indicate repelling branches.} \label{fig:melnikov} \end{figure} \subsection{Direct return maps and multipliers} To integrate the system efficiently, put $q=y-1$ and $p=s\sqrt y$. Here $p$ denotes an oscillator coordinate. With $d\sigma/d\tau=\sqrt{1+q}$, the exact equations are \begin{equation}\label{eq:numerical-oscillator} \frac{dq}{d\sigma}=p,\qquad \frac{dp}{d\sigma}=-q+ \delta\frac{p^3}{(1+q)^{3/2}} P_N\!\left(\frac{p^2}{1+q}\right). \end{equation} Starting at $q=0$, $p=\sqrt{2h}>0$, use the next positive crossing of $q=0$ to define the return map in the coordinate $h=p^2/2$. Its fixed points are found by bracketed scalar root searches. The transverse multiplier is computed by integrating the divergence of \eqref{eq:numerical-oscillator} over the periodic orbit. An independent calculation in the original variables uses \[ \log\mu=\int_\Gamma x\alpha'(x)\,dt. \] This formula follows by integrating the full divergence and using $\int(\alpha(x)-y)dt=\int(x-1)dt=0$ on a periodic orbit. For $\delta=0.1$, the results are in \cref{tab:cycles}. These examples keep the interpolating coefficients fixed as $\delta$ changes; their actual section positions therefore differ slightly from the prescribed first-order energies. The exact placement construction of \cref{sec:moments} would additionally adjust those coefficients. \begin{table}[htbp] \centering \small \begin{tabular}{ccrll} \toprule $N$ & Target $h_i$ & Actual section energy & Multiplier $\mu$ & Stability\\ \midrule 3 & 0.05 & 0.0499999109 & 0.999918785905 & Attracting\\ 3 & 0.10 & 0.1000002816 & 1.000183252618 & Repelling\\ 3 & 0.15 & 0.1499988200 & 0.999059900849 & Attracting\\ \midrule 4 & 0.04 & 0.0400000066 & 1.000007729473 & Repelling\\ 4 & 0.08 & 0.0799999857 & 0.999988415660 & Attracting\\ 4 & 0.12 & 0.1200000385 & 1.000029541335 & Repelling\\ 4 & 0.16 & 0.1599998569 & 0.999819716489 & Attracting\\ \bottomrule \end{tabular} \caption{Numerically located fixed points of the exact return map for $\delta=0.1$. Energies are rounded to ten decimal places. The computed multipliers indicate alternating stability. The calculations use floating-point arithmetic without interval bounds.} \label{tab:cycles} \end{table} The DOP853 integrations were repeated with two controls: $(\mathrm{rtol},\mathrm{atol},\mathrm{max\_step}) =(2\cdot10^{-12},2\cdot10^{-14},0.075)$ and $(5\cdot10^{-14},2\cdot10^{-15},0.02)$. The largest change in a section energy was below $10^{-11}$; the largest multiplier change was below $10^{-13}$. The two divergence calculations of $\log\mu$ agreed to within $2\times10^{-15}$. To recover the original prey coordinate, the code uses \[ x=-W_0(-e^{-1-s^2/2})\quad(s<0),\qquad x=-W_{-1}(-e^{-1-s^2/2})\quad(s>0), \] with a local inverse series near the branch point. The original predator coordinate is $y=1+q$. The three selected cycles in these coordinates are shown in \cref{fig:three-cycles}. \begin{figure}[htbp] \centering \includegraphics[width=.78\linewidth]{figures/three_cycles_original_coordinates.pdf} \caption{Three numerically located limit cycles in the original positive quadrant for the local profile \eqref{eq:numerical-profile}. Solid blue curves are attracting, the dashed orange curve is repelling, and arrows show the flow direction. All surround $E=(1,1)$.} \label{fig:three-cycles} \end{figure} \subsection{Shape checks and reproduction} The analytic expression for $\alpha''$ was evaluated on intervals slightly larger than the cycles' prey projections. A grid of 50,001 points, together with the peak point and local extremum refinement, gave the estimates \[ \begin{array}{c|c|c} N&\text{checked prey interval}&\text{estimated maximum of }\alpha''\\\hline 3 &[0.5086701260,\,1.7083293776]&-0.3636065093\\ 4 &[0.4945100775,\,1.7384221971]&-0.3383078145. \end{array} \] These are numerical checks of local strict concavity. The global smooth and polynomial shape conclusions rest on the extension and approximation proofs, not on this grid. In particular, \eqref{eq:numerical-profile} contains $\log x$ and is not itself claimed to be the final growth law at $x=0$ or for arbitrarily large $x$. The local profiles and their second derivatives appear in \cref{fig:concavity}. \begin{figure}[htbp] \centering \includegraphics[width=\linewidth]{figures/local_growth_and_concavity.pdf} \caption{The local growth profiles in the three-cycle and four-cycle examples, together with their analytically evaluated second derivatives on the checked intervals. The two growth curves are visually close. This figure illustrates strict concavity locally; the global realization is supplied by the analytic construction.} \label{fig:concavity} \end{figure} The separate computational supplement contains \texttt{computations/reproduce\_cycles.py}, instructions, all coefficients and cycle samples, and editable figure sources. Running the script with Python~3.10 or later, NumPy, SciPy, mpmath, Matplotlib, and a LaTeX installation recreates the numerical outputs and figures. The supplied results record the software versions, solver settings, and all sensitivity comparisons. No numerical assertion is used to prove the existence, density, cyclicity, or infinite-cycle theorems. \clearpage \section{Further questions} \label{sec:questions} The construction gives the sufficient degree $D=2N+2\lceil K/2\rceil+2$ for its chosen center jet. The jet-only lower bound is $2\lceil K/2\rceil$, and the endpoint budget and the sizes of higher derivatives give additional independent obstructions. Determining the smallest degree compatible with all these data, or a sharp global cycle bound at a fixed polynomial degree, requires further information beyond the three-scale construction. A separate problem is to impose independently specified low-density growth, carrying capacity, and focus jets while requiring a global exact count. The parity law is necessary when the focus stability is fixed. The arbitrary-data approximation theorem realizes any finite lower bound on the number of cycles, whereas the exact-count theorems use center jets and compatible global reference profiles. A classification of exact realization for arbitrary prescribed data would connect these two forms of flexibility. The explicit coefficients and concavity inequalities also permit quantitative work on the small-parameter range. The present proof gives $\rho_0>0$ by compactness and uniform return-map estimates. Certified bounds on $\rho_0$, return-map remainders, and multiplier separation would turn the symbolic families into finite numerical certificates at selected parameter values. The high-order focus makes these estimates more demanding as $N$ or $K$ increases. \phantomsection \addcontentsline{toc}{section}{References} \begin{thebibliography}{8} \providecommand{\natexlab}[1]{#1} \providecommand{\url}[1]{\texttt{#1}} \expandafter\ifx\csname urlstyle\endcsname\relax \providecommand{\doi}[1]{doi: #1}\else \providecommand{\doi}{doi: \begingroup \urlstyle{rm}\Url}\fi \bibitem[Gasull and Guillamon(1995)]{GasullGuillamon1995} Armengol Gasull and Antoni Guillamon. \newblock Non-existence, uniqueness of limit cycles and center problem in a system that includes predator--prey systems and generalized {Li\'{e}nard} equations. \newblock \emph{Differential Equations and Dynamical Systems}, 3\penalty0 (4):\penalty0 345--366, 1995. \newblock URL \url{https://web.mat.upc.edu/antoni.guillamon/papers/GasullGuillamon-DiffEqDynSyst1995a.pdf}. \bibitem[Gasull et~al.(1996)Gasull, Guillamon, Li, and Zhang]{GasullGuillamonLiZhang1996} Armengol Gasull, Antoni Guillamon, Chengzhi Li, and Zhifen Zhang. \newblock Study of perturbed {Lotka--Volterra} systems via {Abelian} integrals. \newblock \emph{Journal of Mathematical Analysis and Applications}, 198\penalty0 (3):\penalty0 703--728, 1996. \newblock \doi{10.1006/jmaa.1996.0109}. \newblock URL \url{https://web.mat.upc.edu/antoni.guillamon/papers/GasullGuillamonLiZhang-JMAA1996.pdf}. \bibitem[Gavrilov and Iliev(2003)]{GavrilovIliev2003} Lubomir Gavrilov and Iliya~D. Iliev. \newblock Two-dimensional {Fuchsian} systems and the {Chebyshev} property. \newblock \emph{Journal of Differential Equations}, 191\penalty0 (1):\penalty0 105--120, 2003. \newblock \doi{10.1016/S0022-0396(02)00116-X}. \newblock URL \url{https://www.math.univ-toulouse.fr/~gavrilov/publications/34.pdf}. \bibitem[Kuang(1988)]{Kuang1988} Yang Kuang. \newblock Nonuniqueness of limit cycles of {Gause}-type predator--prey systems. \newblock \emph{Applicable Analysis}, 29\penalty0 (3-4):\penalty0 269--287, 1988. \newblock \doi{10.1080/00036818808839785}. \bibitem[Sroka(2022)]{Sroka2022} Grzegorz Sroka. \newblock Constants in {Markov}'s and {Bernstein} inequality on a finite interval in $\mathbb{R}$. \newblock \emph{Analysis and Mathematical Physics}, 12:\penalty0 128, 2022. \newblock \doi{10.1007/s13324-022-00711-8}. \newblock URL \url{https://link.springer.com/article/10.1007/s13324-022-00711-8}. \bibitem[Svirezhev(1984)]{Svirezhev1984} Yu.~M. Svirezhev. \newblock Modern problems of mathematical ecology. \newblock In \emph{Proceedings of the International Congress of Mathematicians, Warszawa 1983}, volume~2, pages 1677--1693. PWN--North-Holland, 1984. \newblock URL \url{https://www.mathunion.org/fileadmin/ICM/Proceedings/ICM1983.2/ICM1983.2.ocr.pdf}. \newblock The limit-cycle question occurs on p. 1682. Conference held in 1983; proceedings published in 1984. \bibitem[Wrzosek(1990)]{Wrzosek1990} Dariusz~M. Wrzosek. \newblock Limit cycles in predator--prey models. \newblock \emph{Mathematical Biosciences}, 98\penalty0 (1):\penalty0 1--12, 1990. \newblock \doi{10.1016/0025-5564(90)90009-N}. \bibitem[Zegeling and Kooij(2021)]{ZegelingKooij2021} Andr\'{e} Zegeling and Robert~E. Kooij. \newblock Several bifurcation mechanisms for limit cycles in a predator--prey system. \newblock \emph{Qualitative Theory of Dynamical Systems}, 20\penalty0 (3):\penalty0 65, 2021. \newblock \doi{10.1007/s12346-021-00501-w}. \newblock Article 65, 48 pages. \end{thebibliography} \end{document}