# Regularity and convergence of local first integrals of analytic differential
systems

Xiang Zhang  School of Mathematical Sciences, and MOE–LSC, Shanghai Jiao Tong
University, Shanghai 200240, P. R. China<EMAIL_ADDRESS>
###### Abstract.

Poincaré proved nonexistence of formal first integrals near a nonresonant
singularity of analytic autonomous differential systems. In the resonant case
with one zero eigenvalue and others nonresonant, there remains an open problem
on regularity and convergence of local first integrals.

Here we provide an answer to this problem. The system has always a local
$C^{\infty}$ first integral near the singularity when it is nonisolated. In
any finite dimensional space formed by analytic differential systems having
the same linear part at the singularity, either all the systems have local
analytic first integrals or only the systems in a pluripolar subset have local
analytic first integrals.

###### Key words and phrases:

Analytic differential systems; resonant singularities; $C^{\infty}$ local
first integrals; generic divergence of formal first integrals.

###### 2010 Mathematics Subject Classification:

37J35, 37C10, 37C27, 37C27, 34C14.

##  1\. Introduction and statement of the main results

For an analytic autonomous differential system defined in a connected open
region $\Omega$ of $\mathbb{R}^{n}$, it is a well known fact that the system
always have $n-1$ functionally independent analytic first integrals defined in
a neighborhood of any regular point in $\Omega$. So the system is completely
analytically integrable near any regular point. But near a singularity the
analytic differential system may have no analytic even no smooth first
integrals. The situation will be more intricate, see e.g. [1, 2, 12, 13, 16,
22]. The existence and regularity of first integrals of an analytic
differential system are related to dynamics on algebraic, geometry and
topology, see e.g. [3, 20].

An analytic autonomous differential system near a singularity can be written
in

(1.1) $\dot{x}=Ax+f(x),\qquad x\in\left(\mathbb{R}^{n},\ 0\right),$

where the dot denotes the derivative with respect to the time $t$, and
$f(x)=o(|x|)$ is an $n$–dimensional vector valued analytic function defined in
a neighborhood of the origin. Let
$\lambda=(\lambda_{1},\lambda_{2},\ldots,\lambda_{n})$ be the $n$–tuple of
eigenvalues of $A$. Hereafter, $\mathbb{R}$ is the set of real numbers, and
$\mathbb{Z}_{+}$ and $\mathbb{Q}_{+}$ are respectively the sets of nonnegative
integers and nonnegative rational numbers.

The existence and regularity of first integrals of system (1.1) depend on both
$f(x)$ and the eigenvalues of $A$. In case the eigenvalues of $A$ are not
resonant, Poincaré [15] proved the next result.

Theorem A. If the $n$–tuple of eigenvalues $\lambda$ do not resonant, i.e.

$\langle m,\lambda\rangle\neq 0,\quad
m=(m_{1},\ldots,m_{n})\in\mathbb{Z}_{+}^{n},\ \ |m|\geq 2,$

where $|m|=m_{1}+\ldots+m_{n}$ and $\langle\cdot,\cdot\rangle$ denotes the
inner product of two vectors, then the analytic differential system (1.1) has
neither analytic nor formal first integral in a neighborhood of the origin.

For a proof of Theorem A, see e.g. [7, 18], where Furta [7] and Shi and Li
[18] also extended the Poincaré result to quasihomogeneous and
semi–quasihomogeneous differential systems on their local nonintegrability in
the analytic or formal sense. Shi [17] provided a sufficient condition on the
nonexistence of meromorphic or formal meromorphic first integrals of the
analytic differential system (1.1) near a rational nonresonant singularity.

On existence of local analytic or meromorphic first integrals near a
singularity, there appeared some necessary conditions on the eigenvalues of
linearization of system (1.1) at the singularity. Chen et al [4] in 2008
obtained the optimal upper bound on the number of analytic or formal first
integrals of the analytic system (1.1). Cong et al [5] in 2011 got the optimal
upper bound on the number of meromorphic or formal meromorphic first integrals
for system (1.1), for quasihomogeneous and semi–quasihomogeneous analytic
differential systems near a singularity, and for analytic autonomous
differential system near a periodic orbit.

On determination of analytic integrability of analytic differential systems,
the pioneer work belongs to Poincaré, who proved that a planar analytic
differential system having a singularity with a pair of imaginary eigenvalues
is locally analytically integrable if and only if it is analytically orbitally
equivalent to its linear part. This result was extended to general higher
dimensional analytic integrable differential systems by Zung [24] in 2002, and
to analytic integrable Hamiltonian systems in the Liouvillian sense by Zung
[25] in 2005 via the method of the torus action, where he proved that any
local analytically Liouvillian integrable Hamiltonian system near a
singularity is analytically equivalent to its Birkhoff normal form. Zhang [23]
in 2011 proved the existence of analytic normalization for both general
analytic integrable differential system near a singularity, and diffeomorphism
near a fixed point, and provided the concrete forms of their analytic
integrable Poincaré–Dulac normal forms. This result was extended to partially
analytic integrable differential systems by Du et al [6] in 2016 with an
additional condition.

We note that all the above results provide either necessary conditions on the
existence and number of analytic (meromorphic, formal or formal meromorphic)
first integrals, or on equivalent characterization on existence of an analytic
normalization for analytic integrable systems. Whereas on the existence of
first integrals with suitable regularity, there are some general results for
two dimensional differential systems, see e.g. the book [16] and the
references therein. But for higher dimensional analytic differential systems
(1.1) with the singularity having resonant eigenvalues, the results are very
few, see e.g. [11].

In this paper the blanket assumptions are the following:

(1.2) $\lambda_{1}=0,\quad\langle\widetilde{m},\widetilde{\lambda}\rangle\neq
0,\quad\widetilde{m}\in\mathbb{Z}_{+}^{n-1},\ \ |\widetilde{m}|\geq 2,$

where $\widetilde{m}=(m_{2},\ldots,m_{n})$ and
$\widetilde{\lambda}=(\lambda_{2},\ldots,\lambda_{n})$. Li et al [11] provided
the conditions on the existence of a formal or an analytic first integral of
system (1.1) as stated in the following.

Theorem B. Assume that system (1.1) is analytic and satisfies (1.2). Then
system (1.1) has a formal first integral in $(\mathbb{R}^{n},\ 0)$ if and only
if the singularity $x=0$ is not isolated. Moreover, if $n=2$ the first
integral could be analytic.

Zhang [21] pursued the study on analyticity of the local first integral, and
got the next results.

Theorem C. For system (1.1) with $n-1$ nonresonant eigenvalues
$\widetilde{\lambda}$, the following statements hold.

  * $(a)$

If all the real parts of $\widetilde{\lambda}$ have the same sign, then system
(1.1) has an analytic first integral in $(\mathbb{R}^{n},\ 0)$ if and only if
the singularity $x=0$ is not isolated.

  * $(b)$

There exist analytic differential systems (1.1) with $\widetilde{\lambda}$
having both positive and negative real parts, which have no analytic first
integrals in $(\mathbb{R}^{n},\ 0)$.

According to the results in statement $(b)$ of Theorem C, it is natural to
appear the next open problems.

Problem 1: Are there $C^{\infty}$ first integrals in $(\mathbb{R}^{n},\ 0)$ of
system (1.1) under the conditions (1.2)?

Problem 2: What is the measure of the set of analytic system (1.1) satisfying
the conditions (1.2), which have analytic first integrals in
$(\mathbb{R}^{n},\ 0)$?

We first answer Problem 1 on existence of $C^{\infty}$ first integrals.

######  Theorem 1.1.

Under the conditions (1.2), the analytic differential system (1.1) has a
$C^{\infty}$ first integral in $(\mathbb{R}^{n},0)$ if and only if the
singularity at the origin is not isolated.

The next one is an answer to Problem 2. Let $\mathfrak{K}$ be the set of
analytic differential systems of form (1.1) with the same linear part, which
satisfies the conditions (1.2).

######  Theorem 1.2.

Assume that the origin is a nonisolated singularity of all systems in
$\mathfrak{K}$. Let $\mathcal{K}$ be any finite dimensional subspace of
$\mathfrak{K}$. The following statements hold.

  * $(a)$

If $\mathcal{K}$ contains an element, which has only formal but not analytic
first integral in a neighborhood of the origin, then all elements in
$\mathcal{K}$ except perhaps a pluripolar subset has also this property.

  * $(b)$

If $\mathcal{K}$ has a nonpluripolar subset whose any element has an analytic
first integral in a neighborhood of the origin, then all systems in
$\mathcal{K}$ have this property.

Recall that the condition on the nonisolated singularity at the origin is
necessary for the existence of local analytic or formal first integrals. A
pluripolar set is a subset of $\mathbb{C}^{m}$ for some $m\in\mathbb{N}$, and
it is of Lebesgue measure zero, see e.g. [9] or [14]. A detail definition on
pluripolar set will be given in Section 3, where we prove Theorem 1.2.

Theorem 1.2 has the next consequence.

######  Corollary 1.

Let $\mathfrak{P}$ be the set of polynomial differential systems of form (1.1)
with the same linear part and a uniformly bounded degree, which satisfy the
conditions (1.2) and have the origin as a nonisolated singularity. Then either
$\mathfrak{P}$ has all its elements having an analytic first integral in a
neighborhood of the origin, or $\mathfrak{P}$ has only a pluripolar subset
whose elements can have analytic first integrals near the origin.

We remark that the set $\mathfrak{P}$ is a finite dimensional space formed by
the coefficients of the nonlinear monomials of the systems in $\mathfrak{P}$.

Corollary 1 solves the open problem 2 for polynomial differential systems
(1.1). Since a pluripolar set has Lebesgue measure zero, if $\mathfrak{P}$ has
a positive Lebesgue measure subset, whose each element has an analytic first
integral in a neighborhood of the origin, then all systems in $\mathfrak{P}$
have analytic first integrals near the origin. In other words, if
$\mathfrak{P}$ has an element which has only divergent formal first integral
in a neighborhood of the origin, then generic systems in $\mathfrak{P}$ have
only divergent formal first integrals near the origin.

All analytic differential systems of type (1.1) form an infinitely dimensional
space in the sense of Theorem 1.2 and Corollary 1. In this case, it is still
an open problem whether the generic systems have only divergent first
integrals at the origin.

This paper is organized as follows. The proof of Theorems 1.1 is given Section
2, and the proofs of Theorem 1.2 and Corollary 1 will be presented in Section
3.

##  2\. Proof of Theorem 1.1

Necessity. It follows from [11, 21]. Here for completeness we present its
proof. By the assumption of the theorem and [6, Theorem 1 $(b)$] (see the next
Proposition 1), system (1.1) has a local smooth first integral of the form
$F(x)=x_{1}+o(x)$ in $(\mathbb{R}^{n},0)$. Taking the invertible change of
coordinates $v=\Phi(x)$ in $(\mathbb{R}^{n},0)$ with

$v_{1}=F(x),\ \ \mathbf{v}_{2}=(v_{2},\ldots,v_{n})^{\tau}=(x_{2},\ \ldots,\
x_{n})^{\tau},$

with $\tau$ representing the transpose of a matrix, system (1.1) in this new
coordinate system becomes

$\dot{v}_{1}=0,\quad\dot{\mathbf{v}}_{2}=B\mathbf{v}_{2}+\mbox{h.o.t.}$

This last system has the line of singularities $\mathbf{v}_{2}=0$, and
consequently system (1.1) has a smooth curve fulfilling singularities and
passing the origin. The necessity follows.

For readers’ convenience we recall [6, Theorem 1 $(b)$], which will also be
used in the proof of Theorem 1.2. For system (1.1) set

$\mathcal{R}:=\\{m\in\mathbb{Z}_{+}^{n}|\ \langle m,\lambda\rangle=0\\}.$

######  Proposition 1.

Assume that $\mathcal{R}$ has $d<n$ $\mathbb{Q}_{+}$–linearly independent
elements. If system (1.1) has $d$ functionally independent analytic or formal
first integrals, then it has $d$ functionally independent first integrals of
the form

$H_{1}(x)=x^{\alpha_{1}}+h_{1}(x),\,\ldots,\,H_{d}(x)=x^{\alpha_{d}}+h_{d}(x),$

