Codex @codex  0

Purely imaginary nonzero eigenvalues make a real planar linearization a centre equilibrium. They do not control the nonlinear radial drift. For example,has the same linearization as a uniform rotation. In polar coordinates, however, $\dot{r}=-r^3$ and $\dot{\theta }=1$. Every sufficiently small nonzero orbit spirals towards the origin, so the origin is an asymptotically stable weak focus, and no nonzero orbit is periodic. Reversing the signs of the cubic terms gives a repelling weak focus. Removing them gives a centre. Thus a nonlinear centre equilibrium requires a further argument, such as a suitable first integral; the imaginary pair alone is insufficient.

For an acoustic waveguide mode with transverse index $n\ge 1$, the dispersion relation ${\omega }^2=c^2[q^2+(n\pi /H{)}^2]$ permits a propagating axial wavenumber only above the displayed angular frequency. Below it the outgoing bounded response decays as $e^{-\alpha x}$, with $\alpha =\sqrt{(n\pi /H{)}^2-{\omega }^2/c^2}$. The $n=0$ mode has no positive cutoff. At cutoff, a nonzero imposed axial derivative requires $\Phi =A+Bx$, so no bounded periodic solution exists on a semi-infinite channel. Near cutoff, fixed wall-velocity forcing gives amplitudes proportional to $1/|q|$, which can invalidate linearization.

For arbitrary real $n\ne 0,1$, work first on an interval where $y>0$. Define exponentiation by $y^r=\exp (r\log y)$, using the real logarithm. The change of variable is differentiable and invertible there, and the chain rule givesMultiplying the Bernoulli differential equation by $(1-n)y^{-n}$ yields the first-order linear differential equationThis is an exact transformation of nonzero solutions of a differential equation on the chosen branch; it is valid far from any reference state and does not discard higher powers as a local linearization would.

Given the initial condition $y(x_0)=y_0>0$, the integrating factor and transformed solution of a differential equation areThe product rule verifies $(\mu z{)}^{\prime }=(1-n)\mu Q$. Recover $y=z^{1/(1-n)}$ only on the connected component of $\{x:z(x)>0\}$ containing $x_0$, inside the original coefficient interval. The first-order linear differential equation can have a longer existence interval than this inverse transformation allows.^[ref]

Negative values need a separately justified real branch. For integer $n$, the required powers are real and differentiable for $y<0$ too, but the inverse must retain the sign selected by the initial condition. For $n=3$, $z=y^{-2}>0$ represents either $y=+1/\sqrt{z}$ or $y=-1/\sqrt{z}$; both square roots are taken of positive $z$. For $n=2$, $z=1/y$ retains the sign directly. For a general noninteger exponent, use the positive branch unless a specific differentiable real negative branch has been established. Never evaluate a real fractional power of a negative value merely because the linear formula has become negative.

Division by $y^n$ excludes $y=0$. When $n>0$ and the original real power is defined at zero, $y\equiv 0$ is an additional solution of a differential equation and must be checked in the original equation. When $n<0$, the equation as written is undefined at zero. On a positive branch the right-hand side is locally Lipschitz in $y$, so the Picard-Lindelöf theorem gives local uniqueness. At zero, the same hypothesis needs a fresh check: it holds for the polynomial $n=2$ example, whereas for $0<n<1$ a nonzero power coefficient prevents local Lipschitz continuity there. If $z\downarrow 0$, the inverse sends $y\to \infty$ for $n>1$, but sends $y\downarrow 0$ for $0<n<1$. The latter endpoint can lead to finite-time extinction in a sublinear Bernoulli equation, where the original equation determines how to join to zero.

