CENTERS AND LIMIT CYCLES OF VECTOR FIELDS
DEFINED ON INVARIANT SPHERES

CLAUDIO A. BUZZI, ANA LIVIA RODERO, AND JOAN TORREGROSA∗

Abstract. The aim of this paper is the study of the center-focus and cyclicity
problems inside the class X of 3-dimensional vector fields that admit a first integral
that leaves invariant any sphere centered at the origin. We classify the centers
of linear, quadratic homogeneous and a family of quadratic vector fields F ⊂ X,
restricted to one of these spheres. Moreover, we show the existence of at least 4
limit cycles in family F .

1. Introduction

Differential equations and dynamical systems appear naturally in the description
of many phenomena for which local processes are known. The central problem is
then to obtain global information on these phenomena. Once the local equations
are formulated in a particular context, the next usual step is to solve them. But,
as in general, the evolution of these process is governed by nonlinear differential
equations, it is not always simple to solve them. The basic idea behind the first
works in the 18th and 19th centuries was to seek solutions that are combinations
of known functions. That is why it is imperative to search for new more geometric
methods for a better understanding of the behavior of the solutions of a system of
differential equations. Integrability is one of them.

The integrability is an intrinsic property of a given system that imposes strong
constraints on the way solutions evolve in phase space. The notion of integrability
was introduced to describe the property of equations for which all local and global
information can be obtained either explicitly from the solutions or implicitly from
invariants. The first class of invariants is the constants of motion, conserved quan-
tities, or first integrals. Of course, there are also other invariants such as integral
invariants, integrating factors, Jacobi multipliers, or symmetries which give rise to
different techniques for integrating differential equations, see for instance [1, 4, 12, 19]
and references therein. We have been motivated to consider the existence of first
integrals.

The importance of the existence of a non-constant first integral lies in the fact
that the trajectories of the vector field remain in the level sets of the function that
defines the first integral, and hence this is a strong constraint on the dynamical
behavior. In the theory of ordinary differential equations, the existence of first
integrals is important not only because they allow decreasing the dimension where
the differential system is defined but also because they simplify the characterization
of the phase portrait. It is important to mention that if we are working in a space

2020 Mathematics Subject Classification. Primary 34C07; Secondary 34C23, 37C27.
Key words and phrases. Vector fields on invariant spheres, integrability, center-focus problem,

local cyclicity.
∗ Corresponding author.

1



2 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

of dimension n and the system of equations has n − 1 independent first integrals
then we say that the system is completely integrable. The complete integrability
means that we can obtain the trajectories just intersecting the level sets of the first
integrals.

In this work we consider the class X of vector fields on R3 that admit the first
integral H(x, y, z) = x2 + y2 + z2, and we denote by Xn when we restrict to the
polynomial of degree n class. It means that 〈x, X(x)〉 = 0, for all x ∈ R3, and
that any sphere centered at the origin is invariant by the flow of X ∈ X. We are
interested on the center focus and cyclicity problems inside this class of systems.
These problems are strongly related to the study and definition of limit cycle (i.e.
isolated periodic orbit) for planar polynomials vector fields. This concept is due
to H. Poincaré, in 1880s. At the end of the XIX century, D Hilbert presented a
list of 23 problems at the International Congress of Mathematicians, in Paris. The
question on the second part of the sixteenth problem about an estimation of the
maximal number and relative positions of the limit cycles of a planar polynomial
vector field remains unsolved. For more details on the history of this important and
open problem, we refer to the interesting survey written by Y. Ilyashenko in [15].

Throughout the work, our strategy is to consider the restriction of the system of
differential equations to a 2-dimensional sphere and to use a stereographic projection
to consider a planar vector field. From there we can use all the tools that are
normally used to study the dynamics of planar differential systems.

Let X ∈ X1 be a linear vector field, we will see that it is always homogeneous and
it writes in the form

ẋ = −a1y − a2z,
ẏ = a1x− a3z,
ż = a2x+ a3y.

(1)

The following result provides a qualitative classification of the equilibrium points
of the above differential system.

Theorem 1. Let p ∈ S2
ρ = {(x, y, z) : x2 + y2 + z2 = ρ2} be an equilibrium point of

system (1) which is isolated on S2
ρ. Then p is of center type. Moreover, the system

is completely integrable.

As we have commented above, if X ∈ X1 is linear then X is homogeneous and
hence XH

1 = X1, where the superscript denotes the homogeneous property and the
subscript the degree of the vector field. Inspired by this fact, we study and classify
the center equilibrium points in the class of homogeneous quadratic vector fields,
that is in XH

2 .
In Section 3 we prove that, without loss of generality, a vector field X ∈ XH

2 writes
in the following canonical form

ẋ = −a4xy − a5xz − (a6 + a7)yz − a8z2,
ẏ = a4x

2 + a6xz − a9z2,
ż = a5x

2 + a7xy + a8xz + a9yz.

(2)

We notice that the equilibrium point is located at (0, 1, 0). On the following result,
we classify the equilibrium points of center type.

Theorem 2. The equilibrium point (0, 1, 0) of system (2) is a nondegenerate center
if, and only if, a7 6= 0, a4 = a9, and a4a5a8a9+a5a6a7a8+a

2
5a7a9+a5a8a

2
9−a7a28a9 = 0.



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 3

Also in Section 3, we will see that the behavior of linear and quadratic homo-
geneous vector fields is the same on all spheres. But this special property can not
be extended for all quadratic vector fields X2. Because of the difficulty of doing a
general study, we restrict our analysis to the unit sphere S2

1 = {(x, y, z) ∈ R3 :
x2 + y2 + z2 = 1}. In this case, we will show that, generically, any X ∈ X2 writes in
its canonical form as

ẋ = −a1y − a2z − a4xy − a5xz − a10y2 − (a6 + a7)yz − a8z2,
ẏ = a1x− a3z + a4x

2 + a10xy + a6xz − a11yz − a9z2,
ż = a2x+ a3y + a5x

2 + a7xy + a8xz + a11y
2 + a9yz.

(3)

In the following result we provide some center families for the above system when
a1 + a10 = 0, a2 + a7 = 1, a3 + a11 = 0, and a9 = 0. In Section 4 we will justify
why we have restricted our analysis to this special family. Although the study of
local cyclicity presented in the next Theorem 4 has been done for this family, when
the centers of Theorem 3 are perturbed without these conditions, no more limit
cycles of small amplitude are obtained. This is explained in Section 4.2. Finally, the
monodromic property around the equilibrium point is guaranteed by w2 = a6+1 > 0
and it is not restrictive to study only w > 0.

Theorem 3. The differential system

ẋ = −a1y − (1− a7)z − a4xy − a5xz + a1y
2 + (1− a7 − w2)yz − a8z2,

ẏ = a1x+ a11z + a4x
2 − a1xy + (w2 − 1)xz − a11yz,

ż = (1− a7)x− a11y + a5x
2 + a7xy + a8xz + a11y

2,

(4)

has a center at the equilibrium point (0, 1, 0) if a4 = 0 and one of the following
conditions is satisfied:

(a) w = 1, a1a5 + a8a11 = 0;
(b) a1 = 0, a8 = 0;
(c) a5 = 0, a11 = 0;
(d) a1 = a8, a5 = −a11;

(e) w 6= 1, a1 =
w2 − 1

w2 + 1
a8, a5 =

w2 + 1

w2 − 1
a11, a7 =

1

w2 + 1
− 1

(w2 + 1)
a28−

w2 + 1

(w2 − 1)2
a211.

We think that the above result provides (generically) a complete center classifica-
tion for the considered family. In Section 4 we have checked it for w ∈ {1/2, 2, 3},
but the expressions to be manipulated are too big to get the proof for every w.

Finally, we have also analyzed the local cyclicity of the families in Theorem 3,
studying bifurcations of small amplitude limit cycles from a weak focus on S2

1. Next
result provides the highest number of limit cycles surrounding a monodromic equi-
librium point in the quadratic family X2 that we have found.

Theorem 4. Consider the system

ẋ = 2αy +
9

20
z − xz − 2αy2 − 89

20
yz − αz2,

ẏ = −2αx+ 2z + 2αxy + 3xz − 2yz,

ż = − 9

20
x− 2y + x2 +

29

20
xy + αxz + 2y2.

(5)



4 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

The equilibrium point p = (0, 1, 0) of (5) is of center type if α = 0. Otherwise, it

is a weak focus of order 4 when α = ±
√

857/488, or of order 3 when α /∈ A :=

{0,±
√

857/488}. Moreover, if α = ±
√

857/488 (α /∈ A) there exist 4 (resp. 3)
small amplitude limit cycles, on S2

1, bifurcating from p considering a perturbation of
(5) inside family (4).

Remark 5. We notice that, near the centers of (5) for α = 0 only weak foci of
order three exist. Consequently, the cyclicity of the center should be also three.

Theorem 1 shows that there exists some similarity between the linear case on R2

and the linear case on S2, since for both we do not have limit cycles, but instead of
the many different types of equilibrium points for linear vector fields on R2, in X1