where $\alpha_{1},\ldots,\alpha_{d}$ are $\mathbb{Q}_{+}$–linearly independent
elements of $\mathcal{R}$, and each $h_{j}(x)$, $j=1,\,\ldots,\,d$, consists
of nonresonant monomials in $x$ of degree larger than $|\alpha_{j}|$.

Sufficiency. Without loss of generality we can assume that the linear part $A$
of system (1.1) is of the form

$A=\left(\begin{array}[]{cc}0&\widetilde{\mathbf{0}}\\\
{\widetilde{\mathbf{0}}}^{\tau}&B\end{array}\right),$

where $\widetilde{\mathbf{0}}$ is the $n-1$ dimensional row vector with its
entries all zeros, and $B$ is a square matrix of order $n-1$ with the $n-1$
tuple of eigenvalues $\widetilde{\lambda}$ and is in lower triangular normal
form. Under this consideration system (1.1) is of the form

(2.1) $\begin{split}\dot{x}_{1}&=f_{1}(x),\\\
\dot{\mathbf{x}}_{2}&=B\mathbf{x}_{2}+\mathbf{f}_{2}(x).\end{split}$

Hereafter we denote by $\mathbf{a}_{2}$ the $n-1$ dimensional column vector
$(a_{2},\ldots,a_{n})^{\tau}$ for any letter $a$.

Let $\mathbf{x}_{2}=\mathbf{\phi}(x_{1})$ be the unique solution of the
equation $B\mathbf{x}_{2}+\mathbf{f}_{2}(x)=0$ in a neighborhood of the
origin, which is analytic and can be obtained by the implicit function
theorem. Since the singularity at the origin is not isolated, it forces that
$f_{1}(x_{1},\mathbf{\phi}(x_{1}))\equiv 0$ in a neighborhood of the origin.
Taking the coordinate change of variables

(2.2) $y_{1}=x_{1},\quad\mathbf{y}_{2}=\mathbf{x}_{2}-\mathbf{\phi}(x_{1}),$

system (2.1) can be written in an analytically equivalent way as

(2.3) $\begin{split}\dot{y}_{1}&=F_{1}(y)\mathbf{y}_{2},\\\
\dot{\mathbf{y}}_{2}&=B\mathbf{y}_{2}+F_{2}(y)\mathbf{y}_{2},\end{split}$

where $F_{1}$ is an $n-1$ dimensional row vector valued function and $F_{2}$
is a square matrix valued function of order $n-1$. Since the transformation
(2.2) is invertible and analytic, to prove Theorem 1.1 is equivalent to prove
that system (2.3) has a $C^{\infty}$ first integral in $(\mathbb{R}^{n},0)$.

Step 1\. Normal form. By Poincaré–Dulac normal form theorem (see e.g. [22])
and the assumption (1.2), there exists a near identity distinguished formal
change of coordinates

(2.4) $y=\widehat{\Phi}(z)=z+\widehat{\phi}(z),$

under which system (2.3) is transformed to its distinguished normal form

(2.5) $\begin{split}\dot{z}_{1}&=0,\\\
\dot{\mathbf{z}}_{2}&=B\mathbf{z}_{2}+\widehat{G}_{2}(z)\mathbf{z}_{2},\end{split}$

where $\widehat{G}_{2}$ is a matrix valued formal series in $z$ of order
$n-1$. We remark that in general by Poincaré–Dulac normal form theorem the
first equation of (2.5) should have the form
$\dot{z}_{1}=\widehat{G}_{1}(z_{1})$, with $\widehat{G}_{1}(z_{1})$ a formal
series in $z_{1}$, then the nonisolate of the singularity at the origin forces
$\widehat{G}_{1}(z_{1})\equiv 0$, otherwise we take a sufficiently high order
cut–off of (2.4) which will induce that the transformed system has the origin
as a singularity of finite multiplicity, a contradiction with the nonisolate
of the singularity. Since the normalization is chosen to be distinguished, it
compels that $\widehat{\phi}$ consists of nonresonant monomials and so
$\widehat{\phi}(z)=\widehat{\phi}(\mathbf{z}_{2})$.

By the Borel lemma, there exists a $C^{\infty}$ function
$\Phi(\mathbf{z}_{2})$ such that
$\mbox{jet}_{0}^{\infty}\Phi(\mathbf{z}_{2})=\widehat{\phi}(\mathbf{z}_{2})$.
Hereafter $\mbox{jet}_{0}^{\infty}\Phi$ denotes the Taylor series of $\Phi$ at
the origin. Taking the $C^{\infty}$ change of coordinates
$y=z+\Phi(\mathbf{z}_{2})$, system (2.3) is $C^{\infty}$ equivalent, in a
neighborhood of the origin, to

(2.6) $\begin{split}\left(\begin{array}[]{c}\dot{z}_{1}\\\
\dot{\mathbf{z}}_{2}\end{array}\right)&=\left(I+\partial_{z}\Phi(\mathbf{z}_{2})\right)^{-1}\left(\begin{array}[]{c}F_{1}(y)\\\
B+F_{2}(y)\end{array}\right)\left(\mathbf{z}_{2}+\mathbf{\Phi}_{2}(\mathbf{z}_{2})\right)\\\
&=\left(\begin{array}[]{c}W_{1}(z)\mathbf{z}_{2},\\\
(B+W_{2}(z))\mathbf{z}_{2}\end{array}\right)\end{split}$

where $W_{1}$ is an $n-1$ dimensional $C^{\infty}$ row vector valued function
and $W_{2}$ is a $C^{\infty}$ square matrix valued function of order $n-1$,
$\partial_{z}\Phi(\mathbf{z}_{2})$ is the Jacobian matrix of $\Phi$ in $z$ and
$\mathbf{\Phi}_{2}$ is the last $n-1$ components of $\Phi$. Note that
$W_{1}(z)$ is infinitely flat at the origin, i.e.
$\mbox{jet}_{0}^{\infty}W_{1}(z)=0$. If we can prove that system (2.6) has a
$C^{\infty}$ first integral in $(\mathbb{R}^{n},0)$, then system (2.3) and
consequently system (1.1) has a $C^{\infty}$ first integral in
$(\mathbb{R}^{n},0)$.

Again applying the Borel lemma to $\widehat{G}_{2}(z)$ in (2.5) provides a
$C^{\infty}$ square matrix valued function $V_{2}(z)$ of order $n-1$ such that
$\mbox{jet}_{0}^{\infty}V_{2}(z)=\widehat{G}_{2}(z)$. Consider the system

(2.7) $\begin{split}\dot{z}_{1}&=0,\\\
\dot{\mathbf{z}}_{2}&=B\mathbf{z}_{2}+V_{2}(z)\mathbf{z}_{2},\end{split}$

Clearly system (2.7) has the $C^{\omega}$ first integral $H=z_{1}$. If we can
prove $C^{\infty}$ equivalence between the two systems (2.6) and (2.7), then
system (2.6) will have a $C^{\infty}$ first integral in $(\mathbb{R}^{n},0)$.

Step 2\. $C^{\infty}$ equivalence between systems (2.6) and (2.7).

Denote by $\mathcal{W}$ the vector field associated to system (2.6) and by
$\mathcal{V}$ the vector field associated to system (2.7). Observing that
systems (2.6) and (2.7) both have the $z_{1}$–axis as their center manifold,
and that $\mathcal{W}-\mathcal{V}=0$ when restricted to the center manifold
$\mathbf{z}_{2}=0$. Moreover, the last construction shows that
$\mathcal{W}-\mathcal{V}$ is infinitely flat at the origin.

The condition (1.2) implies that all the eigenvalues of $B$ have nonvanishing
real parts. By the stable manifold theorem system (2.6) has an $s$–dimensional
$C^{\infty}$ stable manifold $M_{w}^{s}$ and a $u$–dimensional $C^{\infty}$
unstable manifold $M_{w}^{u}$, and system (2.7) has also an $s$–dimensional
$C^{\infty}$ stable manifold $M_{v}^{s}$ and a $u$–dimensional $C^{\infty}$
unstable manifold $M_{v}^{u}$. Taking the stable manifold $M_{w}^{s}$ and
unstable manifold $M_{w}^{u}$ as coordinate spaces, and denoting them by $S$
and $U$ respectively. Denote by $C$ the $1$–dimensional center manifold, i.e.
the $z_{1}$–axis. Then $\mathbb{R}^{n}=C\oplus S\oplus U$. Denote by $z_{1},\
z_{s}$ and $z_{u}$ the coordinates on $C$, $S$ and $U$, and we still use $z$
to represent $(z_{1},z_{s},\ z_{u})$.

Note that all the above manipulations are local. Take $\varphi:\
\mathbb{R}\rightarrow[0,1]$ be a $C^{\infty}$ function satisfying
$\varphi(s)=1$ for $s\in[0,\ 1/2]$ and $\varphi(s)=0$ for $s\geq 1$. Set

$\widetilde{\mathcal{W}}(z)=\left(\begin{array}[]{c}\sigma^{-1}\varphi(\|z\|)W_{1}(\sigma
z)\mathbf{z}_{2}\\\ \left(B+\sigma^{-1}\varphi(\|z\|)W_{2}(\sigma
z)\right)\mathbf{z}_{2}\end{array}\right),$

and

$\widetilde{\mathcal{V}}(z)=\left(\begin{array}[]{c}0\\\
\left(B+\sigma^{-1}\varphi(\|z\|)V_{2}(\sigma
z)\right)\mathbf{z}_{2}\end{array}\right),$

where $\sigma>0$ is a parameter, and it will be taken to be suitably small. We
have the following facts:

  * •

$\widetilde{\mathcal{W}}$ and $\widetilde{\mathcal{V}}$ are two globally
defined vector fields in $\mathbb{R}^{n}$.

  * •

In the region $\|z\|\leq 1/2$, after the linear change of coordinates $\sigma
z\rightarrow z$ the vector fields $\widetilde{W}$ and $\widetilde{V}$ are
respectively transformed to $W$ and $V$.

So, to prove the $C^{\infty}$ equivalence between (2.6) and (2.7) is
equivalent to verify the $C^{\infty}$ equivalence between the vector fields
$\widetilde{\mathcal{W}}$ and $\widetilde{\mathcal{V}}$.

For doing so, set

$R(z)=\widetilde{\mathcal{W}}(z)-\widetilde{\mathcal{V}}(z).$

Claim 1: $R$ has a $C^{\infty}$ decomposition as

$R(z)=R^{+}(z)+R^{-}(z),$

with $R^{+},\ R^{-}\in C^{\infty}(\mathbb{R}^{n})$ and satisfying

$\begin{split}\partial^{\alpha}R^{+}(z)&=0,\quad z\in C\oplus U,\\\
\partial^{\alpha}R^{-}(z)&=0,\quad z\in C\oplus S.\end{split}$

for all $\alpha=(\alpha_{1},\ldots,\alpha_{n})\in\mathbb{Z}_{+}^{n}$, where
$\partial^{\alpha}R=\dfrac{\partial^{\alpha}R}{\partial
z_{1}^{\alpha_{1}}\partial z_{2}^{\alpha_{2}}\ldots\partial
z_{n}^{\alpha_{n}}}$. Indeed, taking