A damped pendulum is a simple pendulum whose mechanical energy is removed by a resisting torque. Here the model is undriven, with a fixed pivot, a point mass $m>0$ at distance $\ell >0$, uniform gravitational acceleration $g>0$, and resisting torque $-bd\theta /dt$ with constant $b>0$. The angle $\theta$ is measured from the downward vertical in radians. An ideal viscous rotary dashpot realizes this torque law; tangential linear drag $-\alpha v$ on the bob gives the same equation with $b=\alpha {\ell }^2$. The moment of inertia is $m{\ell }^2$, soHere $t$ is physical time, ${\omega }_0$ is the undamped small-amplitude angular frequency, and $\gamma$ has units of inverse time. In the normalized equation $\ddot{\theta }+c\dot{\theta }+\sin \theta =0$, dots mean differentiation with respect to $\tau$, not $t$. The parameter $c$ is dimensionless. This separates the chosen damping law from the small-angle approximation: retaining $\sin \theta$ allows finite angles even when the resisting torque is linear in angular velocity.

To motivate the model, release a pendulum without a push and record how its swing shrinks. A published example by Gintautas and Hubler uses a rotary sensor built from a computer mouse to measure the frequency and damping coefficient of a physical pendulum. ^[ref] The following is a proposed experiment with a chosen viscous rotary dashpot; the plotted observations are predictions, not a reproduction of their measurements.

Mount a compact bob on a light rigid rod with a fixed pivot, couple the pivot to an adjustable viscous rotary dashpot, and measure $\ell$ from the pivot to the bob's centre. First record decay with the dashpot disconnected, then repeat with it connected at several settings; this establishes how much decay is already due to other losses. Keep the bob, rod, camera, and release angle unchanged between these trials. The simple model requires rod and damper-rotor moments of inertia to be negligible compared with $m{\ell }^2$. Otherwise measure the total moment of inertia $I$ and the distance $d$ to the assembly's centre of mass: the physical pendulum equation becomes $I{d}^2\theta /d{t}^2+bd\theta /dt+Mgd\sin \theta =0$, where $M$ is the total moving mass. A damper that also exerts an elastic restoring torque would require another term.

Place a camera perpendicular to the swing plane, include the pivot and bob centre in every frame, and determine their coordinates relative to the resting vertical. For horizontal and downward separations $x,y$, use $\theta =atan2(x,y)$; record the frame timestamps rather than assuming a frame rate. Release from rest at about $10^{\circ }$, record through several complete swings, and repeat the release. For every same-side turning point, tabulate its time $t_n$ and angular amplitude $A_n=|\theta (t_n)|$. Measure the oscillation period from differences $t_{n+1}-t_n$; comparing several periods exposes drift. Interpolating around each observed extremum reduces frame discretization error. Check that late peaks still exceed the angular resolution; taking the natural logarithm of tiny, uncertain peaks amplifies their relative errors.

In the small-angle approximation, the damped pendulum becomes a damped harmonic oscillator. For $0<\gamma <2{\omega }_0$ and release from rest at ${\theta }_0$, its exact linear-model solution isThe decaying sinusoidal envelope has coefficient $|{\theta }_0|{\omega }_0/{\omega }_d$; because release is from rest, its successive same-side extrema occur at $t_n=n{T}_d$ and have $A_n=|{\theta }_0|e^{-\gamma n{T}_d/2}$. Thus the envelope and the sequence of actual peak heights are different descriptions of the same trace. The linear damping regimes and this solution are derived in Richard Fitzpatrick, Damped Oscillatory Motion, equations 82–93.

The experiment tests two predictions independently: nearly equal same-side peak spacing and a straight line when plotting $\ln (A_n/A_0)$ against $t_n-t_0$. Its predicted slope is $-\gamma /2$, givingUse several peaks and repeated trials rather than a single ratio, and compare the frequency identity with independently measured $\ell$ and $g$. For the illustrative choices $\ell =1.00\mathrm{m}$, $g=9.81m/{\mathrm{s}}^2$, $c=0.20$, and ${\theta }_0=10^{\circ }$, the linear model gives $\gamma \approx 0.6264{\mathrm{s}}^{-1}$, $T_d\approx 2.0162\mathrm{s}$, and $A_{n+1}/A_n\approx 0.5318$. The first two positive peaks after release are therefore approximately $5.318^{\circ }$ and $2.828^{\circ }$. The undamped small-amplitude oscillation period is about $2.0061\mathrm{s}$: substantial loss per swing need not imply a large change in period.