we only have centers. For planar quadratic polynomials vector fields, N.N. Bautin
proves (see [3]) that there exists at most 3 limit cycles bifurcating from a monodromic
equilibrium point. But our last main theorem exhibits a weak focus of order four
that unfolds 4 limit cycles in X2. That is, we have at least one limit cycle more than
in the classical planar case. Another final difference is obtained, considering also
the case with complex coefficients. In this case, we prove that the local (complex)
cyclicity of centers (e) in Theorem 3 is at least 6 perturbing inside family (4). The
fact that the cyclicity is higher when complex coefficients are considered, can explain
the difficulties that we have found to look for vector fields with high cyclicity values
in the real coefficients study, as well as in the center classification.

This paper is structured as follows. Section 2 is devoted to recall the necessary
classical results for proving our main theorems. In particular, the projection of a
3-dimensional vector field with invariant spheres and the definitions of Lyapunov
constants, first integrals, inverse integrating factors, and integrability. In Section 3
we study the linear and quadratic homogeneous case and we also prove Theorems 1
and 2. In Section 4 we study the quadratic case and we prove Theorems 3 and 4.
We finish studying the family considering that the coefficients are complex numbers.

2. Preliminary results

In this section, we recall some classical concepts and bifurcation techniques that
are necessary for the proofs of the results stated in the paper. Firstly, we introduce
the general vector fields X : R3 → R3 having H(x, y, z) = x2 + y2 + z2 as a first
integral together with some properties that they satisfy. Secondly, as the main
results will be proved projecting each 3-dimensional vector field to a planar one, we
recall some usual notions and definitions for planar vector fields. We mainly study
the center-focus and local cyclicity problems for polynomial systems. So we need
to introduce briefly the Darboux integrability concept and the computation of the
Lyapunov constants. This last notion is the usual planar mechanism to distinguish
when a monodromic equilibrium point is of center or of focus type.

2.1. Setting the problem. Consider a vector field X : R3 → R3 and its associated
differential system ẋ = X(x). Assume that 〈x, X(x)〉 = 0, for all x ∈ R3. That means
that X admits H(x, y, z) = x2 + y2 + z2 as a first integral. In other words all the
spheres of center at (0, 0, 0) and radius ρ, S2

ρ = {(x, y, z) : x2 + y2 + z2 = ρ2}, are
invariant by the flow of X. We recall that in the previous section we have denoted
by X this class of vector fields and by Xn when the components are polynomials of
degree n.



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 5

Next results show that an orthogonal change of coordinates keeps all the spheres
invariant and that it is not restrictive to assume that the equilibrium point, that
always exists, can be located at (0, ρ, 0).

Lemma 6. Let M be an orthogonal matrix. If X ∈ X then M ·X(M t) ∈ X.

Proof. The differential equation ẋ = X(x), with the orthogonal change of coordi-
nates y = M · x, moves to ẏ = M · ẋ = M ·X(x) = M ·X(M t · y) = Y (y) because
M−1 = M t. The proof follows just checking that 〈y, Y (y)〉 = 〈M · x,M ·X(M t · y)〉
= (M · x)t ·M ·X(M t · y) = xt ·X(x) = 〈x, X(x)〉 = 0 for all y ∈ R3. �

Lemma 7. The equilibrium point of X ∈ X can be always located at (0, ρ, 0).

Proof. Let p = (x0, y0, z0) ∈ S2
ρ an equilibrium point of X. Then X(p) = 0 and

‖p‖ = ρ > 0. Consider the unit vector v = (x0, y0, z0) /ρ and the plane

P =

{
(x, y, z) :

xx0
ρ

+
yy0
ρ

+
zz0
ρ

= 0

}
,

passing through the origin, which is perpendicular to the vector v. Note that the
intersection of P with the sphere x2 + y2 + z2 = 1 is a circumference and we can
take u = (0,−z0, y0)/

√
y20 + z20 as a unit vector on this circumference. So, we have

a new orthogonal basis B = {u, v,w}, where w is obtained by the cross product of u
and v. The matrix M that changes the canonical bases to B is orthogonal and this
change of coordinates sends the equilibrium point p = (x0, y0, z0) to (0, ρ, 0). We
notice that Lemma 6 ensures that the new vector field also is in class X. �

The stereographic projection with respect to the antipodal point of the equilibrium
point located at (0, ρ, 0), allow us to consider planar vector fields instead of 3-
dimensional vector fields restricted to spheres.

Let π : S2
ρ \ {(0,−ρ, 0)} → R2 be the stereographic projection on the plane

{(x, y, z) ∈ R3 : y = ρ} given by π(x, y, z) = 2ρ(x, z)/(y + ρ). Then, the projection
Y : R2 → R2 of the vector field X writes as

Y (x) = dππ−1(x) ◦X ◦ π−1(x), (6)

where X = X|S2ρ . Note that π preserves closed curves and contact between curves

contained on its domain of definition. We say that p ∈ S2
ρ is an equilibrium point

of center type of X|S2ρ if π(p) = q and q is an equilibrium point of center type of

Y (x) defined in (6). Moreover, as π(0, ρ, 0) = (0, 0) we can assume that (0, 0) is an
equilibrium point of the planar projected system (6).

There are other transformations to project a 3-dimensional vector field to a planar
one. For example the one used in [17, 18]. Although these works study vector fields
defined on the sphere, the objectives are completely different. In particular, they
study when the maximal circles are invariant for the considered vector fields. Some
of the above properties are used but without proving them.

2.2. Darboux Integrability. Consider the planar differential system

ẋ = P (x, y),

ẏ = Q(x, y),
(7)

where P and Q are polynomials in the variables x and y with coefficients on F,
where F is the field of reals or complex numbers. As usual, we also denote by F[x, y]



6 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

the ring of polynomials in the variables x and y and coefficients in F. Let m be
the maximum between the degree of P and Q. We say that (7) is integrable or
completely integrable on an open subset U ⊂ F2 if there exists an analytic function
H : U → F, non-constant, such that

P
∂H

∂x
+Q

∂H

∂y
= 0,

on U . In this case, H is called a first integral. In addition, a vector field X : R3 → R3

is completely integrable if it has two independent first integrals. Note that, if X ∈ X
then, by definition, it has at least one first integral. We say that an analytic (and
not identically zero) function V : U → F is an inverse integrating factor of (7) on
U if

div(P/V,Q/V ) =
∂(P/V )

∂x
+
∂(Q/V )

∂y
≡ 0,

and, consequently, a first integral can be obtained by direct integration.
Let f ∈ C[x, y] being not identically zero. We say that f(x, y) = 0 is an invariant

algebraic curve of (7) if there exists a cofactor K ∈ C[x, y] such that

Xf = P
∂f

∂x
+Q

∂f

∂y
= Kf. (8)

If the degree of the polynomial differential system is m then, the cofactor K has
degree at most m − 1. By (8) we could see that the gradient of f is orthogonal
to the vector field (7) and, because of that, the flow of (7) is tangent to the curve
f(x, y) = 0. Thus this curve is formed by trajectories of the vector field.

An object that also satisfies (8) is the exponential factor, that we define on the
following. Let g, h ∈ C[x, y] such that g, h are relative prime in C[x, y] or h ≡ 1.
The function exp(g/h) is called an exponential factor of (7) if there exists a cofactor
K ∈ C[x, y] of degree at most m− 1 such that

X
(

exp
(g
h

))
= K exp

(g
h

)
.

Note that although the function exp(g/h) satisfies (8), as it is always nonzero, an
exponential factor do not define invariant curves.

The Darboux Integrability Theory for complex polynomial systems gives us some
conditions in which the existence of invariant algebraic curves and exponential fac-
tors ensure that the vector field is integrable. In this sense, suppose that the vec-
tor field (7), of degree m, admits p irreducible invariant algebraic curves fi = 0,
with cofactors Ki, i = 1, . . . , p and q exponential factors exp(gj/hj) with cofac-
tors Kp+j, j = 1, . . . , q. Then, there exist λi, µj ∈ C not all zero such that∑p

i=1 λiKi +
∑q

j=1 µjKp+j = 0 if and only if the (multivalued) function

fλ11 · · · fλpp
(

exp

(
g1
h1

))µ1
· · ·
(

exp

(
gq
hq

))µq
(9)

is a first integral of system (7). Beside of that, (9) is an inverse integrating factor for
system (7) if, and only if, there exist λi, µj ∈ C not all zero such that

∑p
i=1 λiKi +∑q

j=1 µjKp+j = div(P,Q). We emphasize that fi and exp(gj/hj) could be complex

but when system (7) is real, the Darbouxian function (9) is also real. See [10] for
more details.



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 7

2.3. Lyapunov constants and local cyclicity. We will recall the stability algo-
rithm for equilibrium points of nondegenerate center-focus type having its linear
part in Jordan’s normal form:

ẋ = αx− βy +
n∑
k=2

Pk(x, y),

ẏ = βx+ αy +
n∑
k=2

Qk(x, y),

(10)

where Pk and Qk are homogeneous polynomials of degree k in the variables x and
y. The nondegeneracity condition is β 6= 0. The above system writes, in usual polar
coordinates, (x, y) = (r cos θ, r sin θ), as

ṙ = R(r, θ),

θ̇ = β + Θ(r, θ),

where R and Θ are polynomials in r and trigonometric polynomials in cos θ and
sin θ. Removing the time dependence we consider the 1-dimensional differential
equation

dr

dθ
=

R(r, θ)

β + Θ(r, θ)
= R(r, θ), (11)