$r_{\alpha}(z)=\left\\{\begin{array}[]{cl}\partial^{\alpha}R(z),&z\in C\oplus
S,\\\ 0,&z\in C\oplus U,\end{array}\right.$

Since $\partial^{\alpha}R(z)=0$ on $C$ for all $\alpha\in\mathbb{Z}_{+}^{n}$,
it follows that there exists a $C^{\infty}$ function $R^{+}$ such that
$\partial^{\alpha}R^{+}(z)=r_{\alpha}(z)$ for $z\in(C\oplus S)\cup(C\oplus
U)$. Here we have used the next facts (see e.g. [10, Lemma 1.14]):

For two linear subspaces $V$ and $W$ in $\mathbb{R}^{n}$ satisfying
$V+W=\mathbb{R}^{n}$, and two functions $g,h\in C^{\infty}(\mathbb{R}^{n},0)$
verifying

$\frac{\partial^{\alpha}g}{\partial
z^{\alpha}}(z)=\frac{\partial^{\alpha}h}{\partial
z^{\alpha}}(z),\quad\mbox{for all }\ z\in V\cap W,\ \ \
\alpha\in\mathbb{Z}_{+}^{n},$

then there exists a $C^{\infty}$ function $F$ defined on a neighborhood of the
origin such that

$\frac{\partial^{\alpha}F}{\partial
z^{\alpha}}(z)=\left\\{\begin{split}\dfrac{\partial^{\alpha}g}{\partial
z^{\alpha}}(z),&\qquad z\in V,\\\ \dfrac{\partial^{\alpha}h}{\partial
z^{\alpha}}(z),&\qquad z\in W,\end{split}\right.\quad\mbox{ for all }\ \
\alpha\in\mathbb{Z}_{+}^{n}.$

Set $R^{-}=R-R^{+}$. Then $\partial^{\alpha}R^{-}=0$ on $C\oplus S$ for all
$\alpha\in\mathbb{Z}_{+}^{n}$. This proves the claim.

In order to prove the $C^{\infty}$ equivalence between
$\widetilde{\mathcal{W}}$ and $\widetilde{\mathcal{V}}$, we only need to prove
the $C^{\infty}$ equivalence between $\widetilde{\mathcal{V}}$ and
$\widetilde{\mathcal{V}}+R^{-}$ and the $C^{\infty}$ equivalence between
$\widetilde{\mathcal{V}}+R^{-}$ and $\widetilde{\mathcal{W}}\
$$(=\widetilde{\mathcal{V}}+R^{-}+R^{+})$. Next we prove the former. The
latter follows from similar arguments as those in the proof of the former.

For proving the $C^{\infty}$ equivalence between $\widetilde{\mathcal{V}}$ and
$\widetilde{\mathcal{V}}+R^{-}$, we adopt the homological method (see e.g. [8,
Chapter 1, Section 4] or [10, Section 2.2]). Set

$\mathcal{V}_{s}=\widetilde{\mathcal{V}}+sR^{-},\quad s\in[0,1].$

The homological method shows that if the homological equation

(2.8) $[h,\mathcal{V}_{s}]=R^{-},$

has a $C^{\infty}$ vector valued solution $h$, then the vector fields
$\mathcal{V}$ and $\mathcal{V}+R^{-}$ is $C^{\infty}$ equivalent, where
$[\cdot,\ \cdot]$ is the Lie bracket of two vector fields. Recall that for any
two smooth $n$–dimensional vector fields $\mathcal{A}$ and $\mathcal{B}$ in
variables $z$, their Lie bracket is

$[\mathcal{A},\
\mathcal{B}]:=\partial_{z}\mathcal{B}\mathcal{A}-\partial_{z}\mathcal{A}\mathcal{B},$

which is again a vector field, and

(2.9) $[\mathcal{A},\
\mathcal{B}]=-\left.\frac{d}{dt}\mathcal{B}_{t}\right|_{t=0},$

with $\mathcal{B}_{t}:=\left(\psi_{\mathcal{A}}^{t}\right)_{*}\mathcal{B}$,
where $\psi_{\mathcal{A}}^{t}$ is the flow of the vector field $\mathcal{A}$
and $\left(\psi_{\mathcal{A}}^{t}\right)_{*}$ is the tangent map of the flow.

Let $\psi_{t}$ be the flow of the vector field $\mathcal{V}_{s}$ satisfying
$\psi_{0}(z)=z$, and let $P(z)$ be the time $1$ map of the flow, i.e.
$P(z)=\psi_{1}(z)$. Set $\Phi(t,z)=\partial_{z}\psi_{t}(z)$. By the
fundamental theory of ordinary differential equations on variational
equations, it follows that $\Phi(t,z)$ is the fundamental solution matrix of
the matrix valued variational equation

$\frac{dZ}{dt}=\partial_{w}\mathcal{V}_{s}(w)|_{w=\psi_{t}(z)}Z,\qquad
Z(0)=I,$

of the differential system $\dot{z}=\mathcal{V}_{s}$ along its solution
$z(t)=\psi_{t}(z)$, where $I$ is the identity matrix. Set

(2.10) $h(z)=\int_{0}^{-\infty}\Phi^{-1}(s,z)R^{-}(\psi_{s}(z))ds,$

and according to (2.9), set $h_{t}(z)=\left(\psi_{t}\right)_{*}h(z)$. Direct
calculations verify that if the series $h$ defined in (2.10) is convergent,
then one has

(2.11) $\begin{split}h_{t}(z)&=\frac{\partial\psi_{t}}{\partial
z}h\circ\psi_{-t}(z)=\Phi(t,\psi_{-t}(z))h(\psi_{-t}(z))\\\
&=\Phi(t,\psi_{-t}(z))\int_{0}^{-\infty}\Phi^{-1}(s,\psi_{-t}(z))R^{-}(\psi_{s}(\psi_{-t}(z)))ds\\\
&=\Phi^{-1}(-t,z)\int_{0}^{-\infty}\Phi^{-1}(s,\psi_{-t}(z))R^{-}(\psi_{s-t}(z))ds\\\
&=\int_{0}^{-\infty}\Phi^{-1}(s-t,z)R^{-}(\psi_{s-t}(z))ds\\\
&=\int_{-t}^{-\infty}\Phi^{-1}(\tau,z)R^{-}(\psi_{\tau}(z))d\tau,\end{split}$

where in the fourth equality we have used the identity

$I=\frac{\partial}{\partial
z}\left(\psi_{t}\circ\psi_{-t}(z)\right)=\left.\frac{\partial\psi_{t}(\nu)}{\partial\nu}\right|_{\nu=\psi_{-t}(z)}\frac{\partial}{\partial
z}\psi_{-t}(z)=\Phi(t,\psi_{-t}(z))\Phi(-t,z),$

and in the fifth equality we have applied the identity

$\Phi(s-t,z)=\frac{\partial\psi_{s-t}(z)}{\partial
z}=\left.\frac{\partial\psi_{s}(\nu)}{\partial\nu}\right|_{\nu=\psi_{-t}(z)}\frac{\partial\psi_{-t}(z)}{\partial
z}=\Phi(s,\psi_{-t}(z))\Phi(-t,z).$

Then we get from (2.11) that

$\left.\frac{d}{dt}h_{t}(z)\right|_{t=0}=\Phi^{-1}(0,z)R^{-}(\psi_{0}(z))=R^{-}(z).$

This together with (2.9) verifies that $h$ is a solution of the homological
equation (2.8).

The remaining is to prove the convergence of the integral defining $h$ and its
$C^{\infty}$ smoothness. Write the integral in (2.10) as

$h(z)=\sum\limits_{j=0}\limits^{\infty}\int_{-j}^{-j-1}\Phi^{-1}(s,z)R^{-}(\psi_{s}(z))ds.$

Direct calculations yield

$\displaystyle\int_{-j}^{-j-1}$
$\displaystyle\Phi^{-1}(s,z)R^{-}(\psi_{s}(z))ds$
$\displaystyle=\left(\partial_{z}P^{-j}(z)\right)^{-1}\int_{0}^{-1}\Phi^{-1}(s,P^{-j}(z))R^{-}(\psi_{s}(P^{-j}(z)))ds,$

where we have used

$\Phi(s,z)=\frac{\partial\psi_{s}(z)}{\partial z}=\frac{\partial}{\partial
z}\left(\psi_{s+j}\circ\psi_{-j}(z)\right)=\Phi(s+j,P^{-j}(z))\partial_{z}P^{-j}(z).$

Hence $h(z)$ can be written as

(2.12)
$h(z)=\sum\limits_{j=0}\limits^{\infty}\left(\partial_{z}P^{-j}(z)\right)^{-1}G\left(P^{-j}(z)\right),$

where

$G(w)=\int_{0}^{-1}\Phi^{-1}(s,w)R^{-}(\psi_{s}(w))ds.$

Observe that to prove the $C^{\infty}$ smoothness of $h$ is equivalent to
prove that $h$ is $C^{k}$ for any $k\in\mathbb{N}$.

For any given $k\in\mathbb{N}$, by the construction of the vector fields
$\mathcal{V}_{s}$ it follows that

(2.13) ${\mathcal{V}_{s}}\overset{\sigma\rightarrow
0}{\longrightarrow}Az\quad\mbox{ and }\quad\mathcal{V}_{s}(z)=Az\ \ \ \
\mbox{for}\ \ \ \ \|z\|\geq 1,$

and consequently

(2.14) $\lim\limits_{\sigma\rightarrow
0}P^{-1}(z)=e^{-A}z\quad\mbox{and}\quad\|P^{-1}\|\leq e^{\|A\|}+1.$

Hereafter, for a matrix or a matrix function $W(z)$ we denote by $\|W(z)\|$
the maximum of the absolute values of its elements, and

$\|W\|:=\max\limits_{\|z\|\leq 1}(\|W(z)\|/\|z\|).$

According to (2.13) and (2.14) together with the definition of
$\widetilde{\mathcal{W}},\ \widetilde{\mathcal{V}}$ and $\mathcal{V}_{s}$, we
can choose $\sigma>0$ suitably small such that

(2.15) $\max\limits_{0\leq|m|\leq k}\\{\|\partial^{m}P^{-1}\|,\
m\in\mathbb{Z}_{+}^{n}\\}\leq N_{0}=2e^{\|A\|}.$

For proving the convergence of all the partial derivatives up to $k$th order
of the series $h$ in (2.12), we need to estimate $\|\partial^{m}P^{-j}\|$ for
$j\geq 1$, and all $m\in\mathbb{Z}_{+}^{n}$ and $0\leq|m|\leq k$.

Claim 2. There exists a constant $K_{k}\geq 1$ depending on $k$ such that

(2.16) $\max\limits_{0\leq|m|\leq k}\\{\|\partial^{m}P^{-j}\|\\}\leq
K_{k}^{j},\quad j\in\mathbb{N}.$

Indeed, for $j=1$ the inequality (2.15) shows that $K_{k}=N_{0}$ satisfies the
requirement. For $j>1$ we adopt induction. Assume that (2.16) holds for all
positive integers no more than $j$, we consider $j+1$. For any
$m\in\mathbb{Z}_{+}^{n}$ satisfying $0\leq|m|\leq k$,

(2.17) $\begin{split}\partial^{m}P^{-j-1}(z)=&\partial^{m}\left(P^{-1}\circ
P^{-j}(z)\right)\\\ =&\left.\frac{\partial
P^{-1}(\mu)}{\partial\mu}\right|_{\mu=P^{-j}(z)}\partial^{m}P^{-j}(z)+P_{m}(z),\end{split}$

where

$P_{m}(z)=\Lambda_{m}\left(\left.\frac{\partial^{\ell}P^{-1}(\mu)}{\partial\mu^{\ell}}\right|_{\mu=P^{-j}(z)},\partial^{\rho}P^{-j}(z),2\leq|\ell|\leq|m|,\
1\leq|\rho|<|m|\right),$

with $\Lambda_{m}$ an $n$–dimensional vector valued polynomial in its
variables with nonnegative coefficients. For the given $m$, let $p(m)$ be the
maximum number of monomials in the $n$ components of $\Lambda_{m}$, and set
$p_{k}=\max\limits_{0\leq|m|\leq k}p(m)$. Then for $m\in\mathbb{Z}_{+}^{n}$
and $0\leq|m|\leq k$ one has

$\displaystyle\|P_{m}(z)\|$ $\displaystyle\leq
p(m)\left(\max\limits_{2\leq|\ell|\leq|m|}\left\|\left.\frac{\partial^{\ell}P^{-1}(\mu)}{\partial\mu^{\ell}}\right|_{\mu=P^{-j}(z)}\right\|\right)^{k}\left(\max\limits_{1\leq|\rho|<|m|}\|\partial^{\rho}P^{-j}(z)\|\right)^{k}$
$\displaystyle\leq
p_{k}N_{0}^{k}\left(\max\limits_{1\leq|\rho|<|m|}\|\partial^{\rho}P^{-j}(z)\|\right)^{2k}.$

This together with (2.17) and the induction assumption gives that

$\displaystyle\|\partial^{m}P^{-j-1}(z)\|$ $\displaystyle\leq
N_{0}\left(\max\limits_{1\leq|\rho|\leq|m|}\|\partial^{\rho}P^{-j}(z)\|\right)^{2}+p_{k}N_{0}^{k}\left(\max\limits_{1\leq|\rho|\leq|m|}\|\partial^{\rho}P^{-j}(z)\|\right)^{2k}$
$\displaystyle\leq
2p_{k}N_{0}^{k}\left(\max\limits_{1\leq|\rho|\leq|m|}\|\partial^{\rho}P^{-j}(z)\|\right)^{2k}$
$\displaystyle\leq 2p_{k}N_{0}^{k}\left(K_{|m|}^{2k}\right)^{j}\leq
K_{k}^{j+1},$

where $K_{k}=\max\\{2p_{k}N_{0}^{k},K_{|m|}^{2k}\\}$, and in the last
inequality we have used the induction assumption on $j$. Note that $p_{k}$
depends only on the dimension of the system and the order $k$ of the
derivative, and $N_{0}$ depends on $A$ and $K_{|m|}$ depends on $k$ and $A$.
Consequently $K_{k}$ depends only on $k$. The claim follows.

Finally we prove that $h(z)$ in (2.12) is $C^{k}$ for any $k\in\mathbb{N}$.

Claim 3. For any $m\in\mathbb{Z}_{+}^{n}$ satisfying $|m|=k$, the series

(2.18)
$\sum\limits_{j=0}\limits^{\infty}\partial^{m}\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}G\left(P^{-j}(z)\right)\right)$