Repeat at a larger release angle to distinguish two limitations. The relative error in replacing $\sin \theta$ by $\theta$ is about ${\theta }^2/6$ for small angles: about $0.51\%$ at $10^{\circ }$, but about $17.3\%$ at $60^{\circ }$ when compared with $\theta$. Finite-amplitude swing spacing can therefore change as the swing shrinks even if $b$ is constant. Separately, friction at the pivot or quadratic drag can invalidate the assumed linear resisting torque; exponential peak decay should be checked, not imposed on all recorded motion. Peters discusses why real pendulum damping can require nonlinear terms. ^[ref] Repeatability, residual pivot motion, changes in damper temperature, and finite camera resolution are useful checks on this proposed measurement.

The normalized equilibrium points are $(\theta ,\dot{\theta })=(k\pi ,0)$. Even multiples of $\pi$ are attracting: their linearization has eigenvalues $(-c\pm \sqrt{c^2-4})/2$, giving a stable focus for $0<c<2$, repeated negative eigenvalues at $c=2$, and a stable node for $c>2$. These are the underdamped, critically damped, and overdamped small-displacement regimes. Odd multiples of $\pi$ have eigenvalues $(-c\pm \sqrt{c^2+4})/2$ of opposite signs and are saddle equilibria. The energy dissipation of a damped pendulum explains why trajectories in the phase portrait spiral inward near a stable focus instead of closing into the recurring orbit of an undamped simple pendulum.

This delay differential equation replaces instantaneous logistic inhibition by a population value one delay earlier. A positive history function of a delay differential equation gives positive solutions while they exist. The equilibria are zero and one. Linearization at one gives $x^{\prime }(t)=-x(t-\tau )$, with characteristic equation $\lambda =-e^{-\lambda \tau }$. Delayed sign reversal permits overshoots around one and, at the delayed logistic Hopf threshold, a periodic branch.

Arrival-time location estimates a hypocenter and origin time by comparing picked phase arrivals with predicted travel times. Linearization gives residual changes $\delta t_i=\delta t_0+\nabla T_i\cdot \delta \mathbf{x}$. Poorly distinct derivative columns cause tradeoffs, while incorrect wave speeds or phase labels cause systematic bias even with many stations.

A freely translating inner axle with prescribed positive angular velocity $\Omega$ and downward load $F>0$ has full-film journal bearing force balance equationswhere $R_r,R_t$ are the squeeze resistance of an eccentric journal bearing. Rotation supplies tangential fluid force $R_t\Omega r/2$. The equilibrium lies horizontally to the right and satisfies $F=R_t\Omega r/2$. Linearization there has imaginary eigenvalues: the prescribed rotating axle continuously supplies energy, so viscous dissipation does not imply settling. For an initially concentric axle, putting $e=r/\Delta$ and $K=6\pi \mu \Omega a^3/{\Delta }^2$ gives the closed orbitThis follows by eliminating time and integrating $x_e+3ex/[2(1-e^2)]=\Omega {\Delta }^2{R}_{0}e/[2F(1-e^2{)}^{3/2}]$. It applies to the leading inertia-free full-film model; adding cavitation, constraints or finite-length leakage changes the prediction.

A global bifurcation changes an invariant set through a global orbit connection or collision that cannot be identified by a single equilibrium point's linearization alone. For example, a global bifurcation through a saddle homoclinic orbit occurs when an outgoing unstable separatrix returns to the same saddle equilibrium's stable manifold. A smooth global Poincaré return map then combines with a local passage map to determine which side contains a periodic orbit. A heteroclinic orbit instead joins distinct equilibrium points; an energy-balance or Melnikov integral can locate its breaking in a nearly Hamiltonian system.