which is well defined in a small neighborhood of the origin.
Let r(θ, r0) be the solution of (11) satisfying r(0, r0) = r0. If α 6= 0 by the

Grobman-Hartman Theorem we conclude that (0, 0) is a hyperbolic focus and the
sign of α provides the stability of the origin. When α = 0, the stability depends on
the higher-order terms. For r0 sufficiently small we expand the solution in Taylor’s
series and write

r(θ, r0) = r0 +
∞∑
k=2

rk(θ)r
k
0 ,

with rk(0) = 0, for all k ≥ 2. The Poincaré first return map is defined evaluating
the above solution at 2π:

Π(r0) = r(2π, r0).

The corresponding solution of (10) turns around the origin. The stability of the
origin depends on the sign of the displacement function

∆(r) = Π(r0)− r0,
for r0 is small enough. If there exists k such that ∆′(0) = · · · = ∆(k−1)(0) = 0 and
∆(k)(0) 6= 0 then k is always an odd integer number and we write k = 2K+1, forK ≥
0. The K-th Lyapunov constant is defined as LK = r2K+1(2π) = ∆(2K+1)(r0)/(2K+
1)! when L1 = · · · = LK−1 = 0 and LK 6= 0. Then, we say that the origin of system
(10) is a weak focus of order K if there exists K ≥ 1 such that LK 6= 0 otherwise
we say that the origin is a center.

The above described method to study the center-focus problem comes from Lya-
punov. There is an alternative method due to Poincaré, also when α = 0. It consists
in looking for a function H(x, y) = x2 + y2 +O3(x, y) which satisfies

Ḣ =
dH

dt
= P

∂H

∂x
+Q

∂H

∂y
=
∞∑
k=2

h2k r
2k,



8 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

being r2 = x2 + y2. The first nonvanishing coefficient of Ḣ, which always has an
even subscript, h2K , is called a focal value and determines the stability of the origin
of (10). In this case H acts as a Lyapunov function. It is well known that, both
coefficients, h2K+2 and L2K+1 differs only by a multiplicative nonzero constant and
by the Poincaré–Lyapunov Theorem the origin is a center of (10) if, and only if,
h2K = 0 for all K, in others words, if and only if (10) admits an analytical first
integral. We notice that the Taylor series of H(x, y) = x2 + y2 +O3(x, y) converges
to an analytic first integral because the differential system also is. For more details
we refer the reader to [2, 10, 21].

From the computational point of view, it is better to consider the last method in
complex variables z = x+ i y. In this case, if α = 0 the system (10) writes as

ż = R(z, z̄) = i z +
n∑
k=2

Rk(z, z̄),

being Rk(z, z̄) homogeneous polynomials of degree k in z and z̄. We only write the
first equation because, as (10) has real coefficients, the second one is its conjugate.
In this case we have H(z, z̄) = zz̄ +O3(z, z̄) and, consequently,

dH

dt
= Ḣ = ż

∂H

∂z
+ ˙̄z

∂H

∂z̄
=
∞∑
k=2

gk(zz̄)2k.

Thus the focal values are the coefficients gk.

It is well known that the Lyapunov constants and the focal values differ in multi-
plicative constants and that they are polynomials in the coefficients of Pk and Qk.
Moreover, each LK is always defined modulus the previous vanish. So, independently
of the used mechanism, we will denote them by Lk. The last approach provides some
good algebraic properties about degree and weighted-homogeneity with respect to
the perturbation parameters. See them in [9]. For more details on the center-focus
problem and related problems, we refer the reader to [2, 20].

The classical Hopf bifurcation occurs near α = 0 and when the first Lyapunov
constant (when α = 0) is non-vanishing. It is not restrictive if we assume L1 > 0.
In this case, when α = 0, the origin is unstable and a small amplitude stable limit
cycle bifurcates from the origin when the trace parameter α becomes negative but
small enough. The degenerate Hopf bifurcation occurs when the limit cycles of small
amplitude bifurcate from a weak focus of higher order. It is known that at most
K limit cycles bifurcate from an order K weak focus under analytic perturbations
(see [21]), but the unfolding is not always complete when the perturbation is re-
stricted to a polynomial family of fixed degree. As we are interested in polynomial
perturbations we need a simple condition for proving the existence of K limit cycles
bifurcating from the origin. Instead of looking for weak focus of higher order and its
unfolding, we will study degenerate Hopf bifurcations from centers. The key point
is an interesting application of the Implicit Function Theorem due to Chicone and
Jacobs in [7]. Also Han ([14]) uses it. We will present the approach of Christopher
([8]) that uses the Taylor developments of first-order of the Lyapunov constants with
respect to the pertubation parameters.

Theorem 8 ([8]). Suppose that p is a point on the center variety and that the first
k-Lyapunov constants, L1, . . . , Lk have independent linear parts (with respect to the



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 9

Taylor expansion of Li about p), then p lies on a component of the center variety
of codimension at least k and there are bifurcations which produce k limit cycles
locally from the center corresponding to the parameter value p. If, furthermore, we
know that p lies on a component of the center variety of codimension k, then p is a
smooth point of the variety, and the cyclicity of the center for the parameter value
p is exactly k. In the latter case, k is also the cyclicity of a generic point on this
component of the center variety.

We notice that in the above result we are moving also the trace parameter α. In
[8] we can also find another result for bifurcating limit cycles of small amplitude
using higher order developments of the Lyapunov constants. These higher order
studies are better explained and developed in [11, 13].

3. Centers for linear and quadratic homogeneous vector fields

In this section we classify the type of equilibrium points that a linear vector field
X ∈ X can have proving Theorem 1. As we have already commented before, when
X ∈ X is linear then X is also homogeneous, so X ∈ XH

1 . We have used this
motivation for studying the centers in XH

2 . They are classified in Theorem 2. Before
proving these results, we give some technical lemmas about homogeneous vector
fields on X.

Lemma 9. The homogeneity property is invariant by an orthogonal change of co-
ordinates.

Proof. The proof follows directly from the proof of Lemma 6, just checking that the
change of variables does not break the homogeneity. �

Lemma 10. Let X ∈ XH . The phase portrait in each sphere is topologically equiv-
alent to the one in the sphere of radius 1. Moreover, X has a straight line passing
through the origin filled of equilibrium points.

Proof. The proof follows just doing the change of coordinates y = x/ρ and a time
rescaling, if necessary. �

In the above lemma the existence of an equilibrium point in one sphere is trans-
formed by continuity using a dilation when the vector field is homogeneous. We
notice that this is not the case for a general vector field in X nor for the quadratic
family X2.

Lemma 11. Let X ∈ XH . Then the projected system Y defined in (6) is homoge-
neous if, and only if, X is linear.

Proof. First we will study the projection of X ∈ XH
1 defined in (1). By Lemma 7

we can suppose that any equilibrium point can be located at (0, ρ, 0). This implies
that a1 = a3 = 0. So, system (1) writes as

(X1, X2, X3) = (−a2z, 0, a2x)

and the corresponding projected system (6) is

(Y1, Y2) = (−a2v, a2u). (12)

Which is also homogeneous and linear.
The proof finishes just checking that for degree two, for example, the projection

breaks the homogeneity property. Although we will see in the proof of Theorem 2



10 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

the general case, here we consider a particular system (2) choosing, for example
a4 = a7 = 1 and a5 = a6 = a8 = a9 = 0. Then, the corresponding projected
system (6) writes as the nonhomogeneous cubic system

u̇ = −4u− 4 v − u3 + u2v + uv2 + v3,

v̇ = 4u− u3 − 2u2v − uv2.
�

3.1. Linear case. We recall that if X ∈ X1 then it writes in the canonical form
(1). On the following we prove Theorem 1.

Proof of Theorem 1. By Lemma 10 all the spheres are equivalent, so we can restrict
to the sphere of radius ρ = 1. The proof follows directly from Lemma 11, because the
projected system (12) has a center at the origin when it is an isolated equilibrium
point. That is, when a2 6= 0.

Note that, H2(x, y, z) = a3x−a2y+a1z−k is a first integral for the linear system
(1). As, by definition, H1(x, y, z) = x2 + y2 + z2 is also a first integral, system (1) is
completely integrable. �

3.2. Quadratic homogeneous case. The next technical result ensures that we
can assume that the equilibrium point is located at (0, 1, 0) and so we only work
with X defined on S2

1. Then, we prove our second main result.

Lemma 12. The canonical form of X ∈ XH
2 is system (2).

Proof. The proof follows straightforward using Lemmas 7, 9, and 10. �

Proof of Theorem 2. Using Lemma 12 we can consider system (2) restricted to the
sphere S2

1 and its projection Y defined in (6). We restrict our attention to the
equilibrium point of (2) which is located at the origin after projection. It will
be of nondegenerate center-focus type if the Jacobian matrix J associated to the
projected vector field Y has zero trace and positive determinant. Straightforward
computations ensure that it occurs if, and only if, a4−a9 = 0 and a6a7+a27−a29 > 0.
So, under these conditions we have a weak focus at the origin and, writing w2 =
a6a7 + a27 − a29, the projected system Y can be written as

u̇ = −4a4u−4ξv−4a5uv−4a8v
2−a4u3−(ξ − 2a7)u

2v+(a4 + 2a9)uv
2 + ξv3,

v̇ = 4a7u+4a9v+4a5u
2+4a8uv−a7u3−(2a4 + a9)u