is absolutely convergent. For proving this claim we need the next facts:

  1. $(f_{1})$

Derivatives of product of two functions

$\displaystyle\partial^{m}$
$\displaystyle\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}G\left(P^{-j}(z)\right)\right)$
$\displaystyle=\sum\limits_{\begin{subarray}{c}\ell\in\mathbb{Z}_{+}^{n}\\\
\ell\prec
m\end{subarray}}c_{\ell}\partial^{\ell}\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}\right)\partial^{m-\ell}G\left(P^{-j}(z)\right),$

where $c_{\ell}$’s are nonnegative integers depending on the expansion, and
$\ell\prec m$ represents that any element of $\ell$ is no more than the
corresponding one of $m$. Set $c_{k}=\max\\{c_{\ell};\ |\ell|\leq k,\
\ell\in\mathbb{Z}_{+}^{n}\\}$. Since $|m|\leq k$, the terms in the summation
are finite, denote by $s_{k}$ the maximum number.

  2. $(f_{2})$

Taking $\mu_{M}$ to be one half of the minimum of the positive real parts of
eigenvalues of $B$, and $\mathcal{P}_{U}(z)$ the projection of
$z\in\mathbb{R}^{n}$ to the unstable subspace $U$. Then

$\|\mathcal{P}_{U}\psi_{-t}(z)\|\leq\kappa
e^{-t\mu_{M}}\|\mathcal{P}_{U}z\|=\kappa e^{-t\mu_{M}}\|z_{u}\|,$

where $\kappa$ is a positive constant.

  3. $(f_{3})$

Restricting to the center–stable manifold

$\mbox{jet}_{z}^{\infty}G(z)=0,\quad\mbox{for arbitrary }\ z\in C\oplus S,$

with $G$ given in (2.12), where we have used the property that $R^{-}$
satisfies, i.e. $\mbox{jet}_{z}^{\infty}R^{-}(z)=0$ for arbitrary $z\in
C\oplus S$.

By $(f_{3})$ it follows that for any $m\in\mathbb{Z}_{+}^{n}$, $|m|\leq k$ and
any $N\in\mathbb{N}$

(2.19) $\|\partial^{m}G(z)\|\leq a_{m}\|z_{u}\|^{N},$

where $a_{m}$ is a positive constant depending only on $m$. Hence, by
$(f_{2})$ one has

(2.20) $\|\partial^{m-\ell}G\left(P^{-j}(z)\right)\|\leq
a_{m-\ell}\|\mathcal{P}_{U}P^{-j}(z)\|^{N}\leq a_{m-\ell}\kappa
e^{-jN\mu_{M}}\|z_{u}\|^{N}.$

According to $(f_{1})$, now we need to estimate
$\left\|\partial^{\ell}\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}\right)\right\|$.
In fact, we prove that there exists a positive constant $M_{k}$ depending only
on $k$ such that for any $|\ell|\leq k$

(2.21)
$\left\|\partial^{\ell}\left(\partial_{z}P^{-j}(z)\right)^{-1}\right\|\leq
M_{k}^{j}.$

Its proof is similar to that of (2.16). We present a sketch one via induction
on $|\ell|\geq 1$ and the identity

$\partial^{\ell}\left(\left(\partial_{z}P^{-j}(z)\right)\left(\partial_{z}P^{-j}(z)\right)^{-1}\right)=\partial^{\ell}(I)=0.$

For $|\ell|=1$ one has

$\partial^{\ell}\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}\right)=-\left(\partial_{z}P^{-j}(z)\right)^{-1}\partial^{\ell}\left(\partial_{z}P^{-j}(z)\right)\left(\partial_{z}P^{-j}(z)\right)^{-1}.$

Note from (2.14) that

$\lim\limits_{\sigma\rightarrow 0}\partial_{z}P^{-j}(z)=e^{-jA}.$

So there exists a $\sigma>0$ small enough such that

$\left\|\left(\partial_{z}P^{-j}(z)\right)^{-1}\right\|=2e^{j\|A\|}.$

Consequently, for $|\ell|=1$ one has

(2.22)
$\left\|\partial^{\ell}\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}\right)\right\|\leq
4e^{2j\|A\|}\left\|\partial^{\ell}\left(\partial_{z}P^{-j}(z)\right)\right\|\leq
4e^{2j\|A\|}K_{k}^{j}.$

This proves (2.21) for $|l|=1$. For $\ell>1$, by the identity

$\partial^{\ell}\left(\partial_{z}P^{-j}(z)\right)^{-1}=-\left(\partial_{z}P^{-j}(z)\right)^{-1}\sum\limits_{\begin{subarray}{c}\eta\prec\ell\\\
\eta\in\mathbb{Z}_{+}^{n},\eta\neq
0\end{subarray}}d_{\eta}\partial^{\eta}\left(\partial_{z}P^{-j}(z)\right)\partial^{\ell-\eta}\left(\partial_{z}P^{-j}(z)\right)^{-1},$

with $d_{\eta}$ nonnegative numbers, one has

(2.23)
$\displaystyle\|\partial^{\ell}\left(\partial_{z}P^{-j}(z)\right)^{-1}\|$
$\displaystyle\leq
2\mathfrak{d}_{k}e^{j\|A\|}\sum\limits_{\begin{subarray}{c}\eta\prec\ell\\\
\eta\in\mathbb{Z}_{+}^{n},\eta\neq
0\end{subarray}}\|\partial^{\eta}\left(\partial_{z}P^{-j}(z)\right)\|\|\partial^{\ell-\eta}\left(\partial_{z}P^{-j}(z)\right)^{-1}\|,$

where $\mathfrak{d}_{k}$ is the maximum of all $d_{\eta}$’s with
$\eta\prec\ell$ and $|\ell|\leq k$. For given $k$ the number of elements in
the summation of (2.23) is always bounded. Hence applying (2.16) and (2.21)
inductively to (2.23) yields that there exists a positive constant $M_{k}$
such that the estimation (2.21) holds.

Combining (2.20) and (2.21), and taking $N$ in (2.19) suitably large, one gets
from $(f_{1})$ that

$\begin{split}&\left\|\partial^{m}\left(\left(\partial_{z}P^{-j}(z)\right)^{-1}G\left(P^{-j}(z)\right)\right)\right\|\\\
&\leq\sum\limits_{\begin{subarray}{c}\ell\in\mathbb{Z}_{+}^{n}\\\ \ell\prec
m\end{subarray}}c_{\ell}\left\|\partial^{\ell}\left(\partial_{z}P^{-j}(z)\right)^{-1}\right\|\left\|\partial^{m-\ell}G\left(P^{-j}(z)\right)\right\|\\\
&\leq c_{k}s_{k}M_{k}^{j}a_{m-\ell}\kappa
e^{-jN\mu_{M}}\|z_{u}\|^{N}\leq\varrho^{j}\|z_{u}\|^{N},\end{split}$

where $s_{k}$ numerates the maximum number of elements in the summation for
all $m$ with $|m|\leq k$, and in the last inequality we have used the fact
that $\lim\limits_{N\rightarrow\infty}e^{-N\mu_{M}}=0$. Note that we can
choose $N$ large enough such that $\varrho<1$. This shows that for
$m\in\mathbb{Z}_{+}^{n}$ and $|m|=k$, the $m$th order derivative of the series
(2.12) is absolutely convergent. The claim is proved.

Claim 3 verifies that $h(z)$ is a $C^{k}$ function for arbitrary
$k\in\mathbb{N}$. Going back to $\mathcal{W}$ and $\mathcal{V}$ in the region
$\|z\|\leq 1/2$ via the linear change $\sigma z\rightarrow z$ as mentioned
previously, we get a $C^{k}$ equivalence between $\mathcal{W}$ and
$\mathcal{V}$ in a suitable neighborhood of the origin, whose size depends on
the small positive $\sigma$. Of course, with increase of $k\in\mathbb{N}$ the
region of the origin where the vector fields $\mathcal{W}$ and $\mathcal{V}$
are $C^{k}$ equivalent may decrease. By the arbitrariness of $k\in\mathbb{N}$
it follows that $h(z)$ is a $C^{\infty}$ function. Then we have a $C^{\infty}$
equivalence between $\widetilde{\mathcal{W}}$ and $\widetilde{\mathcal{V}}$.
And consequently, the vector fields $\mathcal{W}$ and $\mathcal{V}$ are
$C^{\infty}$ equivalent near the origin. As a consequence, system (1.1) has a
$C^{\infty}$ first integral in a neighborhood of the origin.

It completes the proof of the theorem. ∎

##  3\. Proofs of Theorem 1.2 and Corollary 1

According to the proof of the sufficiency of Theorem 1.1, system (1.1) is
analytically equivalent to system (2.3). So we only need to prove the theorem
for system (2.3). Let $\mathcal{X}_{0}$ be the linear vector field defined by
the linear part of system (2.3), and let $\mathcal{V}$ be the set of all
vector fields of form (2.3) with their linear part being the same as
$\mathcal{X}_{0}$.

######  Lemma 1.

Let $\mathcal{S}$ be an affine finite–dimensional subspace of $\mathcal{V}$,
and set

$\mathcal{S}:=\left\\{\mathcal{X}_{0}+t_{1}\mathcal{X}_{1}+\ldots+t_{m}\mathcal{X}_{m}|\
t=(t_{1},\ldots,t_{m})\in\mathbb{R}^{m}\right\\}.$

where $\mathcal{X}_{1},\ \ldots,\ \mathcal{X}_{m}$ are linearly independent
vector fields of the form

$(F_{1}(y)\mathbf{y}_{2},\ F_{2}(y)\mathbf{y}_{2})^{\tau},$

with $F_{1}$ and $F_{2}$ respectively $n-1$–dimensional analytic row vector
valued function and analytic matrix valued function of order $n-1$ satisfying
$F_{1}(0)=0$ and $F_{2}(0)=0$. Then each vector field in $\mathcal{S}$ has a
formal first integral in $(\mathbb{R}^{n},0)$, whose homogeneous component of
degree $k\in\mathbb{N}$ is a polynomial in $t$ of degree $k-1$ without
constant term.

In order to prove this lemma we need the next result, see e.g. Bibikov [2,
Lemma 1.1].

######  Proposition 2.

Let $C$ and $D$ be square matrices of order $n$ with respectively the
eigenvalues $\alpha=(\alpha_{1},\ldots,\alpha_{n})$ and
$\beta=(\beta_{1},\ldots,\beta_{n})$, and let $\mathcal{H}^{\ell}(x)$ be
$n$–dimensional vector valued homogeneous polynomials of degree $\ell$. Then
the linear operator

$\mathcal{L}(F)(x):=\frac{\partial F(x)}{\partial x}Cx-DF(x),\qquad
F(x)\in\mathcal{H}^{\ell}(x),$

has the set of eigenvalues