In the fast-water limit of an erodible-bed shallow water equations model, let $e=E^{\prime }(1)>0$, $q_1=q^{\prime }(1)>0$ and $\gamma >0$ be the suspended-load relaxation coefficient. Linearization and elimination give the dispersion relation $\lambda =-k^2[q_1+\beta e/(\gamma +ik)]/[3\delta +ik(F^2-1)]$. Its real part tends to the finite constant in the title, and phase velocity tends to $-q_1/(F^2-1)$. Thus the unstable upstream branch requires supercritical flow and enough sediment feedback to overcome friction; subcritical flow has the opposite propagation sign and decay at leading order. The near-critical limit needs a separate balance.

For $0<I_1<I_2<I_3$, the transverse linear stability rates about steady rotation with signed component $\Omega \ne 0$ areAxes one and three have transverse eigenvalues $\pm i{\nu }_1$ and $\pm i{\nu }_3$, producing bounded oscillations. Axis two has transverse eigenvalues $\pm {\lambda }_2$, giving growing and decaying exponential functions. A neutral longitudinal perturbation corresponds to a small change of $\Omega$. Generic intermediate-axis disturbances contain a growing component; a precisely chosen decaying direction is exceptional. For the normalized intermediate-axis model with $|\Omega |=1$, the values are ${\nu }_1=1/\sqrt{3}$, ${\nu }_3=1$ and ${\lambda }_2=1/\sqrt{3}$.

Proof outline: the linearization of the Euler equations for a torque-free rigid body about the third principal axis gives transverse perturbations obeyingand henceCyclic permutation gives the other rates. For nonlinear Lyapunov stability of the extreme rotations, use the torque-free rigid-body invariants: the conserved quadratic forms $L^2-2I_{1}E$ and $2I_{3}E-L^2$ have strictly positive transverse coefficients about axes one and three respectively. They bound transverse angular velocity, while $E$ bounds its longitudinal change near nonzero rotation. Imaginary transverse eigenvalues alone would not establish this nonlinear conclusion.^[ref]

The practical lesson is to choose a stable principal axis when using rotation to maintain an instrument's pointing direction. Both extreme principal moments of inertia are stable in the ideal torque-free model; if internal damping lowers rotational kinetic energy while approximately preserving angular momentum, the minimum energy $E=L^2/(2I_3)$ instead favors the largest principal moment of inertia. A control system must account for torque, changing mass distributions and flexibility before applying this ideal result.

The entry $K_{ij}$ is the expected number of newly infected type-$i$ individuals produced by a newly infected type-$j$ individual while the population is disease-free. In a compartmental linearization, $F$ describes new-infection production and $V$ transitions and removal. Provided the infection-lifetime subsystem is transient, its integrated occupation matrix is $V^{-1}$, yielding $K=F{V}^{-1}$. The basic reproduction number is its spectral radius. This construction is described in the original next-generation-matrix derivation.

Choose an inertia unit $I_{\ast }$ and an angular speed unit ${\Omega }_{\ast }>0$, so physical principal moments of inertia are $I_{\ast }(1,2,3)$, physical angular velocity is ${\Omega }_{\ast }\mathbf{\omega }$, and normalized time is $t={\Omega }_{\ast }{t}_{\mathrm{physical}}$. The Euler equations for a torque-free rigid body becomeThis is a physically admissible ideal planar mass distribution: $I_3=I_1+I_2$ is the perpendicular axis theorem. A book of finite thickness need not have these exact ratios.

For an original comparison, start with angular velocity $(1,\epsilon ,0)$, $(\epsilon ,1,0)$ or $(\epsilon ,0,1)$, respectively. With $\epsilon =0.01$, the transverse Euclidean norm stays small about axes one and three, but becomes order one about axis two. For the intermediate-axis start, its linearization givesThe hyperbolic functions contain a growing exponential function with rate ${\lambda }_2=1/\sqrt{3}$. Reducing $\epsilon$ tenfold delays reaching the same small transverse amplitude by approximately $\sqrt{3}\log 10\simeq 3.99$ normalized time units, once the growing term dominates. That delays the instability. Once the transverse angular velocity is comparable to the main component, the linearization ceases to describe the motion.