2v−(2ξ − a7)uv2+a9v
3,

(13)

where ξ = (w2 + a4a9)/a7. It is easy to check that the trace and determinant of J
are −4(a4 − a9) and 16w2, respectively. So, when a4 6= a9 the origin is a hyperbolic
focus for system (13) and a4 − a9 = 0 provides the first condition in the statement
to have a weak focus. In this last case, as we have explained in Section 2.3, the
stability depends on the computation of the Lyapunov constants corresponding to
the origin of system (13). Here we only need to compute the first one, that is

L1 =
16(a27 + a29 + w2)C

3(a29 + w2)2
,

where C = −a25a7a9 + a5a
2
7a8 − a5a8a29 − a5a8w2 + a7a

2
8a9. As w 6= 0, the condition

C = 0 is also necessary to have a center at the origin. The second condition in the
statement follows substituting the value of w2 in the above expression of C.



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 11

The proof finishes just showing that under the two conditions in the statement,
system (2) is time-reversible with respect to an straight line passing trough the
origin. We will show this property proving that there exists ϕ such that the trans-
formation

(û, v̂) = (cosϕu− sinϕv, sinϕu+ cosϕv) (14)

changes the vector field (13) to a new one Yϕ(û, v̂) that is time-reversible with respect
to the new v̂-axis, i.e., its phase portrait is invariant under reflection with respect
to the v̂-axis in the direction of time. In other words, it is invariant with respect to
the change (û, v̂, t) 7→ (−û, v̂,−t), being t the time variable. Note that if a system
is time-reversible with respect to a line, then the equilibrium points on it are not
attractors or repellers. See [16] for more details.

Using the rational parameterization sinϕ = 2τ/(1+ τ 2) and cosϕ = (1− τ 2)/(1+
τ 2), the vector field Yϕ(û, v̂) is time-reversible with respect to the v̂-axis if, and only
if,

1

a7(τ 2 + 1)2
(
(4 + û2 − 3v̂2)f1(τ)û+ 4a7(τ

2 + 1)f2(τ)ûv̂
)
≡ 0,

1

a7(τ 2 + 1)2
(
(4− 3û2 + v̂2)f1(τ)v̂ + 4a7(τ

2 + 1)f2(τ)û2
)
≡ 0,

where

f1(τ) = −a7a9τ 4 + 2(a27 − a29 − w2)τ 3 + 6a7a9τ
2 − 2(a27 − a29 − w2)τ − a7a9,

f2(τ) = a5τ
2 + 2a8τ − a5.

Therefore, we only need to check that there exists a common root of the polynomials
f1(τ) and f2(τ). We notice that f2(τ) always has simple real solutions or it vanishes
identically (a5 = a8 = 0). The resultant between f1(τ) and f2(τ) with respect to
τ is res(f1, f2, τ) = 16C2. Hence, there exists a solution for {f1(τ) = f2(τ) = 0}
if, and only if, our second condition C = 0 is satisfied. When a9 6= 0, as a7 6= 0,
there exists a real solution because f1(0)f1(1) = −(a7a9)

2. While when a9 = 0 we
have C = a5a8(a

2
7 − w2) and f1(τ) = 2(a27 − w2)τ(τ 2 − 1) and some cases must be

distinguished: a5 = 0, a8 = 0, or a7 = ±w. When a5 = 0 we take τ = 0 and when
a8 = 0 we take τ = 1. If a7 = ±w, then f1(τ) ≡ 0 and the proof is finished because
f2(τ) always has real solutions. �

As we have proved in Lemmas 7 and 10, when X ∈ XH a straight line of equilib-
rium points always exists and it passes through (0, 1, 0) and (0,−1, 0). The projected
vector field (6) has an equilibrium at the origin and another at infinity. Moreover,
when the degree n is odd, the system is invariant by the change of coordinates
(x, y, z) 7→ (−x,−y,−z) and the point (0, 1, 0) moves to (0,−1, 0). On the other
hand, when n is even, the invariance needs an inversion of the time. So, if the origin
is a center of (13), the infinity is also a center that rotates in the opposite (resp.
same) direction when n is even (resp. odd). This property is exhibited in Theo-
rem 1 for vector fields in XH

1 because the unique phase portrait is a global center.
For vector fields in XH

2 the centers in Theorem 2 located at the origin and at infinity
turns in opposite directions and we can have different global phase portraits. We
can see two of them in Figure 1, where we have drawn the phase portraits of the
one-parameter cubic family (15) studied in next Proposition 13. We notice that the
systems are time-reversible with respect to the x-axis, so the origin and infinity are
simultaneously centers. The other symmetric points can be centers or antisaddles,



12 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

in this last case, with opposite stability. In [5, 6], we can find the classification of the

Figure 1. Phrase portrait of (15) for a = 0 (left) and a > 0 (right).

global phase portraits of reversible cubic centers where the infinity and the origin
rotate in the same direction.

Proposition 13. The system

u̇ = −8v − 4auv − u2v + 2v3,

v̇ = 2u+ 4au2 − 7

2
uv2 − 1

2
u3,

(15)

has centers at the origin and at the infinity simultaneously for all a. These centers
rotate in opposite directions and the only possible phase portraits are those shown in
Figure 1.

Proof. We notice that system (15) is (13) choosing a5 = a, a7 = 1/2, a8 = a9 = 0,
and w = 1. So, by Theorem 2 the origin is a center. We have already explained that,
by the symmetry of the corresponding vector field in R3, the infinity is also a center
which rotates in an opposite direction with respect to the origin. With the change
of coordinates (u, v) 7→ (−u,−v), if necessary, we can assume that a ≥ 0. Further-
more, system (15) is time-reversible with respect to the u-axis and if a = 0 it is also
time-reversible with respect to the v-axis. In particular, as the line of infinity has
no equilibrium points, we do not need to use the Poincaré compactification to study
the dynamics near the infinity. Only the finite real equilibrium points are necessary
to be analyzed, and the nonexistence of limit cycles property. Straightforward com-
putations show that they are (0, 0), A± = (4a±2

√
4a2 + 1, 0) and B± = (0,±2). Let

J be the Jacobian matrix of the vector field associated to (15). Then A± are saddles
because of the time-reversibility and the fact that the determinant of J at A± is
−48(2a

√
4a2 + 1 ± 4a2 ± 1)2 < 0. The matrix J at the equilibrium points B± has

positive determinant. So, the local stability is determined by the sign of the trace of
J, which is ∓8a. So, as a > 0, B+ is an attractor and B− is a repeller. When a = 0
the symmetry with respect to the v-axis proves the existence of a center at each of
those points, obtaining the phase portrait depicted in Figure 1 (left). In this case
we have also the first integral

H(u, v) =
u2 + 4v2

(u2 + v2 + 4)2
.



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 13

When a 6= 0, the Lie derivative of H with respect to (15) is

u̇
∂H

∂u
+ v̇

∂H

∂v
=

24au4v

(u2 + v2 + 4)2

and it is always positive (resp. negative) when v > 0 (resp. v < 0). So, we can use
H as a Lyapunov function on {(u, v) ∈ R2 : v > 0} (resp. on {(u, v) ∈ R2 : v < 0}).
Hence, if a 6= 0 we have no periodic orbits completely contained in the half planes
v > 0 or v < 0. Considering this last property, together with the time-symmetry
and the local phase portraits of all equilibrium points, we have that the global phase
portrait is the one depicted in Figure 1 (right) when a > 0. �

4. Centers and cyclicity for quadratic vector fields

In this section we fix our attention to the quadratic vector fields (3), proving our
main Theorems 3 and 4. Firstly we will prove some center characterization and
secondly some results on limit cycle bifurcation near centers. Before the proofs,
we will show with the next example that the qualitative behavior of a vector field
in X2 on a sphere of radius ρ1 centered at the origin can be totally different from
the behavior on another sphere of radius ρ2 6= ρ1. In particular, the number of
equilibrium points can change.

Example 14. All the spheres x2 + y2 + z2 = ρ2 are invariant for the quadratic
system

(ẋ, ẏ, ż) = (−xz − yz − z2 − z,−z2, x2 + xy + xz + yz + x).

There are two straight lines full of equilibrium points: {x = z = 0} and {x+y+ 1 =
z = 0}. So, for every fixed invariant sphere of radius ρ centered at the origin, we
have four equilibrium points when ρ > 1/

√
2, three when ρ = 1/

√
2, and only two

when ρ < 1/
√

2.

In the following, we will study the existence of periodic orbits and limit cycles on
the invariant spheres. Due to the above example, we will concentrate our efforts to
study them fixing one sphere, describing the centers and the existence of limit cycles
of small amplitude. So, we fix our attention to the unit sphere S2

1 = {(x, y, z) ∈
R3 : x2 + y2 + z2 = 1} and, as we have explained previously, we will assume that
(0, 1, 0) is an equilibrium point of a quadratic vector field X ∈ X2. This fact forces
the conditions a10 + a1 = 0 and a3 + a11 = 0 in system (3), that we can write as

ẋ = −a1y − a2z − a4xy − a5xz + a1y
2 − (a6 + a7)yz − a8z2,

ẏ = a1x+ a11z + a4x
2 − a1xy + a6xz − a11yz − a9z2,

ż = a2x− a11y + a5x
2 + a7xy + a8xz + a11y