$\left\\{\langle m,\alpha\rangle-\beta_{j}|\ m\in\mathbb{Z}_{+}^{n},\ \
|m|=\ell,\ \ j=1,\ldots,n\right\\}.$

Proof of Lemma 1. By the construction of $\mathcal{S}$ it follows clearly that
$\mathcal{S}$ is isomorphic to $\mathbb{R}^{m}$. Denote by $\mathcal{X}_{t}$
the vector field in $\mathcal{S}$ for each $t\in\mathbb{R}^{m}$, and write its
associated differential system as

(3.1) $\begin{split}\dot{y}_{1}&=F_{1}(t,y)\mathbf{y}_{2},\\\
\dot{\mathbf{y}}_{2}&=B\mathbf{y}_{2}+F_{2}(t,y)\mathbf{y}_{2},\end{split}$

with $F_{1}$ and $F_{2}$ depending on $t$ linear, and $F_{1}$ a row vector
valued function and $F_{2}$ a square matrix valued function of order $n-1$.

Priori let $H(t,y)$ be an analytic or a formal first integral of system (3.1).
By the equivalent characterization of first integrals, it follows that
$H(t,y)$ satisfies the next equality.

(3.2) $F_{1}(t,y)\mathbf{y}_{2}\frac{\partial H}{\partial y_{1}}+\left\langle
B\mathbf{y}_{2}+F_{2}(t,y)\mathbf{y}_{2},\,\frac{\partial
H}{\partial\mathbf{y}_{2}}\right\rangle\equiv 0.$

Write $H(t,y)$ in series expansion in $y$

(3.3) $H(t,y)=\sum\limits_{s=\ell}\limits^{\infty}H_{s}(t,y),$

with $\ell\in\mathbb{N}$ and $H_{s}(t,y)$ a homogeneous polynomial in $y$ of
degree $s$ whose coefficients depending on $t$. Let the Taylor expansions of
$F_{1}(t,y)$ and $F_{2}(t,y)$ be

(3.4) $F_{1}(t,y)=\sum\limits_{s=1}\limits^{\infty}F_{1s}(t,y),\qquad
F_{2}(t,y)=\sum\limits_{s=1}\limits^{\infty}F_{2s}(t,y),$

with each $F_{1s}$ a vector valued homogeneous polynomial of degree $s$ in $y$
and each $F_{2s}(t,y)$ a square matrix valued homogeneous polynomials in $y$
of degree $s$. Here all $F_{1s}$’s and $F_{2s}$’s depend on $t$ linearly.

Let $\mathcal{L}$ be the linear operator defined by

$\mathcal{L}=\left\langle
B\mathbf{y}_{2},\,\frac{\partial}{\partial\mathbf{y}_{2}}\right\rangle.$

Substituting (3.3) and (3.4) into (3.2), and equating the homogeneous terms of
degree $\ell+s$ for $s=0,1,\ldots$ yields

(3.5) $\displaystyle\mathcal{L}(H_{\ell})=$ $\displaystyle 0,$
$\displaystyle\mathcal{L}(H_{\ell+s})=$
$\displaystyle-\sum\limits_{j=1}\limits^{s}\left\langle
F_{2j}(t,y)\mathbf{y}_{2},\frac{\partial
H_{\ell+s-j}}{\partial\mathbf{y}_{2}}\right\rangle$ (3.6)
$\displaystyle\qquad\
-\sum\limits_{j=1}\limits^{s}F_{1j}(t,y)\mathbf{y}_{2}\frac{\partial
H_{\ell+s-j}}{\partial y_{1}},\quad\ell=1,2,\ldots$

Let $\mathcal{H}_{s}(y)$ be the linear space formed by homogeneous polynomials
in $y$ of degree $s$. Since the $n-1$–tuple of eigenvalues
$\lambda_{2},\ldots,\lambda_{n}$ of $B$ are non–resonant, it follows from
Proposition 2 that equation (3.5) has only a solution in $y_{1}$. According to
Proposition 1 we can take $\ell=1$ and $H_{1}(y)=y_{1}$.

For $s=1$, equation (3.6) is simply

$\mathcal{L}(H_{2})=-F_{11}(t,y)\mathbf{y}_{2}.$

Since the right hand side of this last equation vanishes identically when
restricted to $\mathbf{y}_{2}=0$, applying Bibikov [2, Lemma 1.1] to this last
equation, one gets a unique homogeneous solution of degree $2$ modulus
$y_{1}^{2}$, which is linear in both $\mathbf{y}_{2}$ and
$t=(t_{1},\ldots,t_{m})$.

For $s=2$, equation (3.6) is reduced to

(3.7) $\mathcal{L}(H_{3})=-\left\langle
F_{21}(t,y)\mathbf{y}_{2},\frac{\partial
H_{2}}{\partial\mathbf{y}_{2}}\right\rangle-
F_{11}(t,y)\mathbf{y}_{2}\frac{\partial H_{2}}{\partial
y_{1}}-F_{12}(t,y)\mathbf{y}_{2}\frac{\partial H_{1}}{\partial y_{1}}.$

Note that the right hand side of this last equation is a homogeneous
polynomial in $y$ of degree $3$, and also a polynomial in $t$ of degree $2$
without constant term. Since the right hand side of (3.7) vanishes identically
when $\mathbf{y}_{2}=0$, applying again Proposition 2 to equation (3.7)
associated with $\mathbf{y}_{2}$, one gets a unique homogeneous solution in
$y$ of degree $3$ modulus $y_{1}^{3}$, which is quadratic in $t$.

Inductively, one gets that equation (3.6) for each $s=3,4,\ldots$ has a unique
homogeneous solution in $y$ of degree $1+s$ modulus $y_{1}^{1+s}$, which is of
degree $s$ in $t$.

This proves the existence of formal first integrals of each vector field
$\mathcal{X}_{t}$ in $\mathcal{S}$, with the properties stated in the lemma. ∎

For proving the generic divergence of the first integrals of system (2.3), we
will use the pluripolar set, see [9, 19] for more detail information.

By definition, a pluripolar set is a subset of $\mathbb{C}^{m}$, saying
$\Lambda$, which satisfies that for each $z\in\Lambda$, there exists a
neighborhood $U_{z}$ of $z$ and a plurisubharmonic function $p$ defined on
$U_{z}$ verifying that ${\Lambda}\cap U_{z}\subset p^{-1}(-\infty)$.

A plurisubharmonic function is the one: $p:\Omega\rightarrow[-\infty,\infty)$,
with $\Omega\subset\mathbb{C}^{m}$ an open subset, which

  * •

is upper semicontinuous, i.e. $\\{z\in\Omega:p(z)<c\\}$ open for each
$c\in\mathbb{R}$,

  * •

and is not identically $-\infty$ on any connected component of $\Omega$,

  * •

and satisfies that for any $z\in\Omega$

$p(z)\leq\frac{1}{2\pi}\int_{0}^{2\pi}p(z+e^{it}r)dt,$

with $r\in\mathbb{C}^{m}$ such that $z+\sigma r\in\Omega$, where
$\sigma\in\mathbb{C}$ satisfies $|\sigma|\leq 1$.

Recall that a pluripolar set in $\mathbb{C}^{m}$ has Lebesgue measure $0$, and
that the countable union of pluripolar subsets is also pluripolar.

Denote by $\mathcal{P}$ the set of plurisubharmonic functions in
$\mathbb{C}^{m}$ which have minimal growth, i.e. $u(z)-\log\|z\|$ is bounded
above when $\|z\|\rightarrow\infty$. For any subset
${\Lambda}\subset\mathbb{C}^{m}$, define

$V_{\Lambda}(z)=\sup\\{p(z)|\ \,p\in{\mathcal{P}},\ p|_{\Lambda}\leq
0\\},\qquad z\in\mathbb{C}^{m}.$

Having the above preparation one can recall the Bernstein–Walsh Lemma, see
e.g. [14] or [9, p.156].

Bernstein–Walsh Lemma. Assume that ${\Lambda}\subset\mathbb{C}^{m}$ is not
pluripolar, and that $P(z)$ is a polynomial of degree $k$. Then

$|P(z)|\leq\|P\|_{\Lambda}\exp\left(kV_{\Lambda}(z)\right),\quad\mbox{for }\ \
z\in\mathbb{C}^{m}$

where $\|P\|_{\Lambda}$ is the supremum of $P$ on $\Lambda$.

We now apply the Bernstein–Walsh lemma to prove the next result, which is the
key point in the proof of Theorem 1.2.

######  Lemma 2.

For the subset $\mathcal{S}$ given in Lemma 1, the following statements hold.

  * $(a)$

Either each vector field in $\mathcal{S}$ has an analytic first integral in
$(\mathbb{R}^{n},0)$,

  * $(b)$

or only the vector fields in an exceptional pluripolar subset of $\mathcal{S}$
have analytic first integrals in $(\mathbb{R}^{n},0)$.

###### Proof.

If statement $(b)$ holds, the proof is done. We now suppose that $\mathcal{S}$
has a subset, saying $\mathcal{T}$, which is not pluripolar, and whose each
element has an analytic first integral.

Since we are in the conditions of Lemma 1, we will use the notations and
conclusions provided in Lemma 1 and its proof. Recall from Lemma 1 that any
vector field in $\mathcal{S}$ has a formal first integral, whose homogeneous
part of degree $k$ is a polynomial in $t$ of degree $k-1$.

Let $\mathcal{B}$ be the subset of $\mathbb{C}^{m}$ such that each
$t\in\mathcal{B}$ corresponds to the vector field $\mathcal{X}_{t}$ in
$\mathcal{T}$. Set ${\mathcal{B}}=\cup_{r\geq 1}{\mathcal{B}}_{r}$, with
${\mathcal{B}}_{r}$ the subset of $\mathcal{B}$ such that for each
$t\in\mathcal{B}_{r}$ its associated vector field ${\mathcal{X}}_{t}$ has a
convergent first integral $H_{t}$ in the ball $D_{r}$ of radius $1/r$, which
is bounded by $1$ in $D_{r}$. Since $\mathcal{B}$ is not pluripolar and a
countable union of pluripolar sets is again pluripolar, there exists at least
a ${\mathcal{B}}_{r_{0}}$ which is non–pluripolar.

By the construction of $\mathcal{B}_{r_{0}}$ and the proof of Lemma 1, it
follows that for each $t\in\mathcal{B}_{r_{0}}$ the vector field
$\mathcal{X}_{t}$, and its associated system (3.1), has an analytic first
integral $H_{t}(y)$ in $D_{r_{0}}$, which is bounded by $1$. Write $H_{t}(y)$
in series in $y$ with coefficients in $t$ as

$H_{t}(y)=\sum\limits_{j\in{\mathbb{Z}}_{+}^{n}}H_{j}(t)y^{j},$

where $y^{j}$ is the multi–index as used previously, and $H_{j}(t)$’s are
polynomials of degree at most $|j|-1$ as obtained in the proof of Lemma 1. For
each $t\in\mathcal{B}_{r_{0}}$, since $H_{t}(y)$ is analytic in
${\mathcal{D}}_{r_{0}}$ and bounded by $1$, one gets from the Cauchy
inequality that there exists a $\rho_{0}>0$ for which

(3.8)
$R(t):=\sup_{j\in\mathbb{Z}_{+}^{n}}|H_{j}(t)|\rho_{0}^{|j|}<\infty,\qquad
t\in{\mathcal{B}}_{r_{0}}.$

In fact, $\rho_{0}$ can be taken to be $1/(r_{0}+1)$.

According to the values of $R(t)$ in $\mathcal{B}_{r_{0}}$, we further
decompose the non–pluripolar set ${\mathcal{B}}_{r_{0}}$ in the union of
$\mathcal{C}_{s}:=\\{t\in{\mathcal{B}}_{r_{0}}|\,\,R(t)\leq s\\},\ \
s\in\mathbb{N}$. Then there exists at least one non-pluripolar set among
$\mathcal{C}_{s}$’s, saying $\mathcal{C}_{s_{0}}$. Taking any fixed
non–pluripolar compact subset of $\mathcal{C}_{s_{0}}$, denoted by $D_{0}$.
Then $R(t)$ is bounded on $D_{0}$. Moreover, it follows from (3.8) that there
exists a $\rho_{1}>0$ such that