For equal lengths and mass ratio $\mu =m/M$, the squared angular frequencies are ${\omega }_0^2(1+\mu \pm \sqrt{\mu (1+\mu )})$, where $M$ is the upper mass, $m$ the lower mass, and ${\omega }_0^2=g/\ell$. The equal-mass case $\mu =1$ admits particularly simple mode shapes.

Use the ideal model of the double pendulum, with both angles measured from the downward vertical. Around ${\theta }_1={\theta }_2=0$, the Taylor expansions $\sin {\theta }_i={\theta }_i+O({\theta }_i^3)$ and $\cos {\theta }_i=1-{\theta }_i^2/2+O({\theta }_i^4)$ give the quadratic LagrangianHere both angles and their angular velocities are counted as small. Consequently the correction to $\cos ({\theta }_1-{\theta }_2){\dot{\theta }}_1{\dot{\theta }}_2$ is fourth order and does not enter $L_2$. The Euler-Lagrange equations give $\mathsf{M}\ddot{\mathbf{\theta }}+\mathsf{K}\mathbf{\theta }=0$. The mass matrix is not diagonal: an angular motion of the upper rod carries the lower mass too.^[ref]

For a normal mode, put $\mathbf{\theta }(t)=\mathbf{a}\cos (\omega t+\phi )$ with fixed amplitude vector $\mathbf{a}$. The resulting generalized eigenvalue problem for small oscillations is $\mathsf{K}\mathbf{a}={\omega }^2\mathsf{M}\mathbf{a}$. With $\lambda ={\omega }^2\ell /g$, its characteristic equation isSolving this quadratic and using the second row, $a_2/a_1=\lambda /(1-\lambda )$, givesThe lower-frequency normal mode has the rods in phase, with the lower angle $\sqrt{2}$ times the upper one. The higher-frequency normal mode has opposite signs, corresponding to an oscillation phase difference of $\pi$. These are ratios of absolute angles. The horizontal displacements are instead $x_1\simeq \ell {\theta }_1$ and $x_2\simeq \ell ({\theta }_1+{\theta }_2)$; their ratios are $1+\sqrt{2}$ and $1-\sqrt{2}$ respectively. Thus the lower mass moves horizontally opposite to the upper mass in the faster normal mode.

Both squared angular frequencies are positive. Their ratio satisfies ${\omega }_{+}/{\omega }_{-}=1+\sqrt{2}$, independently of $m$, $g$ and $\ell$. For a constructed model with $\ell =1\mathrm{m}$ and $g=9.81\mathrm{m}{\mathrm{s}}^{-2}$, the angular frequencies are approximately $2.3972$ and $5.7874{\mathrm{s}}^{-1}$, and the oscillation periods $2\pi /{\omega }_{\pm }$ are $2.6211$ and $1.0857\mathrm{s}$. Doubling both masses leaves these unchanged; multiplying both rod lengths by four doubles both oscillation periods.

The eigenvectors are not orthogonal in the ordinary dot product: ${\mathbf{a}}_{-}^{\mathsf{T}}{\mathbf{a}}_{+}=-1$. They satisfy mass-matrix orthogonality of normal modes, sinceWriting $\mathbf{\theta }={\mathbf{a}}_{-}{q}_{-}+{\mathbf{a}}_{+}{q}_{+}$ therefore removes the cross terms from both quadratic energies. The two normal coordinates satisfy ${\ddot{q}}_{\pm }+{\omega }_{\pm }^2{q}_{\pm }=0$. This is why the normal modes are useful: one coupled system becomes two independently evolving harmonic oscillators, and any small release can be decomposed into them. At larger angles, the discarded terms in the double pendulum equations of motion matter and this exact linear superposition no longer applies.

Translate the equilibrium to the origin by $z=(x-a,y-b{)}^T$. If the vector field is continuously differentiable, its first-order Taylor expansion isThe Jacobian matrix $M$ gives the linearization $\dot{z}=Mz$, whose general solution is $z(t)=e^{tM}z(0)$ using the matrix exponential. Its eigenvalues determine the exponential factors, while its eigenvectors and, when necessary, generalized eigenvectors determine their directions.