2 + a9yz.

(16)

We notice that (0,−1, 0) is an equilibrium point of the above differential system
if, and only if, a1 = a11 = 0. So, in contrary with the homogeneous case, not
always exists the line of equilibrium points passing through the two antipodal points
(0,±1, 0).

The main objective of the next subsections is to exhibit some differences between
the quadratic vector fields in the plane with respect to the quadratic vector fields
in the sphere S1, emphasizing in the number of limit cycles of small amplitude that
can bifurcate from a monodromic nondegenerate equilibrium point. We will see that



14 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

the weak focus order is higher in S2
1 than in R2 and also the number of limit cycles

of small amplitude.

4.1. Center characterizations. In order to have a monodromic equilibrium point
at (0, 1, 0), we will add some extra conditions in the projected system Y obtained
doing the transformation (6) to the quadratic differential equation (16). As before,
we denote by J the Jacobian matrix associated to Y at an equilibrium point. The
origin is a monodromic nondegenerate equilibrium point of Y if, and only if, the trace
and the determinant of J are zero and positive, respectively. That is, when a4 = a9
and a2a6 +a6a7 +2a2a7 +a22 +a27−a29 > 0. Due to the big number of free parameters
and to simplify a little the computational difficulties, we will restrict our analysis
adding two extra conditions: a9 = 0 and a2 + a7 = 1. In this case, the projected
vector field Y has a weak focus at the origin if, and only if, a4 = 0 and a6 + 1 > 0.
Moreover, when the trace is zero the matrix J is in the real Jordan normal form.
Taking all into account and writing w2 = a6 + 1, with w 6= 0, we obtain system (4).
After a reparametrization of the time, the corresponding projected system for (4) is

u̇ = −a4
w
u− v − a1

2
u2 − a5

w
uv − a1 + 2a8

2w2
v2 − a4w

4
u3 +

2a7 − w2

4
u2v

+
a4
4w

uv2 +
w2 + 2a7 − 2

4w2
v3 − a1w

2

8
u4 − a11w

4
u3v − a11

4w
uv3 +

a1
8w2

v4,

v̇ = u+
(2a5 − a11)w

2
u2 + a8uv −

a11
2w

v2 − (2a7 − 1)w2

4
u3 − a4w

2
u2v

− 2w2 + 2a7 − 3

4
uv2 +

w3a11
8

u4 − w2a1
4

u3v − a1
4
uv3 − a11

8w
v4.

(17)

The following result shows that when we fix the value of the determinant of the
Jacobian matrix, or equivalently w, the center-focus problem can be completely
solved. The proof of case (a) in Theorem 3 follows directly from it.

Proposition 15. Consider system (4) with w = 1. Then, (0, 1, 0) is a center on S2
1

if, and only if, a4 = 0 and a1a5 + a8a11 = 0.

Proof. Instead of working with system (4) we will work with the equivalent projected
planar differential system (17). It is easy to check that if a4 6= 0 then (0, 0) is a
hyperbolic focus. So a4 = 0 is the first center condition in the statement. The
second center condition detailed in the statement, a1a5 + a8a11 = 0, follows using
the algorithm described in Section 2.3 for computing the Lyapunov constants for
providing the stability of the origin. In fact, it appears as a common factor in the
first three:

L1 = −2

3
(a1a5 + a8a11),

L2 =
1

15
(a1a5 + a8a11)(3a

2
1 − 6a1a8 + a25 − 18a5a11 + 51a211 − 7a28 − 6a7),

L3 = − 1

630
(a1a5 + a8a11)(3a

4
1 − 6a31a8 + 195a21a

2
5 − 858a21a5a11 − 186a21a

2
11

+ 543a21a
2
8 + 478a1a

2
5a8 − 832a1a5a8a11 − 7182a1a8a

2
11 + 798a1a

3
8 − 18a45

+ 242a35a11 − 509a25a
2
11 + 72a25a

2
8 − 6242a5a

3
11 + 2146a5a8a

2
11 + 14267a411

− 8041a28a
2
11 + 282a48 − 276a21a7 + 600a1a7a8 − 84a25a7 + 1752a5a7a11



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 15

− 6276a7a
2
11 + 516a7a

2
8 − 108a21 − 324a1a8 + 24a25 − 324a5a11 + 540a211

+ 144a27 − 168a28 − 72a7).

The next step is to prove that under these two conditions the origin is a time-
reversible center with respect to a straight line. It follows using the same idea as
in the proof of Theorem 2. That is, from the change (14) again with the rational
parameterization sinϕ = 2τ/(1 + τ 2) and cosϕ = (1− τ 2)/(1 + τ 2). Straightforward
computations show that the transformed vector field (with a4 = 0) is invariant with
respect to the change (û, v̂, t) 7→ (−û, v̂,−t), being (û, v̂) the new variables, if and
only if

ûv̂((û2 + v̂2)f1(τ)− 4f2(τ)) ≡ 0,

(û2 + v̂2)(û2 − v̂2 − 4)f1(τ)− 8f2(τ)û2 ≡ 0,

where f1(τ) = −a11τ 2 + 2a1τ + a11 and f2(τ) = a5τ
2 + 2a8τ − a5. The proof finishes

because, when the second condition a1a5+a11a8 = 0 holds, f1 and f2 have a common
real root, which provides the symmetry line. �

Now we prove the remaining center families in our third main result.

Proof of Theorem 3. As in the above proof, we will start with the projected system
(17). Moreover, we assume a4 = 0, otherwise we have a hyperbolic focus at the
origin.

The case (a) follows directly from Theorem 15 and cases (b) and (c) because
the corresponding systems are time-reversible with respect to the u-axis and v-axis,
respectively. The system in the fourth case (d) is Darboux integrable, having the
next rational first integral

H(u, v) =
(a11wu− a8v − 2a7 + 1)(w2u2 + v2)− (w2 − 1)v2 − 4(a7 − 1)

(w2u2 + v2 + 4)2
.

The remaining case (e) is also Darboux integrable. But the proof of the existence of
an inverse integrating factor is more intricate. We will prove the existence of four
(complex) invariant straight lines for the projected system (17) and then we will
provide an inverse integrating factor of the form (9).

We denote by (Y1, Y2) the vector field (17) and by F = au + bv + 1 a generic
straight line. We recall that, if F is invariant, its respective cofactor K will be a
polynomial of degree 3. Equating the coefficients in u and v of the identity

Y1
∂F

∂u
+ Y2

∂F

∂v
= F K,

we can obtain four possible pairs of values (a, b) 6= 0 that depend on the (complex)
roots of the polynomial

pA(Z) = Z4 + A3Z
3 + A2Z

2 + A1Z
4 + A0, (18)

where A0 = (w2−1)2
16(w2+1)2

, A1 = w2−1
2(w2+1)2

a8, A2 = 1
(w2−1)2a

2
11 + 1

(w2+1)2
a28 − w2−1

2(w2+1)
,

and A3 = −2
w2+1

a8. Note that, if z1 and z2 are roots of (18) then z3 = 1−w2

4z1(w2+1)

and z4 = 1−w2

4z2(w2+1)
also are. So, writing the parameters a8 and a11 of (17) in

terms of the two roots z1 and z2, that is, a8 = (z1+z2)(4z1z2(w2+1)−w2+1)
8z1z2

and a11 =
i(z1−z2)(w2−1)(4z1z2(w2+1)+w2−1)

8z1z2(w2+1)
, we have a complete factorization of the polynomial



16 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

pA(Z) =
∏4

i=1(Z − zi). We notice that we are assuming that all the coefficients are
real, so we have that z2 = z̄1 and z4 = z̄3. Hence, we can write the four invariant
straight lines and its respective cofactors as

F1 = 1 + i z1wu+ z1v,

F2 = 1− i z2wu+ z2v,

F3 = 1− (w2 − 1)(v + iwu)

4(w2 + 1)z1
,

F4 = 1− (w2 − 1)(v − iwu)

4(w2 + 1)z2
,

and

K1 =
1

32z2w(w2 + 1)
((−w3 + w)u+ i((w2 − 1)v − 4z1(w

2 + 1))

w2(4z22(w2 + 1)− w2 + 1)u2 + (4z22(w2 + 1)− w2 + 1)v2 + 8z2w
2v + 8 i z2wu),

K2 =
1