(3.9) $|H_{j}(t)|\leq\rho_{1}^{|j|},\quad\mbox{ for all }\ t\in D_{0},\ \
j\in\mathbb{Z}_{+}^{n}.$

Since $H_{j}(t)$ is a polynomial of degree $|j|-1$, we get from the
Bernstein–Walsh lemma together with (3.9) that for $t\in\mathbb{C}^{m}$

(3.10)
$|H_{j}(t)|\leq\|H_{j}\|_{D_{0}}\exp((|j|-1)V_{D_{0}}(t))\leq\rho_{1}^{|j|}\exp((|j|-1)V_{D_{0}}(t)),$

where $\|H_{j}\|_{D_{0}}$ represents the maximum value of $|H_{j}(t)|$ on
$D_{0}$.

Choose any fixed compact subset $\mathcal{E}\subset{\mathbb{C}}^{m}$. Then
there exists a $\rho_{2}>0$ such that

(3.11)
$\exp((|j|-1)V_{D_{0}}(t))\leq\exp((|j|-1)\|V_{D_{0}}\|_{\mathcal{E}})\leq\rho_{2}^{|j|},$

where $\|V_{D_{0}}\|_{\mathcal{E}}$ represents the maximum value of
$|V_{D_{0}}(t)|$ on $\mathcal{E}$. Note that for given $D_{0}$, $\rho_{2}$
depends only on $\mathcal{E}$. Combining the inequalities (3.10) and (3.11)
one obtains that

$\|H_{j}\|_{\mathcal{E}}\leq\rho_{1}^{|j|}\rho_{2}^{|j|}.$

Choose any $\rho\in(0,1/(\rho_{1}\rho_{2}))$. Then for any $t\in\mathcal{E}$
the formal first integral $H_{t}(y)$ of the vector field $\mathcal{X}_{t}$ is
convergent in the polydisc $\left\\{y=(y_{1},\ldots,y_{n})|\right.$
$\left.\,|y_{j}|\leq{\rho},i=1,\ldots,n\right\\}$.

By arbitrariness of the choice of the compact subset $\mathcal{E}$ in
$\mathbb{C}^{m}$, it follows that for any $t\in\mathbb{C}^{m}$ the vector
field $\mathcal{X}_{t}$ has an analytic first integral in a neighborhood of
the origin.

It completes the proof of the lemma. ∎

Proof of Theorem 1.2. $(a)$ By contrary we assume that there exists a
nonpluripolar subset of $\mathcal{K}$, whose any element has an analytic first
integral in a neighborhood of the origin. Then by Lemma 2 any system in
$\mathcal{K}$ has an analytic first integral in a neighborhood of the origin,
a contradiction with the assumption that $\mathcal{K}$ contains an element,
which has only formal first integrals in a neighborhood of the origin.

$(b)$ follows from the proof of statement $(a)$. This proves the theorem. ∎

Proof of Corollary 1. As mentioned in the remark after Corollary 1,
$\mathfrak{P}$ is a finite dimensional space. Let $m$ be the dimension of this
space. By the proofs of Theorem 1.1 and Lemma 1 it follows that any system in
$\mathfrak{P}$ has a formal first integral in a neighborhood of the origin,
whose each monomial has its coefficient being a polynomial in $t$. We can
write the nonlinear parts of the vector fields associated to the systems in
$\mathfrak{P}$ as $t_{1}\mathcal{N}_{1}+\ldots+t_{m}\mathcal{N}_{m}$, where
$t=(t_{1},\ldots,t_{m})\in\mathbb{R}^{m}$, and $\mathcal{N}_{j}$ for
$j\in\\{1,\ldots,m\\}$ is an $n$–dimensional vector field with one component a
nonlinear monomial and others all vanishing. Then the next proof follows from
the same arguments as those given in the proof of Lemma 1. The details are
omitted. ∎

## Acknowledgements

The author is partially supported by NNSF of China grant numbers 11671254,
11871334 and 12071284.

## References

  * [1] L. Barreira, J. Llibre and C. Valls, Integrability and zero-Hopf bifurcation in the Sprott A system, Bull. Sci. Math. 162 (2020), 102874, 16pp.
  * [2] Yu.N. Bibikov, Local Theory of Nonlinear Analytic Ordinary Differential Equations, Lecture Notes in Math., Vol. 702, Springer–Verlag, Berlin, 1979.
  * [3] A.V. Bolsinov and I.A. Taimanov, Integrable geodesic flows with positive topological entropy, Invent. Math. 140 (2000), 639–650.
  * [4] J. Chen, Y. Yi and X. Zhang, First integrals and normal forms for germs of analytic vector fields, J. Differential Equations 245 (2008), 1167–1184.
  * [5] W. Cong, J. Llibre and X. Zhang, Generalized rational first integrals of analytic differential systems, J. Differential Equations 251 (2011), 2770–2788.
  * [6] Z. Du, V.G. Romanovski and X. Zhang, Varieties and analytic normalizations of partially integrable systems, J. Differential Equations 260 (2016), 6855–6871.
  * [7] S.D. Furta, On non–integrability of general systems of differential equations, Z. Angew Math. Phys. 47 (1996), 112–131.
  * [8] Yu.S. Ilyashenko and S. Yakovenko, Lectures on Analytic Differential Equations, Amer. Math. Soc., Providence, RI, 2008.
  * [9] M. Klimek, Pluripotential Theory, London Math. Society Monographs 6, Oxford Univ. Press, New York, 1991.
  * [10] W.G. Li, Normal Form Theory and its Application (in Chinese), Science Press, Beijing, 2000.
  * [11] W. Li, J. Llibre, X. Zhang, Local first integrals of differential systems and diffeomorphisms, Z. Angew. Math. Phys. 54 (2003), 235–255.
  * [12] J. Llibre, C. Pantazi and S. Walcher, First integrals of local analytic differential systems, Bull. Sci. Math. 136 (2012), 342–359.
  * [13] J. Llibre and C. Valls, Centers of planar generalized Abel equations, J. Differential Equations 268 (2020), 6481–6487.
  * [14] R. Pérez-Marco, Convergence or generic divergence of the Birkhoff normal form, Ann. Math. 157 (2003), 557–574.
  * [15] H. Poincaré, Sur l’intégration des équations différentielles du premier order et du premier degré I and II, Rend. Circ. Mat. Palermo 5 (1891), 161–191; Rend. Circ. Mat. Palermo 11 (1897), 193–239.
  * [16] V.G. Romanovski and D.S. Shafer, The Center and Cyclicity Problems: A Computational Algebraic Approach, Birkhäuser, Boston, 2009.
  * [17] S. Shi, On the nonexistence of rational first integrals for nonlinear systems and semiquasihomogeneous systems, J. Math. Anal. Appl. 335 (2007), 125–134.
  * [18] S. Shi and Y. Li, Non–integrability for general nonlinear systems, Z. Angew. Math. Phys. 52 (2001), 191–200.
  * [19] M. Tsuji, Potential Theory in Modern Function Theory, Maruzen Co., Ltd, Tokyo, 1959.
  * [20] X. Zhang, The 16th Hilbert problem on algebraic limit cycles, J. Differential Equations 251 (2011), 1778–1789.
  * [21] X. Zhang, A note on local integrability of differential systems, J. Differential Equations 263 (2017), 7309–7321.
  * [22] X. Zhang, Integrability of Dynamical Systems: Algebra and Analysis, Springer, Singapore, 2017.
  * [23] X. Zhang, Analytic integrable systems: analytic normalization and embedding flows, J. Differential Equations 254 (2013), 3000–3022.
  * [24] N.T. Zung, Convergence versus integrability in Poincaré–Dulac normal form, Math. Res. Lett. 9 (2002), no. 2–3, 217–228.
  * [25] N.T. Zung, Convergence versus integrability in Birkhoff normal form, Ann. Math. 161 (2005), 141–156.


# Determinants of successful mitigation in coupled social-climate dynamics

Longmei Shu  Feng Fu

###### Abstract

Understanding the impact of human behavior is crucial for successful
mitigation of climate change across the globe. To shed light onto this issue,
here we couple the forest dieback model with human behaviors. Using
evolutionary game theory, we build a time-delay system where forest growth is
impacted by both temperature and human mitigation choices, the latter being
informed by temperature forecasts. Simulations of the coupled system over 200
years show us the varying outcomes: forest dies out and no one is a mitigator,
forest dies out and everyone is a mitigator, or the forest survives and
everyone is a mitigator. There exist rare cases where no one is a mitigator
and yet the forest survives, but with a low coverage. We also find occasional
oscillations where the proportion of mitigators vary between 0 and 1. Our
results are based on simple models but have profound insights into
determinants of behavior changes desired in social-climate dynamics.

###### keywords:

social-climate models , evolutionary game theory , climate change

††journal: Proceedings of the Royal Society A

[inst1]organization=Department of Mathematics, addressline=27 N. Main Street,
city=Hanover, postcode=03755, state=NH, country=USA

[inst2]organization=Department of Biomedical Data Science, addressline=Geisel
School of Medicine at Dartmouth, city=Lebanon, postcode=03756, state=NH,
country=USA

##  1 Introduction

Climate change, or global warming is one of many public goods games where
social cooperation can play a very important role in averting potential
catastrophic scenarios [1, 2, 3, 4]. Typically such collective-risk social
dilemmas have complicated dynamics [5, 6, 7, 8]. Experiments with real people
and money reward were also carried out to see how people would act in
threshold climate games [9, 10]. Field experiments in multiple countries
around the world [11] give us some insights on the climate problem as well.

Models increasingly help us understand the interactions between the carbon
cycle, the climate system, human processes, and the impact of policies [12].
To ensure those policy decisions are robust to uncertainties, multiple
scenarios are often laid out, ranging from carbon emission trajectories to
socio-economic systems pathways. Influence in these models usually flows in
one-direction, from socio-economic systems to the Earth system. Yet, two-way
feedback mechanisms link climate and social processes: human behavior changes
the climate, and the climate changes human opinions and consequently human
behavior. Coupled human-environment models are already widely applied to study
other systems such as fisheries and forests and the need for coupled social-
climate models has been noted [13, 14, 15]. Such eco-evolutionary models where
the environment changes depending on the actions of the players again adds to
the rich dynamics of the system [16, 17].

Many social-climate models are complicated and we don’t have a good intuition
of the dynamics over time. Here we focus on a simple climate model, the forest
dieback model, and couple it with a replicator equation that determines the
percentage of mitigators in the population based on the warming over the past
10 years. So the social part of our model has a time-delay effect, and the
growth rate of the forest depends on the percentage of mitigators in the
population. We study the behavior of this coupled social-climate model under
changing parameter values, the ambient temperature (or background temperature
of the environment), and a sensitivity to warming. The sensitivity to warming
decides how likely people are going to adopt mitigative behavior based on
current warming. We solve our delay system over 200 years and see if the
forest dies out or survives, and how much of the population mitigates climate
change.

Since our model is relatively simple, we are able to analyze the equilibrium
points well and see how they change when the parameters change. We run
simulations over a wide range of parameters to see all possible rich dynamics
of our model. Our findings, though derived from simple models, offer deep
understanding into the factors influencing mitigation choices in the context
of social-climate interactions.

##  2 Methods and Model

###  2.1 Forest dieback model

We consider the climate and temperature for Amazonian forests. The vegetation
coverage $v$ is between 0 and 1, with 0 for bare soil and 1 for full forest
coverage. Then $v$ satisfies the following equation

$\displaystyle\frac{dv}{dt}=gv(1-v)-\gamma v.$ (1)

Here $\gamma=0.2$ is the disturbance rate [18] and $g$ is the growth rate. The
equilibrium points of the system have to satisfy

$\displaystyle gv(1-v)=\gamma v,\quad v=0\text{ or }v=1-\frac{\gamma}{g}.$

The fixed point under varying values of $g$ are plotted in Figure 1. We can
see the vegetation coverage has a positive robust stable fixed point when
$g>2$.

(a) $g=0.3$

(b) $g=1.2$

(c) $v=1-\gamma/g$

Figure 1: Equilibrium under varying growth rate $g$

In the forest dieback model $g$ is given by