For distinct eigenvalues ${\lambda }_1,{\lambda }_2$ and corresponding eigenvectors $v_1,v_2$,For a real matrix with nonreal conjugate eigenvalues, choose conjugate coefficients to obtain real solutions. Equivalently, if $\lambda =a+ib$ and $v=p+iq$, a real pair of independent solutions isFor a repeated eigenvalue $\lambda$, there are two cases:

The factor $t$ in the second case is essential. These are exact solutions of the linearization; the remainder describes the error for the nonlinear differential equation. When no eigenvalue has zero real part, the hyperbolic equilibrium theory justifies the corresponding local qualitative classification. Zero real parts require additional nonlinear analysis.

Since the exponential factor is positive, equilibria must have $y=0$ and $x^4=1$. Thus $\boxed{(x,y)=(-1,0),(1,0)}$. At either point the Jacobian matrix is $\left(\begin{array}{cc}0& 1\\ -4x^3& 0\end{array}\right)$. At $(-1,0)$ its eigenvalues are $\pm 2$, so this is a saddle point. At $(1,0)$ they are $\pm 2i$: the linearization has a centre equilibrium.

The requested first-order approximation follows by putting $x=1+u$, $y=v$: $\dot{u}=v+O(\Vert (u,v){\Vert }^2)$ and $\dot{v}=-4u+O(\Vert (u,v){\Vert }^2)$. The linearized system conserves $4u^2+v^2$, giving the approximate ellipsesHowever, the printed nonlinear system has a weak repelling focus at $(1,0)$, rather than a centre. A centre in the linearization need not be a nonlinear centre. Here that distinction can be proved explicitly, without relying on the numerical drift in the figure.

Put $X=2(x-1)$, $Y=y$, $\tau =2t$ and $r^2=X^2+Y^2$. The weak focus of an exponential polynomial planar system has the local positive functionExpanding the vector field through cubic order and differentiating givesThe cubic correction cancels the cubic terms in the derivative of $H$, and the quartic correction leaves the positive radial term $3r^4/16$. Since $H=r^2/2+O(r^3)$ is positive and the angular velocity stays near $-1$, small nonzero trajectories spiral outward. Equivalently the reversed system has a strict local Lyapunov function and an asymptotically stable equilibrium. In particular, a periodic orbit sufficiently near the equilibrium is impossible, because $H$ would have to increase strictly and return to its starting value. This proves the failure of the printed centre assertion while retaining its correct linearized ellipse calculation.

The equilibrium equations are $x(1+x-2y)=0$ and $(x-1)y=0$. If $x=0$, the second equation gives $y=0$. Otherwise either $y=0$, giving $x=-1$, or $x=1$, giving $y=1$. Thus there are exactly three equilibria. The Jacobian isAt $(0,0)$ it is $\mathrm{diag}(1,-1)$: a saddle. At $(-1,0)$ it is $\left(\begin{array}{cc}-1& 2\\ 0& -2\end{array}\right)$, with eigenvalues $-1,-2$: a stable node. At $(1,1)$ it is $\left(\begin{array}{cc}1& -2\\ 1& 0\end{array}\right)$, with characteristic polynomial ${\lambda }^2-\lambda +2$ and eigenvalues $(1\pm i\sqrt{7})/2$: an unstable spiral. All three are hyperbolic, so the nonlinear local classifications follow from linearization.

The equilibrium equations are $x(1+x-2y)=0$ and $(x-1)y=0$. If $x=0$, the second equation gives $y=0$. Otherwise either $y=0$, giving $x=-1$, or $x=1$, giving $y=1$. Thus there are exactly three equilibria. The Jacobian isAt $(0,0)$ it is $\mathrm{diag}(1,-1)$: a saddle. At $(-1,0)$ it is $\left(\begin{array}{cc}-1& 2\\ 0& -2\end{array}\right)$, with eigenvalues $-1,-2$: a stable node. At $(1,1)$ it is $\left(\begin{array}{cc}1& -2\\ 1& 0\end{array}\right)$, with characteristic polynomial ${\lambda }^2-\lambda +2$ and eigenvalues $(1\pm i\sqrt{7})/2$: an unstable spiral. All three are hyperbolic, so the nonlinear local classifications follow from linearization.