32z1w(w2 + 1)
((w3 − w)u+ i((w2 − 1)v − 4z2(w

2 + 1))

(w2(−4z21(w2+ 1) + w2− 1)u2− (4z21(w2+ 1)− w2+ 1)v2 − 8z1w
2v + 8 i z1wu),

K3 =
(−z1wu+ i(z1v + 1))(w2 − 1)

z1((−w3 + w)u+ i((w2 − 1)v − 4z1(w2 + 1)))
K1,

K4 =
(z2wu+ i(z2v + 1))(w2 − 1)

z2((w3 − w)u+ i((w2 − 1)v − 4z2(w2 + 1)))
K2.

The proof finishes showing that
∏4

i=1 F
λi
i is an inverse integrating factor of the

system. This last property holds because taking

λ1 =
4(w2 + 1)z21 + 2(w2 − 1)

4(w2 + 1)z21 + w2 − 1
, λ2 =

4(w2 + 1)z22 + 2(w2 − 1)

4(w2 + 1)z22 + w2 − 1
,

λ3 =
8(w2 + 1)z21 + w2 − 1

4(w2 + 1)z21 + w2 − 1
, λ4 =

8(w2 + 1)z22 + w2 − 1

4(w2 + 1)z22 + w2 − 1
,

the condition
4∑
i=1

λiKi = div(Y1, Y2)

is satisfied. Moreover, not all λi, for i = 1, . . . , 4, are zero. �

Remark 16. We notice that the above proof is also valid when all the parameters
are complex. But in this case we will have no relation between the two main roots
z1 and z2.

If we fix the value of the special parameter w we conclude that there are no other
centers cases different from the ones detailed in Theorem 3. But the expressions
that appear are very big to provide a complete proof for an arbitrary w. This is why
in the next result we have chosen some explicit values.

Proposition 17. Consider the differential system (3) with w ∈ {1/2, 2, 3}. Then,
(0, 1, 0) is a center if, and only if, a4 = 0 and one of the conditions (b), (c), (d), or
(e) in Theorem 3 is satisfied.



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 17

Proof. From the proof of Theorem 3 we know that all the families detailed in the
statement are centers. Hence, we only need to check that there are no others. We
present only the proof for the case w = 2. The key point is the computation of
enough center conditions as it was done for Proposition 15, but in this case we need
to compute more Lyapunov constants using the method explained on Section 2.3 for
system (17). We have needed six to finish the proof. In the proof we will denote by
SL = {L1 = L2 = · · · = L6 = 0} the system of equations needed to be solved. The
other cases w = 3 and w = 1/2 are completely analogous.

So, we fix w = 2 and the first Lyapunov constant writes as

L1 = −a1
37a5 − 15a11

48
+ a8

15a5 − 11a11
24

.

Clearly, it vanishes when a1 = a8 = 0 which is Family (b), so we have centers in this
case. So, if necessary, we can assume that a1 and a8 do not vanish simultaneously.
The same conclusion holds for a5 = a11 = 0 (Family (c)). Consequently, we have

two possibilities a1 = 15a5 − 11a11 = 0 or a8 = a1(37a5−15a11)
2(15a5−11a11) .

The first case follows easily computing the next Lyapunov constants using that
a5 and a8 do not vanish, otherwise we obtain the previous studied families. From
L2 = 0 we obtain

a7 = −507a25 + 1452a28 + 15367

13915
,

and then Li = a5a8Li(a5, a8), for i = 3, 4, 5, being Li polynomials of degree 4, 6,
and 8 respectively. We have no other families in this case because the system of
equations {L3 = L4 = L5 = 0} reduces to {225a25− 363 = 225a28 + 9408 = 0}, which
has no real solutions.

For the second case, we can assume that a8 = a1(37a5−15a11)
2(15a5−11a11) and 15a5 − 11a11 6= 0.

Straightforward computations show that L1 = 0 and

L2 =
a1(a5 + a11)

320(15a5 − 11a11)3
(
196(2a5 − a11)(8a5 − 5a11)(42a5 − 23a11)a

2
1

+ (672a35 − 1124a25a11 + 624a5a
2
11 − 115a311 − 543a5 + 576a11)(15a5 − 11a11)

2

+ 15(106a5 − 67a11)(15a5 − 11a11)
2a7
)
.

(19)
When 106a5 − 67a11 = 0, the numerator of the above expression, up to a nonzero
rational factor, writes as a1a

4
5(1149184a21 +33856a25 +59367025) and we have no new

center families when L2 vanishes. When a5 + a11 = 0 we have a1 = a8 which is
Family (d) in Theorem 3. Therefore, on the following we write a7, from (19), as a
rational function of (a1, a5, a11). The next Lyapunov constants write as

Li =
a1(a5 + a11)(3a5 − 5a11)Li(a1, a5, a11)
(106a5 − 67a11)i−1(15a5 − 11a11)2i−1

(20)

for i = 3, 4, 5, 6, being Li polynomials with rational coefficients of degrees 5(i − 1)
and 50, 120, 235, 406 monomials, respectively. We do not write them here because
of their size. The common factor of the above expressions get, respectively, the
Families (b), (d), and (e) in Theorem 3.

The proof finishes checking that the noncommon factors Li do not vanish simulta-
neously in R3 in new families, different from the ones detailed in the statement. We



18 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

will show this fact computing some crossing resultants between them and proving
that the unique new possible intersection points have complex coordinates.

We start removing the parameter a1 computing the crossing resultants with re-
spect to the first factor: Rk = res(L3,Lk+3, a1), for k = 1, 2, 3, are polynomials in
(a5, a11) of degrees 30, 40, and 50, respectively, which decompose in some irreductible
factors with some natural powers (multiplicity). As the powers are not essential for

solving the equations we can remove all, and defining R̂i from Ri, having the same
factors but with multiplicity one. They have the common factor

L̂b = (1200a25 − 1110a5a11 + 225a211 − 169)(2a5 − a11)(11a5 − 9a11)

(8a5 − 5a11)(5a5 + a11)(15a5 − 11a11)(106a5 − 67a11).

The analysis of each new factor in L̂b follows analogously to the ones described
before, 15a5− 11a11 = 0 or 106a5− 67a11 = 0. When 2a5− a11 = 0, 11a5− 9a11 = 0,
8a5−5a11 = 0, or 5a5+a11 = 0, we always obtain Family (c). The condition 1200a25−
1110a5a11 + 225a211− 169 = 0 needs a more accurate analysis. Assuming it, the new

real solutions of SL should satisfy a21 = (±19712a5
√
a25 + 1 − 20288a25 − 7744)/225

and a11 = (37a5 ± 13
√
a25 + 1)/15, but as ±19712a5

√
a25 + 1− 20288a25 − 7744 ≤ 0

for each a5 we have no new center families because, additionally, we are assuming
a1 6= 0.

Hence, from now on we can assume L̂b 6= 0 and we can take R̃i = R̂i/L̂b. Removing
the parameter a5 from the resultants R12 = res(R̃1, R̃2, a5) and R13 = res(R̃1, R̃3, a5),
we obtain two polynomials on a11 with degrees 353 and 479, respectively. As above
we define R̂12 and R̂13 leaving only one factor if there are multiplicity higher than
one after the irreductible factorization in the rationals field. Here, the common
factor is

L̂c = a11(10443a211 + 6875)(48a211 + 1331)(1375a211 − 81)(5125a211 − 507)(15a211 − 1).

As only a11 = 0 gets a real solution of SL, providing Family (c), we can assume that

L̂c 6= 0.
The proof finishes because the last resultant R123 = res(R̂12/L̂c, R̂13/L̂c, a11) is a

non vanishing rational number. �

In an attempt to conclude the general case, i.e, to prove that without fixing w
the only centers families are the ones stated in Theorem 3 we have followed the
same scheme than in the previous proof but adding the next Lyapunov constant,
L7, of (6). The first steps follows easily. From L1 = 0 and L2 = 0, we write a8
and a7 as rational functions on (a1, a5, a11, w) and we can define the corresponding
Li(a1, a5, a11, w) functions, after removing all the common factors as in (20). The
first crossing resultants, Rk = res(L3,Lk+3, a1), for k = 1, . . . , 4, are big polynomials
with rational coefficients on (a5, a11, w) of degrees 119, 163, 207, and 251, and the
next R1k = res(R̃1, R̃k, a5), for k = 2, 3, 4 are even big polynomials of degrees 2454,
3438, and 4422, respectively, on (a11, w). Unfortunately, we are unable to conclude
the proof of the general case because we can not get the next crossing resultants
(removing a11) R123, R124 and the last (removing w) R1234, due to the huge memory
requirements to do all these computations.

4.2. Bifurcation of limit cycles of small amplitude. As in many problems of
degenerate Hopf bifurcation, the cyclicity of centers is usually less than the cyclicity



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 19

of some special systems having a weak focus of very high-order. This is also the case
in our families in X2. Theorem 4 provides a weak focus of order four but we have
no such points near the centers in Theorem 3. We have studied the local cyclicity
of some of them and next proposition provides our highest bifurcation result (near
centers) studying Family (e) in Theorem 3. The other families have even less local
cyclicity. We remark that we have also considered the perturbation of some centers
in Theorem 3 without choosing the conditions detailed in Section 4.1 but the local
cyclicity is not higher than the presented in the following results. So we have decided
to do not to write such results here.

Proposition 18. Consider the coefficients of (4) satisfying (e) in Theorem 3. Then,
there exists a perturbation in X2 such that at least 3 small amplitude limit cycles
bifurcate from the equilibrium point (0, 1, 0) on S2

1.

Proof. We take the parameter values (a1, a5, a7, a8, a11,w) satisfying Theorem 3.(e)
and we consider (a1, a5, a7, a8, a11, w) = (a1+ε1, a5+ε2, a7+ε3, a8+ε4, a11+ε5,w+ε6)
in the projected system (17). We denote by Li(ε), with ε = (ε1, . . . , ε6), the corre-
sponding Lyapunov constants, computed with the method explained on Section 2.3.
Clearly, when ε = 0 the origin is of center type, so Li(0) = 0 for all i. Then we

write the Taylor series of first-order with respect to ε as Li(ε) = L
[1]
i (ε)+O2(ε). The

proof follows because the linear terms of the first three Lyapunov constants have
rank 3 with respect to ε and, adding the trace parameter and using Theorem 8 or the
Implicit Function Theorem, we can get 3 limit cycles of small amplitude bifurcating
from the origin. �

Remark 19. We notice that we have no more limit cycles up to first-order study be-
cause the rank does not increase considering more Lyapunov constants. The second-
order study does not generate more either.

Before proving our last main result, we observe that the system (5) comes taking
a1 = −2α, a4 = 0, a5 = 1, a7 = 29/20, a8 = α, a11 = 2, and w = 2 in (4).

Proof of Theorem 4. As previously, we will work with the projected system (17)
corresponding to system (5). With the algorithm described in Section 2.3, the
Lyapunov constants are L1 = L2 = 0, L3 = 248832α(488α2 − 857)/5 and L4 =
−5971968α(361920α4 + 513328α2 − 3529237)/5. As L3 = 0 if, and only if, α = 0 or

α = ±
√

857/488, the proof of the weak focus order and center statements is clear.
When α = 0 we have a center since we are in Family (b) of Theorem 3. In the

other case, when α = ±
√

857/488, the origin is a weak focus of order 4 because
L4 = 6717957338234880α/3721 6= 0. For other values of α, the weak focus order is
only 3.

The cyclicity statement part follows considering the partial perturbation a1 =
−2α+ε1, a7 = 29/20+ε2, and a8 = α+ε3. The other parameters remain unchanged.
Then we compute the Taylor series of the Lyapunov constants up to first-order,

writing Li(ε) = L
[1]
i (ε) +O2(ε) with

L
[1]
1 (ε) = −7ε1 − 14ε3,

L
[1]
2 (ε) =

48

5
(2920α2 + 3809)ε1 − 51840αε2 +

24

5
(2920α2 + 3809)ε3,

L
[1]
3 (ε) =

248832

5
α(488α2 − 857)− 72

5
(3886080α4 − 1845344α2 − 380679)ε1



20 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

+ 62208α(2272α2 + 1909)ε2 −
36

5
(3886080α4 + 3214240α2 − 3342471)ε3.

When α = ±
√

857/488 we have proved above that L
[1]
3 (0) = 0. So, as the 3 × 3

matrix formed with the coefficients of (L
[1]
1 , L

[1]
2 , L

[1]
3 ) with respect to (ε1, ε2, ε3) has

nonzero determinant, we get (adding the trace parameter and using the Implicit
Function Theorem) four limit cycles bifurcating from the origin. When α /∈ A, as

L
[1]
3 (0) 6= 0, only three limit cycles can bifurcate from the origin. The proof finishes

just taking ε3 = 0 and checking that the matrix formed with the coefficients of

L
[1]
1 , L

[1]
2 with respect to (ε1, ε2) has rank two. �

4.3. Cyclicity on complex systems. Usually the cyclicity problem is considered
when the coefficients of the vector field are real numbers. But as, in fact, it depends
on the study of the zeros of a polynomial, this problem can be studied also when
the values of the parameters are complex numbers. In fact the existence of complex
solutions increase the difficulties in finding the center classification and, sometimes,
we need to go further in the Lyapunov constants because of this fact. Of course, they
do not provide neither real weak foci nor real centers. As we have shown in almost
all the proofs, there are complex solutions such that the orders of the complex weak
foci are higher than the associated to real solutions and some difficulties can appear
in the discussions. Moreover, the unfoldings can have more hyperbolic zeros that
could be considered as complex limit cycles. The aim of the last result is to describe
this phenomenon in our problem, showing why a carefully study is very important to
distinguish the existence or not of (real) centers, as for example in Proposition 17, or
which is the highest (real) weak focus order. In particular, we show the existence of
six (complex) limit cycles near the origin for family (e) in Theorem 3, instead of the
three that appear in Proposition 18 when we consider the problem in the reals. We
notice that we have not found any point in the corresponding (real) center variety
having more than three (real) limit cycles. In fact these three limit cycles have been
appeared generically up to a first-order analysis. Instead of study higher order,
we use the technique detailed in [11] that explain that the local cyclicity changes
moving parameters inside the center variety. This technique allow us to increase the
number of (complex) limit cycles but not the real ones.

Proposition 20. Consider system (4) with a1 =
w2 − 1

w2 + 1
a8, a5 =

w2 + 1

w2 − 1
a11 and

a7 =
1

w2 + 1
− 1

(w2 + 1)
a28−

w2 + 1

(w2 − 1)2
a211 for every fixed (a8, a11, w) such that w2 6= 1,

and a8a11 6= 0. Then it has a center at p = (0, 1, 0) and there exist (â8, â11, ŵ) ∈ C3

and complex perturbations inside the class X2 such that at least 6 (complex) small
amplitude limit cycles bifurcate from p.

Proof. The center property follows directly from Theorem 3.(e) because the proof
does not change if we consider complex coefficients instead of reals.

As the previous proofs we will work with the projected (perturbed) system (17)
corresponding to system (4) but with (a1+ε1, a5+ε2, a7+ε3, a8+ε4, a11+ε5, w+ε6)
instead of (a1, a5, a7, a8, a11, w). The algorithm explained in Section 2.3 allow us to
compute the first 6 Lyapunov constants and the first-order Taylor series write them

as Lk(ε) = L
[1]
k (ε) +O2(ε), for k = 1, . . . , 6, being ε = (ε1, . . . , ε6). We notice that,

from the center property, Lk(0) = L
[1]
k (0) = 0. Straightforward computations show



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 21

that the 3 × 6 matrix generated by the coefficients of (L
[1]
1 (ε), L

[1]
2 (ε), L

[1]
3 (ε)) with

respect to ε has rank 3. In fact the 6 × 6 matrix obtained from the first 6 has also
rank three. Hence, adding the trace parameter, a first-order study only provides
(generically) 3 small limit cycles using Theorem 8. Now we will use in a different
way the Implicit Function Theorem to prove the statement. We will closely follow
the scheme developed in [11].

Considering only the first two, instead of the first three, we can make an analytic
change of coordinates in the parameter space, changing (ε1, ε2) to (u1, u2) in order
that L1 = u1 and L2 = u2. Then we can use these new coordinates to continue our
analysis under the condition u1 = u2 = 0. Under this assumption we can obtain

L
[1]
k =

a8a11Mk(a8, a11, w)

w4k−7(w2 + 1)2k−3(w2 − 1)2k−2N(a8, a11, w)
U3(a8, a11, w, ε3, ε4, ε5, ε6),

for k = 3, . . . , 6, being

U3 = −(w2 + 1)2(w2 − 1)3ε3 − 2a8(w
2 + 1)(w2 − 1)3ε4 − 2a11(w

2 − 1)(w2 + 1)3ε5

+ 2w((w2 + 3)(w2 + 1)2a211 + (w2 − 1)3a28 − (w2 − 1)3)ε6,

N = (2(6w6 + 25w4 + 20w2 + 5))(w2 − 1)3a28 + 2w2(10w8 + 49w6 + 81w4 + 43w2

+ 9)(w2 + 1)2a211 − 3w2(w2 − 1)3(w2 + 1)4

and Mk polynomials with rational coefficients of degrees 30, 46, 62, 78 with 74, 185,
369, 640 monomials, respectively. We only show the first one because of the size of
them.

M3 = −(3w12 + 30w10 − 51w8 + 420w6 + 749w4 + 350w2 + 35)(w2 + 1)4a48

− (35w16 + 350w14 + 752w12 + 450w10 − 102w8 + 450w6 + 752w4 + 350w2

+ 35)(w4 − 1)2a28a
2
11 − w4(35w12 + 350w10 + 749w8 + 420w6 − 51w4

+ 30w2 + 3)(w2 + 1)4a411 − (6w16 + 60w14 − 3w12 + 390w10 + 1799w8

+ 2240w6 + 1267w4 + 350w2 + 35)(w2 − 1)4a28 − w4(35w16 + 350w14

+ 1267w12 + 2240w10 + 1799w8 + 390w6 − 3w4 + 60w2 + 6)(w4 − 1)2a211

− 3w4(w4 + 8w2 + 1)(w2 + 1)4(w2 − 1)6.

We notice that, as a8a11 6= 0 and w2 6= 1, we can change ε3 by u3 isolating from
u3 = U3. Clearly when M3 is nonzero we have rank 3 and no more limit cycles, up
to first-order exist, because the next Lk vanish if L3 does.

The proof ends checking that there exist values η = (â8, â11, ŵ) ∈ C3 such that
M3(η) = M4(η) = M5(η) = 0, M6(η) 6= 0, N(η) 6= 0, and the determinant of
the Jacobian matrix, J , of (M3,M4,M5) does not vanish at η. Hence, for each
transversal solution η, there exists a small neighborhood such that the Implicit
Function Theorem allow us to write, for ε4 = ε5 = ε6 = 0,

Lk = v̂ku3 +
∞∑
i=2

L̃k,i(v̂3, v̂4, v̂5)u
i
3, for k = 3, 4, 5,

and

L6 = ṽ6u3 +
∞∑
i=2

L̃6,i(v̂3, v̂4, v̂5)u
i
3,



22 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

with ṽ6 ∈ C\{0} and L̃k,i analytic functions. The next step needs again the Implicit
Function Theorem, dividing each equation by u3, to write Lk = vku3 for k = 3, 4, 5.
Finally, for u3 small enough, we have a (complex) weak focus of order 6 that unfolds
6 (complex) limit cycles of small amplitude. We must use the trace parameter
together with u1, u2 to have the complete versal unfolding.

We will finish now showing how the transversal intersection together with the
nonzero conditions hold. From the scheme used in the proof of Proposition 17, we
start computing the crossing resultantsR4 = res(M3,M4, a8) andR5 = res(M3,M5, a8)
and after removing the common factors we get that one of the factors of R45 =
res(R̂4, R̂5, a11) is the polynomial

p(w) = 2358125w36 + 21253750w34 + 74700325w32 + 418209680w30

+ 2696915172w28 + 11471299432w26 + 32469933620w24

+ 65296240368w22 + 97631971158w20 + 111336439620w18

+ 97631971158w16 + 65296240368w14 + 32469933620w12

+ 11471299432w10 + 2696915172w8 + 418209680w6

+ 74700325w4 + 21253750w2 + 2358125,

which has degree 18 in w2 and has no real roots. Straightforward computations
provide that the solutions of

{M3(a8, a11, w) = M4(a8, a11, w) = M5(a8, a11, w) = p(w) = 0}
can be written as η = (â8, â11, ŵ) = (γ/393216, β/393216, α) with p(α) = 0 and

γ2 = (−2109840270437897288946799375α34 − 17713873005809414028007091875α32

− 55845814488941782868405995600α30 − 339226787759422100359736546240α28

− 2202132815857798644636420812396α26 − 8895632031878776983832830622252α24

− 23501844430016397280367284455712α22 − 43681378432973931255571646316816α20

− 59799471686092166179474533669890α18 − 61675400652851000448826407685450α16

− 47993652204362508680940846437904α14 − 27604081780746071222684026969568α12

− 11207519567177779180170267397628α10 − 2966440603081008995602043288284α8

− 468074157329937609377913199360α6 − 66069795390963180008678824400α4

− 24129993592921092989446394375α2 − 3181668269212383369855156875)/9396192821225170,

β2 = (40248036412021566334043770625α34 + 348148658927405932769817170625α32

+ 1139390577520977495682159181200α30 + 6652155763250716876027785442080α28

+ 43419594803251066283846291084532α26 + 178700034808562373265475821347460α24

+ 480531845343567179492990058264000α22 + 906984466960916743690778969263248α20

+ 1258033583163550539675725186999550α18 + 1313015179221532924371751690040910α16

+ 1033425110388324368012929243720848α14 + 601374377745001674638591851302720α12

+ 247291049686230857389940010207140α10 + 66539552467687537839070855101012α8

+ 10700521536131995853057121609120α6 + 1483387556836472011860839676560α4

+ 538998995971563271556998588905α2 + 69737060750810475560536625305)/32886674874288095.

Using the above relations we can check that, for each α, simple root of the polyno-
mial p, we have that M3(η),M4(η), and M5(η) vanish and M6(η) = p1(α), N(η) =
p2(α), and J (η) = αβγp3(α) do not, because the polynomials pi(α), for i = 1, 2, 3,



CENTERS AND LIMIT CYCLES ON INVARIANT SPHERES 23

which are of degree 17 in α2 and with rational coefficients, have nonzero resultants
with the polynomial p.

�

We remark that, in the above proof, the coefficients of M3 with respect to (a8, a11)
are negative for all w ∈ R, so the described bifurcation mechanism does not provide
more than the three limit cycles in Proposition 18.

Acknowledgments

This work has been realized thanks to the Catalonia AGAUR 2017 SGR 1617
grant; the Spanish Ministerio de Ciencia, Innovación y Universidades - Agencia Es-
tatal de Investigación PID2019-104658GB-I00 grant; the Brazilian CNPq 304798/2019-
3 and São Paulo Research Foundation (FAPESP) 2017/08779-8, 2019/00440-7 and
2019/10269-3 grants; and the Coordenação de Aperfeiçoamento de Pessoal de Nı́vel
Superior – Brasil (CAPES) – Finance Code 001; and the European Community
H2020-MSCA-RISE-2017-777911 grant.

References

[1] M. J. Ablowitz, A. Ramani, and H. Segur. A connection between nonlinear evolution equations
and ordinary differential equations of P -type. I. J. Math. Phys., 21(4):715–721, 1980.

[2] A. A. Andronov, E. A. Leontovich, I. I. Gordon, and A. G. Măıer. Qualitative theory of second-
order dynamic systems. Halsted Press (A division of John Wiley & Sons), New York-Toronto,
Ont.; Israel Program for Scientific Translations, Jerusalem-London, 1973.

[3] N. N. Bautin. On the number of limit cycles which appear with the variation of coefficients
from an equilibrium position of focus or center type. Number 100. American Mathematical
Society, 1954.

[4] L. R. Berrone and H. Giacomini. Inverse Jacobi multipliers. Rend. Circ. Mat. Palermo (2),
52(1):77–130, 2003.

[5] M. Caubergh, J. Llibre, and J. Torregrosa. Global phase portraits of some reversible cubic
centers with collinear or infinitely many singularities. International Journal of Bifurcation and
Chaos, 22(11):1250273, 2012.

[6] M. Caubergh and J. Torregrosa. Global phase portraits of some reversible cubic centers with
noncollinear singularities. Internat. J. Bifur. Chaos Appl. Sci. Engrg., 23(9):1350161, 30, 2013.

[7] C. Chicone and M. Jacobs. Bifurcation of critical periods for plane vector fields. Trans. Amer.
Math. Soc., 312(2):433–486, 1989.

[8] C. Christopher. Estimating limit cycle bifurcations from centers. In Differential equations with
symbolic computation, pages 23–35. Springer, 2005.

[9] A. Cima, A. Gasull, V. Mañosa, and F. Mañosas. Algebraic properties of the Liapunov and
period constants. Rocky Mountain J. Math., 27(2):471–501, 1997.

[10] F. Dumortier, J. Llibre, and J. C. Artés. Qualitative theory of planar differential systems.
Universitext. Springer-Verlag, Berlin, 2006.

[11] J. Giné, L. F. S. Gouveia, and J. Torregrosa. Lower bounds for the local cyclicity for families
of centers. J. Differential Equations, 275:309–331, 2021.

[12] A. Goriely. Integrability and nonintegrability of dynamical systems, volume 19 of Advanced
Series in Nonlinear Dynamics. World Scientific Publishing Co., Inc., River Edge, NJ, 2001.

[13] L. F. S. Gouveia and J. Torregrosa. Lower bounds for the local cyclicity of centers using high
order developments and parallelization. J. Differential Equations, 271:447–479, 2021.

[14] M. Han. Liapunov constants and Hopf cyclicity of Liénard systems. Ann. Differential Equa-
tions, 15(2):113–126, 1999.

[15] Y. Ilyashenko. Centennial history of Hilbert’s 16th problem. Bulletin of the American Math-
ematical Society, 39(3):301–354, 2002.

[16] J. S. W. Lamb and J. A. G. Roberts. Time-reversal symmetry in dynamical systems: a survey.
Phys. D, 112(1-2):1–39, 1998. Time-reversal symmetry in dynamical systems (Coventry, 1996).



24 C. BUZZI, A. L. RODERO, AND J. TORREGROSA

[17] J. Llibre and C. Pessoa. Homogeneous polynomial vector fields of degree 2 on the 2-dimensional
sphere. Extracta Math., 21(2):167–190, 2006.

[18] J. Llibre and C. Pessoa. Invariant circles for homogeneous polynomial vector fields on the
2-dimensional sphere. Rend. Circ. Mat. Palermo (2), 55(1):63–81, 2006.

[19] P. J. Olver. Applications of Lie groups to differential equations, volume 107 of Graduate Texts
in Mathematics. Springer-Verlag, New York, 1986.

[20] V. G. Romanovski and D. S. Shafer. The center and cyclicity problems: a computational
algebra approach. Birkhäuser Boston, Ltd., Boston, MA, 2009.

[21] R. Roussarie. Bifurcation of planar vector fields and Hilbert’s sixteenth problem, volume 164
of Progress in Mathematics. Birkhäuser Verlag, Basel, 1998.

Universidade Estadual Paulista (UNESP), Instituto de Biociências Letras e Ciências
Exatas, São José do Rio Preto-SP, 15054-000, Brazil

Email address: claudio.buzzi@unesp.br

Universidade Estadual Paulista (UNESP), Instituto de Biociências Letras e Ciências
Exatas, São José do Rio Preto-SP, 15054-000, Brazil

Email address: ana.rodero@unesp.br

Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Be-
llaterra, Barcelona (Spain); Centre de Recerca Matemàtica, Campus de Bella-
terra, 08193 Bellaterra, Barcelona (Spain)

Email address: torre@mat.uab.cat


	1. Introduction
	2. Preliminary results
	2.1. Setting the problem
	2.2. Darboux Integrability
	2.3. Lyapunov constants and local cyclicity

	3. Centers for linear and quadratic homogeneous vector fields
	3.1. Linear case
	3.2. Quadratic homogeneous case

	4. Centers and cyclicity for quadratic vector fields
	4.1. Center characterizations
	4.2. Bifurcation of limit cycles of small amplitude
	4.3. Cyclicity on complex systems

	Acknowledgments
	References