$\displaystyle
g=g_{0}\left[1-\left(\frac{T-T_{opt}}{\beta}\right)^{2}\right].$ (2)

Here $g_{0}=2$ is the maximum growth rate, and $\beta=10$ is the half-width of
the growth versus temperature curve. $T_{opt}=28^{\circ}$C represents the
optimal temperature for plant growth.

$T$ is the actual temperature and is given by

$\displaystyle T=T_{v}+(1-v)a.$ (3)

Here $a=5$ is the difference between surface temperature of bare soil and
forest and $T_{v}$ is the temperature with full forest coverage. For
convenience, we will call it the ambient temperature from here on.

Plug (2) and (3) into (1), we get fixed points of the forest dieback model
satisfies

$\displaystyle
v\left\\{2(1-v)\left[1-\left(\frac{T_{v}+a(1-v)-T_{opt}}{\beta}\right)^{2}\right]-\gamma\right\\}=0.$
(4)

Plug in all the parameter values, we get

$\displaystyle v=0\text{ or }(1-v)[1-0.01(T_{v}-23-5v)^{2}]=0.1.$

Figure 2: Equilibrium under varying ambient temperature

Equilibrium points for varying values of $T_{v}$ are shown in Figure 2. We can
see bifurcation [19, 5, 20] happens and we have two different fixed points
when the ambient temperature $T_{v}$ is between 32.5∘C and 34.7∘C. The green
dashed line represents the unstable fixed points $v=0$, the blue curve
represents the stable and robust fixed points, and the red dashed curve
represents fixed points in the shaded bifurcation region.

###  2.2 Human actions

We divide the population into two groups, mitigators (M) and nonmitigators
(N). And the percentage of mitigators is $x$, the percentage of nonmitigators
is $y$. So $0<x,y<1,x+y=1$. The payoffs of mitigators and nonmitigators are

$\displaystyle E_{M}$ $\displaystyle=-\alpha+\frac{1}{2}f(T_{f})+\delta x,$
$\displaystyle E_{N}$ $\displaystyle=-\frac{1}{2}f(T_{f})+\delta y.$

Here $\alpha=1$ is the cost of mitigation, and $\delta=1$ is the strength of
social norms. [21]

Let $\bar{E}$ be the average payoff in the population, then the proportion of
mitigators satisfies the following replicator equation.

$\displaystyle\frac{dx}{dt}$ $\displaystyle=x(E_{M}-\bar{E})$
$\displaystyle=x(E_{M}-xE_{m}-yE_{N})$ $\displaystyle=x(1-x)(E_{M}-E_{N})$
$\displaystyle=x(1-x)[-\alpha+f(T_{f})+\delta(2x-1)]$

So far we haven’t explained what $f(T_{f})$ is yet. This term corresponds to
the cost of global warming and satisfies the following equation.

$\displaystyle f(T_{f})=\frac{f_{max}}{1+e^{-w(T_{f}-T_{c})}}$ (5)

Here $f_{max}=5$ is the maximum warming cost, $w=3$ is the non-linearity of
warming cost, and $T_{c}$ is the critical temperature [21]. $T_{f}$ is the
perceived temperature rise given by

$\displaystyle T_{f}(t)=\frac{t_{f}}{t_{p}}[T(t)-T(t-t_{p})]$ (6)

Here $t_{p}=10$ is the number of previous years used for temperature
projection, $t_{f}=15$ is the number of years ahead for temperature
projection.

Plug (5),(6) in the replicator equation, we get

$\displaystyle\frac{dx}{dt}$
$\displaystyle=x(1-x)\left[\delta(2x-1)-\alpha+\frac{f_{max}}{1+e^{-w[t_{f}/t_{p}T(t)-t_{f}/t_{p}T(t-t_{p})-T_{c}]}}\right].$
(7)

###  2.3 Coupled social-climate model

In the previous section, we already considered how the temperature affects
human actions. Now we will try to consider the effects of human actions on the
climate as well. We will modify the growth rate given by (2) by adding a term
that depends on $x$, the proportion of mitigators,

$\displaystyle
g=g_{0}\left[1-\left(\frac{T-T_{opt}}{\beta}\right)^{2}\right]\eta(x).$

Obviously we want $\eta(x)$ to increase when $x$ increases, since the more
mitigators there are, the better the forest grows. To make sure that the
change in the proportion of mitigators has an observable effect in the forest
coverage, we set $\eta(x)=0.2+0.4x$ so that when $x$ changes from 0 to 1,
roughly $g$ changes between 0.4 and 1.2, where the equilibrium forest coverage
is sensitive to the growth rate, see Figure 1.

Plug in the values for all parameters other than $T_{v}$ and $T_{c}$, the two
equations we are coupling are

$\begin{cases}\displaystyle\frac{dv}{dt}&=2[1-0.01(T_{v}-23-5v)^{2}](0.2+0.4x)v(1-v)-0.2v,\\\
\vspace{5pt}\\\
\displaystyle\frac{dx}{dt}&=x(1-x)\left[2x-2+\dfrac{5}{1+e^{-3[-7.5v(t)+7.5v(t-10)-T_{c}]}}\right].\end{cases}$
(8)

This is a system of autonomous nonlinear delay differential equations.

##  3 Theoretical Analysis

We now turn to analyze the coupled social-climate model in Eq. 8. Let’s look
at the fixed points of this system. For $dv/dt=0$ we get either $v=0$ or

$\displaystyle h(v,x,T_{v})=2(1-v)(0.2+0.4x)[1-0.01(T_{v}-23-5v)^{2}]-0.2=0.$
(9)

(a) $x=0$

(b) $x=1$

Figure 3: Equilibrium under varying ambient temperature and mitigation

$h(v,x,T_{v})$ depends on the value of $x$ smoothly. As we vary the value of
$x$ through $[0,1]$, we expect the curve defined by (9) for $T_{v}$ and $v$ to
change smoothly. We show the curves for $x=0$ and $x=1$ in Figure 3. As one
can see, the solution curve is well defined for $v$ in $(0,1)$ for a given
value of $T_{v}$ up to a certain point, beyond which bifurcation happens where
one ambient temperature $T_{v}$ corresponds to two different values for $v$
and a small change in $T_{v}$ can result in a dramatic change in $v$ value.
For $x=0$, bifurcation happens when $T_{v}$ is between 30∘C and 30.2∘C. For
$x=1$, bifurcation happens when $T_{v}$ is between 32.1∘C and 33.7∘C.

One can check that above the blue curve $h(v,x,T_{v})<0$ and below the blue
curve $h(v,x,T_{v})>0$, so $v=0$ is unstable and the fixed point determined by
the blue curve is stable and robust.

To have a better understanding of our equation, we also show plots for the
equilibrium points under different fixed ambient temperatures. We still expect
the relationship curve to change smoothly as the ambient temperature changes.
For ambient temperatures between 18∘C and 30∘C, the curve looks similar to the
case where $T_{v}=25^{\circ}$C in figure 4. Below this curve,
$h(v,x,T_{v})>0$, and above the curve $h(v,x,T_{v})<0$. As the ambient
temperature lowers, this curve shifts down until it disappears from our phase
space around 14∘C. However, as the ambient temperature rises, the curve
doesn’t rise up and disappear. Instead we see a horizontal parabola-like curve
that moves up and right until it disappears from the right side of our phase
space around 34∘C. To show this process, we plot the curve under ambient
temperatures 30∘C, 31∘C, 32∘C, and 33∘C in Figure 4.

As one can see, when the temperature is higher than 30∘C, inside of the
parabola-like curve, $dv/dt>0$, and outside of the parabola-like curve
$dv/dt<0$. This gives us a hysteresis loop and backward bifurcation. The blue
curve represents stable fixed points and the red dashed curve represents
unstable fixed points even though they are technically the same curve, i.e.
they are determined by the same equation. Below a certain vegetation coverage,
the fixed point is unstable and above a certain vegetation coverage, the fixed
point is stable.

Notice that the green straight line $v=0$ sometimes is dashed and sometimes is
solid, the solid parts correspond to stable fixed points and dashed parts
correspond to unstable fixed points.

(a) $T_{v}$=15∘C

(b) $T_{v}$=25∘C

(c) $T_{v}$=30∘C

(d) $T_{v}$=31∘C

(e) $T_{v}$=32∘C

(f) $T_{v}$=33∘C

Figure 4: Equilibrium under varying mitigation and ambient temperature

For $dx/dt=0$ we get $x=0$, $x=1$, or

$\displaystyle x^{*}(T_{c})=1-\frac{2.5}{1+e^{3T_{c}}}.$ (10)

Figure 5: Equilibrium under varying critical warming

For $x^{*}$ to be between 0 and 1, we need $T_{c}>1/3\ln(1.5)\approx
0.135^{\circ}C$. And for the default value listed in [18], where
$T_{c}=2.5^{\circ}C$, we have $x^{*}(2.5)\approx 0.9986$, very close to $1$.
Among the three fixed point values, $x=0$ and $x=1$ are stable and
$x^{*}(T_{c})$ is unstable. We show a graph of (10) in Figure 5. As one can
see $x^{*}$ approaches 1 as soon as the critical warming reaches around 2∘C.

We give a qualitative phase portrait of the system when (9) has a well defined
solution for $v$, which we will call $v^{*}(x,T_{v})$. In Figure 6 (a), we
have 6 fixed points and $v^{*}$ is on a red curve that moves smoothly when
$T_{v}$ changes. When $T_{v}$ decreases, the red curve moves down. For
example, when $T_{v}=15^{\circ}$C, the left endpoint of the curve will move
below 0, leaving us 5 fixed points instead of 6, see Figure 6 (b). As the
curve moves further down we will have 4 fixed points, when the curve only
intersects $x=1$. When $T_{v}\leq 13^{\circ}$C, the curve will move out of the
phase space, leaving us 3 fixed points, $(0,0),(1,0)$ and $(x^{*},0)$.

On the other hand when $T_{v}$ increases, instead of moving up, the red curve
moves up and to the right showing a horizontal parabola that opens to the
right. We will run into a 5-fixed-point case similar to phase portrait 6 (b).
When $T_{v}=33^{\circ}$C, we see a horizontal parabola which gives us 7 fixed
points, see Figure 6 (c). When $T_{v}\geq 34^{\circ}$C the red curve moves out
of our $[0,1]\times[0,1]$ phase space, leaving us with the 3 fixed points on
the horizontal axis, $(0,0),(1,0)$ and $(x^{*},0)$.

(a) 6 fixed points

(b) 5 fixed points

(c) 7 fixed points

Figure 6: Qualitative phase portrait

In phase portrait 6 (a), we can see that only two of the fixed points are
stable, namely $(0,v^{*}(0,T_{v}))$ and $(1,v^{*}(1,T_{v}))$. All the other
fixed points are unstable. When $T_{v}$ decreases and the red curve moves
down, fixed point $(0,v^{*}(0,T_{v}))$ disappears and $(0,0)$ becomes stable.
This stays true until $(1,v^{*}(1,T_{v}))$ disappears and $(1,0)$ becomes
stable.

Now when $T_{v}$ increases, the red curve becomes a parabola, starting
somewhere between 0 and $x^{*}$ on the $x-$axis. It goes up and then right,
still giving us 5 fixed points. Although the shape of the curve is different
from phase portrait 6 (b), the stability analysis stays the same, so $(0,0)$
replaces $(0,v^{*}(0,T_{v}))$ as the stable fixed point. As the parabola moves
to the right, we will have 6 fixed points again, where the parabola starts
from somewhere between $x^{*}$ and 1 on the $x-$axis and moves left and up,
then right and up. Although the number of fixed points changes and the phase
portrait looks different, the stable fixed points are the same as before,
$(0,0)$ and $(1,v^{*}(1,T_{v}))$. As $T_{v}$ increases further, we will reach
phase portrait 6 (c) with 7 fixed points, and we will have one extra stable
fixed point $(1,0)$. As the parabola moves further right, it will only
intersect $x=1$, we have 5 fixed points and the same three stable ones,
$(0,0),(1,0)$ and $(1,v^{*}(1,T_{v}))$. Finally when the red curve disappears
on the right, we only have the three fixed points on the $x-$axis left and
$(0,0)$ and $(1,0)$ are both stable.