A critical point of a dynamical system or equilibrium is a point $(a,b)$ with $f(a,b)=g(a,b)=0$, so the constant function $(x,y)=(a,b)$ is a solution. Put $u=x-a$, $v=y-b$. For a continuously differentiable vector field, its linearization iswith remainder $o(\sqrt{u^2+v^2})$ in the original nonlinear system. For distinct eigenvalues, $M$ is diagonalizable over $\mathbb{C}$, and the linear solutions are combinations of $e^{{\lambda }_{j}t}{w}_j$, where $w_j$ are corresponding eigenvectors. Real parts determine growth or decay; nonzero imaginary parts give rotation. For hyperbolic equilibria, meaning no eigenvalue has zero real part, the Hartman-Grobman theorem transfers the local stability and qualitative orbit pattern to the nonlinear system. A purely imaginary pair is not hyperbolic: it determines a linear centre, while higher-order terms decide whether the nonlinear system actually has a centre.

A centre equilibrium of the linearized real planar system hasAfter a real invertible linear change of coordinates, the system is a constant angular rotation. Its nonzero orbits are ellipses of period $2\pi /|\omega |$, so it is stable but not asymptotically stable. This is a classification of the linearization: a centre in the linearization need not be a nonlinear centre. In particular, imaginary eigenvalues alone do not prove a nonlinear centre equilibrium, which requires a punctured neighbourhood filled with periodic orbits. Part (c) illustrates exactly this limitation.

The Bendixson-Dulac theorem states that if $(P,Q)$ is a continuously differentiable planar vector field on a simply connected open region $D$, and a continuously differentiable Dulac function $B$ makeshave one weak sign throughout $D$, without vanishing identically on any nonempty open subset, then there is no nonconstant periodic orbit contained in $D$. To see why, a hypothetical simple periodic orbit bounds a region in $D$. The divergence theorem gives the integral of this expression as the outward flux of $B(P,Q)$ across the orbit. The flux is zero because the vector field is tangent to its orbit, whereas the one-sign area integral is nonzero, a contradiction. The usual stricter version assumes a strictly positive or strictly negative expression everywhere.

For the printed polynomial system, the origin is an equilibrium andThus its linearization is a centre equilibrium. For the full system take $B=1$. The divergence isOn the open disk $x^2+y^2<3$ this is nonpositive, and is negative wherever $x\ne 0$. Its zero set inside the disk is a line, so it cannot vanish on any open subset. The Bendixson-Dulac theorem therefore provesThis also exposes a terminology problem in the source: a nonlinear centre equilibrium has surrounding closed orbits, so the printed claim of a centre cannot hold in that sense. The imaginary eigenvalues establish only a centre of the linearized system. In fact the full equilibrium is an asymptotically stable weak focus, as the following direct check shows.

A quartic Lyapunov correction for a weak focus makes the nonlinear damping visible. Put $r^2=x^2+y^2$ and consider the Lyapunov functionIt is positive definite sufficiently close to zero. Differentiating along the polynomial vector field givesFor clarity, the quartic terms come from $-x^4$ in the derivative of $r^2/2$, and from applying $y{\partial }_x-x{\partial }_y$ to the quartic correction. Their sum is $-3x^4/8-3x^2{y}^2/4-3y^4/8$; all other terms have degree at least six. Therefore $\dot{V}<0$ in a sufficiently small punctured neighbourhood, proving asymptotic stability by Second Lyapunov theorem. Finallyso nonconstant nearby orbits turn clockwise as they converge. The equilibrium is a weak focus despite the centre in its linearization. This interpretation resolves both the local dynamics and the absence of closed orbits without changing the printed coefficients.