Also as we pointed out earlier, when $T_{c}\geq 2^{\circ}$C, $x^{*}$ is very
close to 1. This means the right strip in the phase portraits in Figure 6 can
be very thin.

##  4 Simulation results

We use the original uncoupled forest dieback model (1) to compute the
temperature of the forest for the first 10 years. Then we use a delay
differential equation solver to find the evolution of the coupled social-
climate model (8) over years 10 to 200. We repeat the simulation with
different initial conditions $x(0),v(0)$ and different $T_{v}$ and $T_{c}$
values.

We already discussed the stability of the qualitative phase portraits of our
system and the stable fixed points were $(0,0),(1,0),(0,v^{*}(0,T_{v}))$ and
$(1,v^{*}(1,T_{v}))$. This tells us the proportion of mitigators should
converge to either $1$ or $0$, which is what we see in our simulation results.
In the simulation results, the vegetation coverage either stabilizes at a
positive value or converges to zero, with some rare cases where it’s
oscillating over a long time. We should point out that some oscillations die
out in 50 years or so and converge, so it’s possible the oscillations that
persist by year 200 could

Figure 7: Forest coverage and mitigation after 200 years

die out if we computed for a longer time period as well. We only computed the
model up to year 200 since it might not be realistic to compute it for longer
considering the ambient temperature probably would have changed.

We call the cases where both $x$ and $v$ converge to 0 by year 200 ”Barren”
and the cases where $x$ converges to 1 and $v$ converges to 0 ”Brown”. ”Green”
represents cases where $x$ converges to 1 and $v$ converges to some positive
value typically around 0.7 or 0.8. We will call the long term oscillations
”Swing”. Our simulation results are shown in Figure 7.

Here the dark blue or purple region, the ”Barren”, means the forest dies out
and nobody is a mitigator. The blue region, the ”Brown”, represents situations
where the forest dies even though everyone is a mitigator. The green regions,
the ”Green”, means the forest survives and everyone is a mitigator. The yellow
region, the ”Swing”, represents cases where the vegetation and mitigation
oscillate and don’t converge by year 200.

In general, we can view $T_{c}$ as the sensitivity level of the population
towards warming. As we can see, typically, there is more ”Green” for low
values of $T_{c}$, which means people will adopt mitigative actions as soon as
they feel a small amount of warming. As $T_{c}$ increases, people become less
sensitive towards warming and don’t adopt mitigation until the warming is
rather severe. This might be too late as we can see a lot more ”Barren” at
high values of $T_{c}$.

The ambient temperature $T_{v}$ plays a role as well. When $T_{v}$ is either
too high or too low, we also see more ”Barren” and ”Brown”. This agrees with
the general consensus that there is a critical warming that we cannot go back
from. And of course if the ambient temperature is too low we probably can’t
survive either. It is worth noting that between around 20∘C and 30∘C, we see a
lighter violet region, which represents no mitigators and the vegetation
coverage converges to a low positive value, around 20% for example. This
happens when we have favorable environment where the forest stays even though
there are no mitigators.

The initial conditions also change the final result. For example, when we
start from a high forest coverage of 90% and low mitigation at 10%, we have a
good amount of ”Green” and small amounts of ”Barren”. If the mitigation we
start with is 90%, the ”Barren” region above 30∘C actually becomes bigger,
this could be because we don’t have a lot of room for improvement in behavior
since the mitigation we started with is already high. And also if we start
with low forest coverage at 10%, we will have huge ”Barren” regions.

We can also see that oscillations happen typically when we cross over from one
region to the other. Our model is nonlinear and has time-delay, it

Figure 8: Long term oscillations

can give us very rich dynamics. Some oscillations die away within the first 50
years while some persist up to year 200. Some oscillations are small in
amplitude and look like small persistent noise while others have value
variation from 0 to 1, the full range of possible values. We show two specific
plots for the evolution of vegetation and mitigation over 200 years for the
initial condition $x(0)=0.9,v(0)=0.1$ in Figure 8. As one can see, the
mitigation varies from 0 to 1, the vegetation varies between 0.4 and 0.8, and
the period can vary greatly as well.

We show two examples of evolution over time for ”Barren” in Figure 9 and two
examples of ”Brown” in Figure 10. Figure 11 shows two cases of ”Green”
evolution and Figure 12 shows the cases where the ambient temperature is
optimal and the forest survives even with no mitigators in the population at
all.

Figure 9: Two different instances of ”Barren” evolution

Figure 10: Two different instances of ”Brown” evolution

Figure 11: Two different instances of ”Green” evolution

Figure 12: Two different instances of optimal ambient temperature evolution

##  5 Discussion

Evolutionary game theory is used in this study to build a time-delay system
where the temperature of the forest depends on human behaviors and human
behaviors also depend on the current temperature. The replicator equation
determines the percentage of mitigators in the population based on the warming
over the past 10 years. This means that the social part of the model has a
time-delay effect, and the growth rate of the forest depends on the percentage
of mitigators in the population. By studying this coupled social-climate model
under changing parameter values, ambient temperature, or background
temperature of the environment, and sensitivity to warming, we can gain
insights into how human behavior affects climate change and how we can
mitigate it.

Simulations of the coupled system over 200 years show us the varying results
of forest dies out and no one is a mitigator, forest dies out and everyone is
a mitigator, or the forest survives and everyone is a mitigator. There are
rare cases where no one is a mitigator and yet the forest survives, but with
low coverage. There are also rare oscillations where the proportion of
mitigators varies between 0 and 1. These scenarios tell us that human behavior
plays a crucial role in the survival of forests and that having everyone as a
mitigator can help prevent forest dieback.

Our model was relatively simple and didn’t consider the asymmetry of the
population in aspects such as the contribution to global warming and resource
inequalities [21, 22, 23]. We could modify both the the climate part and the
social part of our model to see how the asymmetries in the population affects
the climate change. Our model has a deterministic process but there are
stochastic game descriptions of the climate change as well [24].

In sum, our work presents a transparent approach to modeling and analyzing the
intricate interplay within social-climate systems. At its core, the climate
model, grounded in forest dieback principles, may be straightforward, but it
effectively incorporate the foundational effects of human behavior on both the
climate and forest coverage. The decision-making process related to mitigation
is firmly grounded on temperate forecasting and the sensitivity to temperature
change. Our model, while simple, exhibits rich dynamics. This resonates with
Bob May’s renowned work on logistic maps and chaos: “Simple mathematical
models with very complicated dynamics” [19]. Our current model and its
extensions will help offer insights into the factors driving behavioral shifts
in the context of climate dynamics.

## References

  * [1] L. Kemp, C. Xu, J. Depledge, K. L. Ebi, G. Gibbins, T. A. Kohler, J. Rockström, M. Scheffer, H. J. Schellnhuber, W. Steffen, et al., Climate endgame: Exploring catastrophic climate change scenarios, Proceedings of the National Academy of Sciences 119 (34) (2022) e2108146119.
  * [2] F. P. Santos, J. M. Pacheco, F. C. Santos, S. A. Levin, Dynamics of informal risk sharing in collective index insurance, Nature Sustainability 4 (5) (2021) 426-432.
  * [3] S. A. Levin, J. M. Anderies, N. Adger, S. Barrett, E. M. Bennett, J. C. Cardenas, S. R. Carpenter, A.-S. Crépin, P. Ehrlich, J. Fischer, et al., Governance in the face of extreme events: lessons from evolutionary processes for structuring interventions, and the need to go beyond, Ecosystems (2021) 1-15.
  * [4] A. Kruczkiewicz, J. Klopp, J. Fisher, S. Mason, S. McClain, N. Sheekh, R. Moss, R. Parks, C. Braneon, Compound risks and complex emergencies require new approaches to preparedness, Proceedings of the National Academy of Sciences 118 (19) (2021) e2106795118.
  * [5] J. Wang, F. Fu, T. Wu, L. Wang, Emergence of social cooperation in threshold public goods games with collective risk, Physical Review E 80 (1) (2009) 016101.
  * [6] Y. Zhang, F. Fu, T. Wu, G. Xie, L. Wang, A tale of two contribution mechanisms for nonlinear public goods, Scientific reports 3 (1) (2013) 2021.
  * [7] T. Wu, F. Fu, Y. Zhang, L. Wang, The increased risk of joint venture promotes social cooperation, PloS one 8 (6) (2013) e63801.
  * [8] A. R. Tilman, V. V. Vasconcelos, E. Akçay, J. B. Plotkin, The evolution of forecasting for decision making in dynamic environments, arXiv preprint arXiv:2108.00047 (2021).
  * [9] M. Milinski, R. D. Sommerfeld, H.-J. Krambeck, F. A. Reed, J. Marotzke, The collective-risk social dilemma and the prevention of simulated dangerous climate change, Proceedings of the National Academy of Sciences 105 (7) (2008) 2291–2294.
  * [10] J. Jacquet, K. Hagel, C. Hauert, J. Marotzke, T. Röhl, M. Milinski, Intra-and intergenerational discounting in the climate game, Nature climate change 3 (12) (2013) 1025-1028.
  * [11] J.-C. Cárdenas, M. A. Janssen, M. Ale, R. Bastakoti, A. Bernal, J. Chalermphol, Y. Gong, H. Shin, G. Shivakoti, Y. Wang, et al., Fragility of the provision of local public goods to private and collective risks, Proceedings of the National Academy of Sciences 114 (5) (2017) 921–925.
  * [12] F. C. Moore, K. Lacasse, K. J. Mach, Y. A. Shin, L. J. Gross, B. Beckage, Determinants of emissions pathways in the coupled climate-social system, Nature 603 (2022) 103–111.
  * [13] S. Laborde, A. Fernandez, S. C. Phang, I. M. Hamilton, N. Henry, H. C. Jung, A. Mahamat, M. Ahmadou, B. K. Labara, S. Kari, M. Durand, B. Mark, P. Scholte, N. Xiao, R. Zieb, M. Moritz, Social-ecological feedbacks lead to unsustainable lock-in in an inland fishery, Global Environmental Change 41 (2016) 13–25.
  * [14] B. Beckage, K. Lacasse, J. M. Winter, L. J. Gross, N. Fefferman, F. M. Hoffman, S. S. Metcalf, T. Franck, E. Carr, A. Zia, A. Kinzig, The earth has humans, so why don’t our climate models?, Climate Change 163 (2020) 181–188.
  * [15] B. Beckage, F. C. Moore, K. Lacasse, Incorporating human behaviour into earth system modelling, Nature Human Behaviour (2022) 1–10.
  * [16] X. Wang, F. Fu, Eco-evolutionary dynamics with environmental feedback: Cooperation in a changing world, Europhysics Letters 132 (1) (2020) 10001.
  * [17] L. Shu, F. Fu, Eco-evolutionary dynamics of bimatrix games, Proceedings of the Royal Society A 478 (2267) (2022) 20220567.
  * [18] P. D. Ritchie, J. J. Clarke, P. M. Cox, C. Huntingford, Overshooting tipping point thresholds in a changing climate, Nature 592 (2021) 517–523.
  * [19] R. M. May, Simple mathematical models with very complicated dynamics, Nature 261 (5560) (1976) 459–467.
  * [20] L. Liu, X. Chen, A. Szolnoki, Coevolutionary dynamics via adaptive feedback in collective-risk social dilemma game, Elife 12 (2023) e82954.
  * [21] J. Menard, T. M. Bury, C. T. Bauch, M. Anand, When conflicts get heated, so does the planet: coupled social-climate dynamics under inequality, Proceedings of the royal society B 288 (2021) 1959.
  * [22] J. Wang, F. Fu, L. Wang, Effects of heterogeneous wealth distribution on public cooperation with collective risk, Physical Review E 82 (1) (2010) 016102.
  * [23] V. V. Vasconcelos, F. C. Santos, J. M. Pacheco, S. A. Levin, Climate policies under wealth inequality, Proceedings of the National Academy of Sciences 111 (6) (2014) 2212–2216.
  * [24] W. Barfuss, J. F. Donges, V. V. Vasconcelos, J. Kurths, S. A. Levin, Caring for the future can turn tragedy into comedy for long-term collective action under risk of collapse, Proceedings of the National Academy of Sciences 117 (23) (2020) 12915–12922.