Cyril Gavoille

Claire Hilaire
email: claire.hilaire@labri.fr

Minor-Universal Graph for Graphs on Surfaces

Keywords: Minor, Universal Graphs, Bounded Genus Graphs

published or not. The documents may come   L'archive ouverte pluridisciplinaire

Introduction

Graph universality plays an important role in graph theory and in algorithmics. Rather than considering individually each graph of a given familly, it is usually simpler to manipulate a reduced set of graphs, or even a single graph, whose properties are close enough from those of the graphs of the family.

For instance, the Strong Product Theorem [DJM + 20], which states that every planar graph is a subgraph of the strong product of a graph of bounded treewidth and a path, has provided many new insights in planar graph theory, as it suffices to study such product to derive properties on all planar graphs. This includes a proof for bounding their queue number and other variant of chromatic number [DEJ + 20]. Another example is the design of NCA-labeling schemes [START_REF] Alstrup | Nearest common ancestors: A survey and a new algorithm for a distributed environment[END_REF], whose goal is to associate with each vertex of each n-vertex tree a short label so that the label of the nearest-common ancestor between any two vertices u, v can be retrieved from the labels of u and of v. Interestingly, it is shown in [GKŁ + 18] how a specific tree with n c vertices containing every n-vertex trees as topological minor can be used for the design of succinct NCA-labeling schemes, namely with c log nbit labels. Similar connections are known between adjacency-labeling schemes and induceduniversal graphs [ADBTK17, DEG + 21]. The last example is the refinement of celebrated Excluded Grid Theorem of Robertson and Seymour [START_REF] Robertson | Quickly excluding a planar graph[END_REF] which states that the graphs excluding as minor a fixed planar graph have bounded treewidth. This result relies on the fact that every planar graph is a minor of some relatively small grid. So excluding such a grid as minor implies excluding as minor any small enough planar graph.

Our results. Motivated by the latter example, we focus in this paper on universality under minor1 containment. More precisely, a graph U is minor-universal for a given graph family F if every graph of F is a minor of U. Keeping in mind that we want the graph U to be representative of some properties of graphs in F, most of the time we will impose that U and the graphs in F share a common property, say, for instance, planarity. Note that the Strong Product Theorem fails for this purpose as it provides for planar graphs a universal graph which is certainly of n2 vertices but far from being planar in general. It can be shown 2 , for instance, that the strong product of two paths of n vertices is already of genus Ω(n).

From [RST94, Theorem (1.3) & (1.4)], the 2n × 2n-grid contains as minor every planar graph with n vertices. In other words, planar n-graphs have a minor-universal grid with at most 4n 2 vertices. This quadratic bound is essentially optimal, as long as grids are concerned, since [BCE + 19, Lemma 5] showed that any minor-universal grid for planar graphs must have Ω(n 2 ) vertices.

We generalize the result of [START_REF] Robertson | Quickly excluding a planar graph[END_REF] to graphs on surface as follows:

Theorem 1 For every n and every surface Σ of Euler genus g, there is a graph U n,Σ embedded on Σ with O(g 2 (n + g) 2 ) vertices that contains as minor every graph embeddable on Σ with n vertices.

In other words, U n,Σ is a minor-universal graph for n-vertex graphs embeddable on Σ. In fact, as we will see in Theorem 2, the graph U n,Σ depends only on the Euler genus and orientability of Σ. Moreover, the embedding of the minor is preserved in U n,Σ .

Our proof is constructive. The minor-universal U n,Σ , as well as the witness of the minor in U n,Σ , can be constructed in polynomial time. Interestingly, along the way, it gives an alternative and self-contained proof of [START_REF] Robertson | Quickly excluding a planar graph[END_REF] for the minor-universal grid of planar graphs.

Related results. For trees, [START_REF] Bodini | On the minimum size of a contraction-universal tree[END_REF] proved that the smallest minor-universal tree for nvertex trees has between Ω(n log n) and O(n 1.985 ) vertices. In fact, the result of [GKŁ + 18] implies that this smallest minor-universal tree has between Ω(n 1.724 ) and O(n 1.895 ) vertices3 . From computational perspectives, the smallest minor-universal tree of two bounded degree trees can be computed in cubic time [START_REF] Nishimura | Finding smallest supertrees under minor containment[END_REF], whereas it is NP-complete to decide whether or not a tree is minor of another one. Minor universality has also been investigated for infinite graphs [DK99, HMŠ + 21]. However, as noticed by [START_REF] Diestel | A universal planar graph under the minor relation[END_REF], the problem already becomes different in the planar case since the infinite grid is no longer minor-universal for infinite planar graphs.

Preliminaries

We consider compact connected surfaces without boundary. A canonical system of loops for a surface Σ of Euler genus g is a one-vertex graph embedded on Σ with g loop-edges such that by cutting Σ along these loops results in a surface homeomorphic to a polygon P with 2g sides. The way the pairs of sides of P are reattached to form Σ can be described by a signature, a word associating clockwise one symbol with each side of P . The polygon P and its signature form a canonical polygonal schema. Note that the loops are either two-sided (the surface Σ is orientable with oriented genus g/2) or one-sided (Σ is non-orientable with non-oriented genus g).

For g > 0, the signature is a 1 a 2 ā1 ā2 . . . a g-1 a g āg-1 āg if Σ is orientable, and a 1 a 1 . . . a g a g if Σ is non-orientable. For the sphere, g = 0, the signature is a 0 ā0 by convention. Sides are reattached if they have the same (or opposite) symbol in the signature, given the convention that positive symbols correspond to a clockwise orientation of the sides, and negative symbols to an anticlockwise orientation.

Consider a graph G embedded4 on some surface Σ with a canonical system of loops L. Denote by G L the graph embedded on Σ obtained from the union of G and L, and by adding a new vertex at each intersection between G and L.

A polygonal embedding for G is the planar embedding Π into a polygon obtained from G L by cutting Σ along the edges of L. (We refer to Fig. 1 for an illustration.) The edges and vertices on L appear duplicated on the boundary of outerface of Π. More precisely, the boundary of the outerface of Π is a cycle that can be cut into a clockwise sequence of paths sharing their extremities. These paths, called sides of Π, correspond to the edges of the polygon P , and their extremities, called corners of Π, corresponds to the vertices of P . Note that corners have degree two in Π. The ordered sequence of corners is called the border of Π. The signature of Π, denoted by σ(Π), describes how the sides of Π are reattached to form G L on Σ. More precisely, the ith symbol of the signature is associated with the ith sides of Π, i.e., the path between the ith and (i + 1)th corners. These sides are merged according to the orientation given by the symbols, and there are called twin sides of Π. Obviously, twin sides must contain the same number of vertices. Futhermore, the vertices and edges that are identified in this process are called respectively twin vertices and twin edges.

Given a polygonal embedding Π, the sewing of Π, denoted by sewing(Π), is the graph obtained by reattaching the sides of Π according to its border and its signature. It is clear that, if Π comes from G L embedded on Σ by cutting along the edges of L, then sewing(Π) is isomorphic to G L. In particular, sewing(Π) embeds on Σ and contains G as minor.

Main result

To prove Theorem 1, we consider a graph G embedded on the surface Σ of Euler genus g and having a canonical systems of loops L with vertex r and signature σ. Thus G has a polygonal embedding Π G such that G is a minor of sewing(Π G ), a graph that embeds on Σ. Note that |V (Π G )| 2|V (G L)| + |σ| since every vertex of G L appears at most on two sides of Π G , and each corner is a copy of r. In fact, Π G has at most n internal vertices due to G, and a

r 1 2 a 1 ā1 a 2 ā2
Figure 1: The graph K 6 embedded on the torus (of Euler genus 2) with the system of loops L = ({r} , { 1 , 2 }), and a polygonal embedding Π of K 6 of signature σ = a 1 a 2 ā1 ā2 . By construction, K 6 is a minor of the graph sewing(Π) = K 6 L.

certain number of vertices lying on the sides of Π G due to r and to the intersections between G and L. We can bound this number thanks to the following result5 : Lemma 1 ( [START_REF] Lazarus | Computing a canonical polygonal schema of an orientable triangulated surface[END_REF][START_REF] Fuladi | Short topological decompositions of non-orientable surfaces[END_REF]) Given an embedding of G on a surface of Euler genus g, there is a polynomial time algorithm that computes a canonical system of g loops such that each loop intersects any edge of G in at most6 30 points.

It follows that |V (Π G )| = O(|V (G)| + g|E(G)|) = O(g(n + g))
since each of the g loops of L crosses O(1) times each of the O(n + g) edges7 of G. Our strategy then is to transform step by step the initial polygonal embedding Π G into some unique polygonal embedding U that preserves signature and minor containment of its sewing. The number of vertices of this final embedding is |V (U)| = O(|V (Π G )| 2 ). Our minor-universal graph U n,Σ as in Theorem 1 is nothing else than sewing(U) that has less than |V (U)| = O(g 2 (n + g) 2 ) vertices as claimed.

To formalize signature and minor containment preserving, we introduce the following relation. A polygonal embedding Π is a p-minor of a polygonal embedding Π if they have the same signature and if sewing(Π) is a minor of sewing(Π ). Moreover, we say that Π has size (m, n) if each side has at most m vertices, the corners excluded, and at most n internal vertices (those that are not lying on the boundary of the outerface of Π). In particular, every polygonal embedding Π of size (m, n) has at most n + |σ(Π)|(m + 1) vertices8 , and sewing(Π) has at most n + |σ(Π)|(m + 1)/2 vertices (this is n + |σ|m/2 + 1 vertices if σ(Π) is canonical, see Section 4.3). A polygonal embedding of size of (m, 0) is a Hamiltonian outerplane graph.

The final polygonal embedding is denoted hereafter by U σ,m , where σ is any signature and m an integral parameter. We will see later in Section 4.3 that U σ,m has signature σ and size (m, O((|σ|m) 2 )). Roughly speaking, it is a square half-grid cut by its diagonal. The vertices on the boundary of the outerface are those of the diagonal, see Fig. 2 for an example. The technical theorem is the following. Note that it applies to any polygonal embedding, not only those of canonical signature.

Theorem 2 Every polygonal embedding Π of size (m, n) is a p-minor of U σ(Π),m+2n .

Before proving Theorem 2, let us show that it implies Theorem 1. Consider a graph G on n vertices enbedded on a surface Σ of Euler genus g. By Lemma 1, G has a polygonal embedding Π G whose signature σ depends only on Σ (as it corresponds to that of the canonical polygonal schema of Σ), and such that

Π G has size (m, n) with m = O(n + g). Moreover, G is a minor 9 of sewing(Π G ). By Theorem 2, Π G is a p-minor of U σ,m+2n
, which by definition of p-minor implies that sewing(Π G ) is a minor of sewing(U σ,m+2n ). It follows that G is a minor of sewing(U σ,m+2n ) = U n,Σ , that has the required properties: it is embedded on Σ and has O(g 2 (n + g) 2 ) vertices.

To prove Theorem 2, we proceed in two steps, which are summarized and formalized by the next two lemmas whose proofs are sketched hereafter and fully presented in Section 4:

Lemma 2 Every polygonal embedding Π of size (m, n) is a p-minor of a polygonal embedding of size (m + 2n, 0).
Intuition is that all internal vertices of Π can all be somehow pushed to the sides of the outerface, without increasing much the number of vertices on each side. It follows we can consider that all the vertices belong to the outerface of the polygonal embedding making it outerplanar. This is done by first cutting along the edges of some spanning forest of the inner graph (rooted in a vertex of the outerface), creating an empty space and transforming the graph into an outerplanar graph. Along this process, each vertex that is not a leaf appears more than once (the number of new vertices is bounded by the number of edges of the forest, thus by n). To keep track of this vertex duplication, we add some edges that could be contracted later to get back to the original polygonal embedding. Other transformations are then needed to embed those new edges and obtain an outerplanar polygonal embedding of the right size and containing Π as a p-minor.

Lemma 3 Every polygonal embedding Π of size (m, 0) is a p-minor of U σ(Π),m .
Note that Theorem 2 is a straightforward combination of Lemma 2 and Lemma 3.

The proof of Lemma 3 is inspired from the embedding into a n × n-grid for hamiltonian planar graphs with n vertices [RST94, Theorem (1.3)]. From this embedding, one can easily show that outerplanar graphs with n vertices are minor of a half-n × n-grid. In our context, the diagonal of this half-grid is the place for the sides of our polygonal embedding, as depicted in Fig. 2. Interestingly, along the way, we give in Proposition 1 an alternative prove of a result of [RST94, Theorem (1.4)] that states that every planar graph with n vertices is minor of a Hamiltonian planar graph with 2n vertices.

Proofs of the lemmas

Along the proofs, if H is a minor of G, it will be sometimes more convenient to say that G is a major of H. Similarly, we say that polygonal embedding

Π G is a p-major of Π G , if Π H is a p-minor of Π G .
Because plane graphs we will consider are actually polygonal embeddings, the boundary of the outerface is a cycle. Whenever we talk about outerface we refer to this cycle. We say that two plane graphs H, G such that H is minor of G have the same outerface if the minor H has been obtained from G without removing vertices or edges of the outerface of G, nor contracting edges with both endpoints on the outerface of G. It follows:

Property 1 If a polygonal embedding Π H is a minor of polygonal embedding Π G with same outerface, border and signature, then Π H is a p-minor of Π G .

Recall that a polygonal embedding is nothing else than a plane graph with a given cycle outerface, border and signature. Property 1 allows us to consider that we work on plane graph, as long as we keep the same outerface and border when constructing a major.

In Section 4.1, we consider plane graphs with a fixed cycle outerface, i.e., with an outerface whose boundary is a cycle. We give tools to construct a plane major that preserves this outerface. These tools will be used later in the construction of the p-major of Π in Section 4.2.

Considering plane graphs

Let G be a planar graph G with a given planar embedding and having a cycle outerface O. We denote by

n = |V (G)| -|V (O)| the number of internal vertices of G.
Because the boundary of the outerface of G is a connected graph (a cycle), we can assume w.l.o.g. that G is connected as well. If not, we can for instance triangulate all the faces except the outerface. By doing this we get a plane major of G without altering its number of vertices nor its cycle outerface.

Let F be a forest of G, composed of k trees denoted by T 1 , . . . , T k , and constructed as follows:

• Start with a spanning tree T of G;

• Remove one edge of any path in T that connects two vertices of O, if any;

• Remove all remaining isolated vertices that are in O; and • Root each T i at the only one vertex of T i ∩ O, denoted by r i .

Observe that F spans all internal vertices of G, thus |V (F )| = n + k and |E(F )| = n. Because we consider planar major of G, it is quite tempting to add edges to force F to be a single tree. However, we cannot do that and suppose that the inner vertices induced a connected subgraph of G, as there could be an edge between two vertices of the outerface separating the internal vertices into two parts.

We now construct a major of G, denoted by G 1 , that consists in blowing up each tree T i into a graph T i = C i ∪ E i , composed of a cycle C i plus a set E i of extra edges. More precisely, blowing up T i consists in traversing the tree according to a plane Euler tour from r i (or in other words, a walk along the boundary of the outerface of T i considered here as a single plane graph). The cycle C i is constructed iteratively by adding a new vertex at each vertex of T i visited along this tour, two consecutive vertices on this tour being connected by an edge. See Fig. 3 for an illustration. Since F has no isolated vertices, T i has at least one edge. If T i has exactly one edge, then C i consists of one single edge. The set E i is composed of all the edges connecting any two vertices u, v of C i if they corresponds to the same visited vertex of T i and that, among them, appear consecutive during the visit. In other words, we add an edge between u, v in C i if they corresponds to same visited vertex w of T i and that none of the vertices between u and v in C i corresponds to w.

This completes the description of blowing up

T i into T i . E i a i r i C i O G T i G 1 T i
Figure 3: Blowing up the tree T i (with black edges) into T i = C i ∪ E i (with red and violet edges). The root r i of T i becomes an anchor

a i of C i in T i .
The graph G 1 is then obtained from G where each T i is replaced by T i which is outerplanar. This is possible because the walk along the boundary of the outerface of T i and T i are isomorphic, so T i can be plugged into T i by keeping the plane embedding of G.

Claim 1 The graph G 1 has the following properties:

• G 1 is a plane major of G with same outerface O; • C 1 , . . . , C k form a partition of internal vertices of G 1 . • Each C i contains exactly one vertex in O, called its anchor; • G 1 has 2n -k internal vertices;
Proof. It is easy to check that T i is a major of T i , as it suffices to contract all edges of E i to obtain a graph isomorphic to T i . By doing this for every T i in G 1 we get exactly G (up to some isomorphism) since these are the only differences.

We have seen that the trees T i 's of F form a partition of the internal vertices of G. It follows that T i forms a partition of the internal vertices of G 1 , and C i as well since

V (T i ) = V (C i ).
Each T i intersects O in exactly one vertex, its anchor, the vertex corresponding to r i in T i . In particular, no edge of T i belongs to O. So contracting edges of E i to get T i cannot affect the outerface, and thus G 1 and G have the same outerface O.

To obtained C i from T i , each edge of T i is traversed twice.

It follows that

|V (C i )| = 2|E(T i )|. We have seen that |E(F )| = k i=1 |E(T i )| = n. Therefore, the number of internal vertices of G 1 is k i=1 (|V (C i )| -1) = 2n
-k as we need to remove the anchors.

We now construct from G 1 a major G 2 , that consists in applying a splitting anchor operation at every anchor of G 1 . More precisely, consider a cycle C i with anchor a i in G 1 . Let u, w 1 , w 2 , v be the neighbors of a i , taken in this cyclically ordered around a i , such that u, v belong to O and w 1 , w 2 to T i . By construction u, a i , v belongs to the boundary of a common face, as well as w 1 , a i , w 2 . The edge a i -w 1 belongs to C i , and a i -w 2 belongs to either C i or E i , depending whether r i has one or several children, and w 1 = w 2 if C i consists in one edge. The splitting operation consists in replacing a i by an edge a i -a i , called edge-anchor of a i , and reconnecting to a i all the a i 's neighbors going from u to w 1 , and reconnecting to a i all the a i 's neighbors going from w 2 to v. We denote by P i the path going from a i to a i in C i , and by F i the new set of edges E i in G 2 . Note that C i and P i differ by a i , a i , a i , and that E i and F i differ by a i in G 2 that plays the role of a i in G 1 . See Fig. 4 for an illustration.

u w 2 w 2 E i F i C i O G 1 G 2 P i a i u a i v a i v w 1 w 1 Figure 4: Splitting anchor a i of C i into the edge-anchor a i -a i . Claim 2 The graph G 2 is a Hamiltonian plane major of G 1 with |V (O)| + 2n vertices.
Proof. Clearly, G 2 is a plane major of G 1 , as it suffices to contract every edge-anchor of G 2 to obtained G 1 .

To construct a Hamiltonian cycle in G 2 , we start with the cycle outerface of G 2 . Then, each edge-anchor a i -a i is replaced by P i . Inherited from Claim 1 and from the fact that C i 's partition all internal vertices of G 1 , all the internal vertices of G 2 are (disjointly) spanned by the P i 's. Thus it forms an Hamiltonian cycle for G 2 .

The graph G 1 has |V (O)| vertices on its outerface O and 2n -k internal vertices (by Claim 1). Splitting each of the k anchors adds k vertices in G 2 . Therefore,

|V (G 2 )| = |V (G 1 )| + k = (|V (O)| + (2n -k)) + k = |V (O)| + 2n.
Note that the splitting anchor operation alters the outerface. So, for the proof of Lemma 2, that needs to consider polygonal embedding whose twin sides, in particular, have to be of same length, we will need to make extra transformation on the outerface as explained later in Section 4.2.

Before, we observe that the two previous transformations (by Claim 1 and Claim 2) provide an alternative proof of [RST94, Theorem (1.4)] that states that every n-vertex planar graph is a minor of a planar Hamiltonian graph with at most 2n vertices.

Recall that a circuit in a graph extends the notion of cycle to subgraph composed of one single vertex or edge. A circuit is separating if its deletion increases the number of connected component of the graph, and it is non-separating otherwise. It is not difficult to see that every graph has a non-separating circuit with at least one vertex.

Proposition 1 Every planar graph with n vertices and with a non-separating circuit of k vertices is minor of a Hamiltonian planar graph with at most 2n -k vertices. In particular, every triangulation with n 4 vertices is minor of a Hamiltonian planar graph with at most 2n -4 vertices.

Proof. Let G be a plane graph with n vertices and with a non-separating circuit of k vertices, denoted by O. W.l.o.g. we can assume that O is a cycle, i.e., k 3, since otherwise we can triangulated G and consider its triangle outerface as non-separating circuit. And, obviously the triangulation is a major of G and we will get a minor of 2n -3 vertices that is less than 2n -k if k < 3.

First, we embed G in the plane such that all the vertices of G \ O are inside O. Note that O is not necessarily the boundary of the outerface of this embedding, since there can be some chords outside of O. Denote by C the set of these chords. Let G be the plane graph obtained from G \ C such that O is the boundary of its outerface which is a cycle.

Since G is a plane graph with cycle outerface O, we can construct the major G 1 as in Claim 1, and the Hamiltonian plane major G 2 from G 1 (by Claim 2). Due to the edge-anchor, the outerface of G 2 is a subdivision of O. Therefore, we can add all the chords of C in G 2 that we can embedded outside the outerface of G 2 , while preserving its planarity. Thus, G 2 ∪ C is the desired major of G.

The graph G has n = n -k internal vertices by construction. So, by Claim 2, G 2 (and G 2 ∪ C as well) has k + 2n = 2n -k vertices as claimed.

If G is a triangulation with at least four vertices, then G has a non-separating circuit of length k + 1. Indeed, the outerface is a non-separating triangle that can be increased by one using any incident internal face. So the above property construction applies with k = 4.

Actually, the proof of [RST94, Theorem (1.4)] gives an upper bound of 2n -4 vertices (assuming n 4), and relies on Whitney's Theorem 10 [Whi31]. In contrast Proposition 1 gives an explicit and direct construction with the same upper bound. We believe that it gives an interesting construction of the major where the Hamiltonian cycle cuts the graph into two sides. Inside the cycle are the copies of edges that where initially in G, and outside of the cycle are the edges introduced by the operations of blowing up trees and splitting anchors: contracting those edges results in G.

Proof of Lemma 2

The goal of this section is to show that: Lemma 2 Every polygonal embedding Π of size (m, n) is a p-minor of a polygonal embedding of size (m + 2n, 0).

Let Π be a polygonal embedding of size (m, n). Up to taking a p-major of Π, we assume that each corner of Π is on the boundary of a inner triangle: if not, we can connect by an edge the two neighbors of each corner. This does not alter the size of Π. This is to prevent against any further edge insertion that could increase the degree of corners which must be exactly two in any polygonal embedding. The boundary of the outerface of Π is a cycle, so the construction of G 2 as in Claim 1 applies. Let Π 1 be the polygonal embedding obtained from Claim 1, with Π in the role of G, Π 1 in the role of G 1 . The embedding Π 1 admits Π as a minor while preserving the cycle outerface (border and signature does not matter for Claim 1). Thus Π is a p-minor of Π 1 . Keeping the same notation as in Claim 1, Π 1 has exactly 2n -k internal vertices that are partitioned by the cycles C 1 , . . . , C k , each C i intersecting the outerface in one vertex, its anchor a i .

Unfortunately, as said previously the transformation of G 1 into G 2 as in Claim 2 does not preserve a polygonal embedding: twin edges in G 1 are not twin anymore in G 2 because of the new created edge-anchor. So it breaks the p-major sequence we are constructing. To preserve a polygonal embedding, we need to modify accordingly twin sides whenever we apply a splitting anchor so that each vertex of one side has a twin vertex on its twin side.

For this purpose we define a twin splitting operation that applies to each anchor of Π 1 . Roughly speaking, we apply a splitting anchor operation while subdividing one of the two twin edges that are incident to the twin vertex of the anchor. The edge to subdivide depends on a total ordering ≺ defined on the vertices of the cycle outerface. It is based on their rank, starting with the very first corner of the border of Π 1 . So u ≺ v if vertex u is visited before v when traversing the cycle outerface from the first corner of the border.

Each step t of this twin process consists in applying one twin splitting operation which results into a new polygonal embedding Π t 1 obtained from Π t-1 1 , starting with Π 0 1 = Π 1 . More precisely, consider u -a -v be three consecutive vertices on a side of Π t-1 1 where a is an anchor and u ≺ a ≺ v. Let x, b, y be the twin vertices respectively of u, a, v. Because an anchor cannot be a corner, vertices u, v, x, y are well defined. The twin splitting operation consists in applying a splitting anchor operation on a followed by (cf. In the former case, the edge b -y is replaced by the path b -s -y and b -s becomes the twin edge of a -a , the edge-anchor of a. In the latter case, the edge x -b is replaced by the path x -s -b and s -b becomes the twin of a -a . All the other twin relations are unchanged. See Fig. 5 for an illustration. It is important to observe that b may be an anchor that will be eventually split at some step t > t. In this representation, the focus on the twin sides creates the illusion that the inside part of the embedding (in gray) is a disconnected region, which is not true. However, in this representation choice, the inside part, that is locally planar, may appear twisted if fully represented in the case of a non-orientable surface.

After this twin process, denote by Π 2 the final embedding obtained from Π 1 , i.e., the embedding obtained from Π 1 by applying successively all twin splitting operations.

Claim 3 The embedding Π 2 has the following properties:

• Π 2 is a polygonal embedding that is p-major of Π 1 ; • Π 2 has size (m + k, 2n -k); • P 1 , . . . , P k form a partition of the internal vertices of Π 2 ; • During the twin process, no edge-anchor is subdivided and no two edge-anchors are twins.

Proof. Each twin splitting operation inserts edges and split anchors on the outerface of Π 1 . So, by contracting these new edges, we get back to Π 1 which shows that Π 2 is a major of Π 1 . Moreover, thanks to the twin splitting, eventually every vertex (and edge) of a side in Π 2 of have a well defined twin. So, Π 2 is a polygonal embedding that is p-major of Π 1 .

At each twin splitting operation, the number of vertices of a side increases by at most one while the number of internal vertices does not change. Since Π 1 contains k anchors, after applying the twin process, each side of Π 2 contains at most m + k vertices (corners excluded) and 2n -k internal vertices. So, the size of

Π 2 is (m + k, 2n -k).
The internal vertices in Π 1 and in Π 2 are the same. And after each splitting anchor operation the paths P i 's span the same internal vertices than C i 's. So by Claim 1, P i 's also form a partition of the internal vertices of Π 2 .

It remains to show that during the twin splitting process no two edge-anchors can be twin or subdivided. So, consider any step t that performs a twin splitting operation at some anchor b of Π t -1 1 . We need to check that indeed the operation does not create twin edge-anchors or subdivide an edge-anchor. If the twin of b does not belong to an edge-anchor, then it is fine since the twin splitting of b will not result in a subdivision or a twin edge of an edge-anchor. So we only need to check the situation where the twin of b belongs to an edge-anchor. This occurs only if b and its twin a were both anchors in Π t-1 1 for some t < t . W.l.o.g. assume that t is the step where a is split into a -a , and let s be the neighbor of b in Π t 1 resulting of the twin splitting of a in Π t-1 1 .

We remark that two edge-anchors cannot be incident, and that on a side with path x-b-y that has been subdivided into x -b -s -y, then neither x -b nor b -s can be subdivided anymore (because the twin of b cannot be an anchor anymore). The path x -b -s can only be replaced by an edge-anchor x -b -b -s if b is an anchor.

We use the notations as above and as in Fig. 5. Let u -a -v be the vertices on the twin side of b in Π t-1 1 , where a is an anchor and twin of b. W.l.o.g. assume u ≺ a ≺ v, and let x, y be the twins of u, v respectively. In Π t 1 , we have u ≺ a ≺ a ≺ v, and there is some vertex s adjacent to b coming from an edge subdivision incident to b (of either b -y or x -b). Let us apply step t on b which produces the edge-anchor b -b in Π t 1 . To simplify notations we will reuse names u, v, x, y to denote the neighbors of a , a , b, s in Π t 1 as well as in Π t -1 1 .

Assume that b ≺ a and b ≺ y, i.e., the case 2 holds for a. In that case, we have the path

x -s -b -y in Π t -1 1
, and also u ≺ a ≺ a . By exchanging the roles of s -b -y with a -a -v the case 1 holds for b. This is because the condition a ≺ b rewrites in b ≺ a which is true (a plays the role of a). Therefore, the twin splitting of b subdivides the edge a -v which cannot be an edge-anchor by the previous remark. The edge-anchor b -b is twin with a non edge-anchor as well. Thus, step t does not twin or subdivide any edge-anchor in this case.

Assume that a ≺ b or y ≺ b, i.e., the case 1 holds for a. In that case, we have the path

x -b -s -y in Π t -1 1
, and also that u ≺ a ≺ a . Thus, in both cases, the twin splitting of b subdivides the edge u -a which cannot be an edge-anchor by the previous remark. The edge-anchor b -b is twin with a non edge-anchor as well. And therefore, in all the cases, step t does not twin or subdivide any edge-anchor.

The last property of Claim 3 (no edge-anchors is subdivided) ensures that the k edgeanchors in Π 1 are still existing in Π 2 after applying the twin process on Π 1 .

In order to transform Π 2 into an outerplanar embedding, we will apply at each of its edgeanchor a swapping operation. So, each step t of this process consists in applying one such operation which results into a new polygonal embedding Π t 2 obtained from Π t-1 2 , starting with Π 0 2 = Π 2 . The swapping operation is defined as follows. (We reuse the same notations as in Fig. 4 describing the splitting anchor operation which has been used to make Π 2 .) Let a i -a i be an edge-anchor of Π t-1 2 , P i be the path connecting a i to a i , and F i be the edges not in P i and embedded inside the cycle P i ∪ {a i -a i }. Finally, let b , b be the twin vertices of respectively a i , a i .

To obtain Π t 2 , we first delete the edge a i -a i and all edges of F i . By this way, the side containing a i -a i is extended by P i . Then, on the outerface of Π t-1 2 , we connect b to b by a copy P i of P i and including a copy F i of all the edges of F i such a way that the direction of P i is preserved. It means that if u is traversed before v when going from a i to a i on P i , then u is traversed before v when going from b in b in P i where u , v are the copies in P i of u, v respectively. Furthermore, the twin of any vertex u in P i is its copy of u in P i . By this way, the side containing b -b is extended by P i . See Fig. 6 for an illustration.

inside Π t-1 2 a i a i u v P i inside Π t 2 b b outerface of Π t-1 2 outerface of Π t 2 u v P i F i b b inside Π t-1 2 a i a i P i u v F i inside Π t 2
Figure 6: A swapping operation for the edge-anchor a i -a i in Π t-1 2 (on the left) leading to a new embedding Π t 2 (on the right). Vertex u is traversed before v when going from a to a on P i , and so for their twin vertices u , v when going from b to b on P i . The vertices of P i that were internal in Π t-1 2 belongs to the outerface of Π t 2 .

The last property of Claim 3 (no edge-anchors are twin) ensures that b -b was not an edge-anchor in Π t-1 2 . So the swapping operation on the edge-anchor a i -a i leaves unchanged all the other edge-anchors.

We denote by Π 3 the final embedding obtained from Π 2 by applying successively a swapping operation on all its edge-anchor. Let us denote by A ∼ = B if A, B are isomorphic graphs or homemorphic embeddings.

Claim 4 The embedding Π 3 is a polygonal embedding p-major of Π 2 with sewing(Π 3 ) ∼ = sewing(Π 2 ) and size (m + 2n, 0).

Proof. Up to a permutation of the indices, we can assume that the swapping operation for a i -a i is performed at step i from Π i-1 2 . It is easy to see that each operation does not create edge crossings, that the outerface remains a cycle, and that each vertex (edge) of each side has a well-defined twin vertex (edge) in its twin side. Therefore, Π i 2 is a polygonal embedding, and also Π 3 , that is Π k 2 , by transitivity. To prove it has the expected size, let n i = |V (P i ) \ {a i , a i } | be the number of internal vertices of P i in Π i-1 2 . Note that this number is the same as the number of vertices of P i in Π 0 2 = Π 2 because P i is not altered by any of the previous step j < i. Therefore, from Claim 3, k j=1 n j = 2n -k. The swapping operation for a i -a i increases by n i the number of vertices of its side and its twin side, and, at the same time, decreases by n i the number of internal vertices of Π i-1 2 . By Claim 3, Π 2 has size (m + k, 2n -k). So after step i, the size of

Π i 2 is (m + k + i j=1 n j , n - 2k -i j=1 n j ). It follows that Π 3 , which is Π k 2 , has size (m + 2n, 0) since k j=1 n j = 2n -k.
It is clear that signature and border, i.e., the sequence of corners, are not altered by any swapping operation. So, it remains to show that sewing(Π 3 ) ∼ = sewing(Π 2 ) (and by this way it will also show that Π 3 is a p-major of Π 2 ). By transitivity, it suffices to show that sewing(Π i

2 ) ∼ = sewing(Π i-1 2 ). Let R i be the plane graph composed of P i ∪ {a i -a i } ∪ F i in Π i-1 2 such that the boundary of its outerface is the cycle P i ∪ {a i -a i }. Similarly, let R i be the plane graph composed of

P i ∪ {b -b } ∪ F i in Π i 2 with cycle outerface P i ∪ {b -b }. Note that V (R i ) = V (P i ) and V (R i ) = V (P i ).
For convenience, denote by G j = sewing(Π j 2 ). We want to show that

G i-1 ∼ = G i .
Observe that in G i-1 and in G i , whenever the twin side are merged, a i = b , a i = b and thus the edge a i -a i = b -b exists in both graphs (whereas a i -a i exists only in Π i-1

2 ). Also, the way the vertices of P i are twin with the vertices P i (by preserving the order when going from a i to a i in P i ), ensures that P i = P i and

F i = F i in G i . It follows that R i ∼ = R i ,
and since the swapping operation alters only R i , V (G i-1 ) = V (G i ). And thus, if there is a difference between G i-1 and G i , it must be an edge.

From previous equalities, if an edge x -y lies in R i then the edge exists in R i and thus in G i . And conversely, if x -y lies in R i , the edge exists in R i and thus in G i-1 . If x -y lies outside R i (x, y may be both in V (R i ) with an embedding not inside R i ), then x -y lies also in Π i 2 since the swapping operation alters only R i (and its inside). So

x -y ∈ E(G i-1 ) implies x -y ∈ E(G i ) (even if x -y belongs a some side of Π i-1 2 ). Conversely, if x -y lies outside R i in Π i 2 , then x -y lies also in Π i-1 2 (and outside of R i ). It follows that edge sets of G i-1 and G i are the same, proving that G i-1 ∼ = G i .
By transitivity of the p-minor relation, every polygonal embedding Π of size (m, n) is a p-minor of Π 2 (by Claim 3) that is a p-minor of Π 3 of size (m + 2n, 0) (by Claim 4), which completes the proof of Lemma 2.

Proof of Lemma 3

The goal of this section is to show that:

Lemma 3 Every polygonal embedding Π of size (m, 0) is a p-minor of U σ(Π),m .
The proof is an adaptation of the proof of [RST94, Theorem (1.3)] to the case of outerplanar polygonal embeddings. The original proof uses grid major whereas we need to use U σ,m .

We first define the polygonal embedding U σ,m , which is defined for all signature σ and m ∈ N. To construct U σ,m , we start from a square half-grid of dimension |σ|m, that is with |σ|m rows and |σ|m columns. The vertices are column-row pairs of integers (i, j) where 0 j i < |σ|m, and (i, j) and (i , j ) are adjacent if and only if |i -i | + |j -j | = 1. Then, the vertices of the diagonal are connected by a cycle, (i, i) being connected to (i -1, i -1) modulo |σ|m. This cycle forms the cycle outerface of U σ,m . Finally, every m edges along this cycle is subdivised by one vertex, a corner. More precisely, for each i ∈ {0, . . . , |σ| -1}, we add a vertex c i between (im, im) and (im -1, im -1) modulo |σ|m. The border of U σ,m is the sequence (c 0 , . . . , c |σ|-1 ), and the ith side of U σ,m , associated with the ith symbol of σ, is the path between c i and c i+1 mod |σ| on the cycle outerface. See Fig. 7 and Fig. 2 for illustrations. Clearly, the number of internal vertices of U σ,m is n = 1+2+• • •+|σ|m-1. By merging the twin sides of U σ,m to obtain sewing(U σ,m ), we destroy at least half of the vertices of the cycle outerface of U σ,m . It follows that sewing(U σ,m ) has no more than n + |σ|(m + 1)/2 vertices. However, depending on the signature, this number can be lower as the set of corners may collaps even more when merging the sides. E.g., if the signature is canonical, sewing(U σ,m ) has only n + |σ|m/2 + 1 vertices since a canonical signature tell us that all the corners are twins, resulting into the single vertex of the canonical system of loops.

To summarize, we have:

Claim 5 The embedding U σ,m is a polygonal embedding of signature σ with 1 2 |σ|m • (|σ|m -1)) internal vertices, each of the |σ| sides having m non-corner vertices. In particular, sewing(U σ,m ) has no more than 1 2 (|σ| 2 m 2 + |σ|) vertices.

Let Π be a polygonal embedding of size (m, 0). Up to taking a p-major of Π, we can suppose that each side of Π is composed of exactly m non-corner vertices. If needed, one can subdivide an edge and its twin of any too small side, resulting in a p-major of Π. So, w.l.o.g., the cycle outerface of Π is composed of |σ(Π)| sides of m non-corner vertices each.

Since Π is outerplane, its cycle outerface contains all its vertices. The embeddings Π and U σ(Π),m have same signature and also m non-corner vertices per side. So both embeddings have isomorphic border and a cycle outerface of |σ(Π)|(m + 1) vertices.

We denote by (u 0 , . . . , u |σ(Π)|-1 ) the border of Π, and by v 0 , . . . , v |σ(Π)|m-1 its non-corners vertices ordered clockwise around its cycle outerface so that this cycle is

u 0 -v 0 -v 1 -• • • -v m-1 -u 1 -v m -• • • -v im-1 -u i -v im -• • •
Up to taking a p-major of Π, we will assume that all inner faces of Π are triangulated such that corners have degree two. For each non-corner v i , we define the two indices a(i), b(i) such that {v k : k ∈ [a(i), b(i)]} is the minimal subset of vertices containing v i and all its non-corner border and signature. By Property 1, Π is a p-minor of U σ(Π),m that completes the proof of Lemma 3.

Conclusion

We have shown that Euler genus-g graphs have an Euler genus-g minor-universal graph of O(g 2 (n + g) 2 ) vertices. It is not clear whether the dependencies in g and in n are tight. Even for the planar case, it is not obvious that the n 2 term is required, noting that the Ω(n 2 ) known lower bound exists only if the minor-universal has to be a grid.

It would be interesting to extend the result to other graph families, in particular bounded treewidth graphs and minor-free graphs. Another direction could be to extend the result to topological minor-universal graphs.

Figure 2 :

 2 Figure 2: The polygonal embedding U σ,m , with m = 3 and signature σ = a 1 a 2 ā1 ā2 .

  Fig. 5): • a subdivision of the edge b -y, if a ≺ b or y ≺ b; or • a subdivision of the edge x -b, if b ≺ a and b ≺ y.

Figure 5 :

 5 Figure5: The twin splitting operation for an anchor a: it combines a splitting anchor operation on a (in brown) and an edge subdivision (in blue) incident to its twin b. Twin relations are with dotted lines. In this representation, the focus on the twin sides creates the illusion that the inside part of the embedding (in gray) is a disconnected region, which is not true. However, in this representation choice, the inside part, that is locally planar, may appear twisted if fully represented in the case of a non-orientable surface.

  If y ≺ b, then we have y ≺ s ≺ b ≺ x. By exchanging the roles of s -b -x with u -a -a the case 1 holds for b. This is because the condition u ≺ a rewrites in y ≺ b which is true. If b ≺ y, then we have x ≺ b ≺ s. By exchanging x -b -s with u -a -a the case 1 holds again for b. This is because the condition a ≺ a rewrites in b ≺ y which is true (s plays the role of y).

Figure 7 :

 7 Figure7: The polygonal embedding U σ,m for m = 2 and σ = a 1 a 2 a 1 a 2 , a non-canonical signature for a non-orientable surface of genus 2. Its size is (2, 28). The cycle outerface isc 0 -(0, 0) -(1, 1) -c 1 -(2, 2) -(3, 3) -c 2 -(4, 4) -(5, 5) -c 3 -(6, 6) -(7, 7) -c 0 .

Recall that H is a minor of G if H can be obtained from G by vertex removals, edge removals, and edge contractions.

Using the Euler's Formula.

For these results, rooted binary trees and topological minors where considered. However, this applies to our setting since, up to a factor of two, the rooted and unrooted versions are equivalent, and minor and topological minor containments are equivalent in binary trees.

We only consider simple graphs and cellular embeddings, i.e., where each face is homeomorphic to a disc.

We refer to[START_REF] De Verdière | Algorithms for embedded graphs[END_REF][Theorem 8.1] for a reformulation of the original statement of[START_REF] Lazarus | Computing a canonical polygonal schema of an orientable triangulated surface[END_REF].

This is only 4 points in the case of orientable surface.

From Euler's formula in simple graphs.

|σ(Π)| denotes the number of symbols in σ(Π).

In fact ΠG is a minor of sewing(ΠG) \ {r} as the unique vertex r of system of loops is not part of G and can be removed.

Every planar triangulation without separating triangle has a Hamiltonian cycle.

Acknowledgement. We are grateful to Arnaud de Mesmay for prelimary discussions and suggestions for this problem.

* gavoilleThis work is partially funded by the French ANR projects ANR-17-CE40-0015 (DISTAN-CIA), ANR-18-CE47-0010 (QUDATA) and ANR-22-CE48-0001 (TEMPOGRAL).

neighbors. Note that a(i) = i or b(i) = i is possible. For technical reasons, we increase b(0) to b(0) = |σ(Π)|m.

In order to show that Π is a minor of U σ(Π),m , we construct a witness of Π in U σ(Π),m with the property that it preserves its outerface and its border. A witness for a minor Π in U σ(Π),m is a collection {W (u)} u∈V (Π) of nonempty pairwise disjoint subsets of vertices of U σ(Π),m such that each set W (u) induces a connected subgraph of U σ(Π),m , and for each edge

The witness of Π in U σ(Π),m is defined by:

Clearly, for corner vertices, W (u i ) = {c i } fulfills all desired properties: non-emptyness, disjointness and connectedness. For non-corner vertices,

Let us check that W (v i ) and W (v j ) are pairwise disjoint (witness for corner and non-corner are clearly disjoint). Assume j < i. A non-empty intersection is only possible between C(v i ) and R(v j ). And, this can only occur at (i, j). Note that the case j = 0 provides disjoint witnesses, since

and thus v i has a neighbor v a(i) with a(i) < j. If (i, j) ∈ R(v j ), then i ∈ [j, b(j)), and thus v j has a neighbor v b(j) with i < b(j). It follows that a(i) < j < i < a(j) and that v a(i) -v i and v j -v b(j) is a pair of crossing edges, which is impossible in the outerplanar embedding Π.

It remains to check that, for each u -v ∈ E(Π), W (u) ∪ W (v) induces a connected graphs in U σ(Π),m . We can restrict our attention to non-corner neighbors, since c i ∈ W (u i ) have exactly two neighbors that are (im -1, im -1) ∈ W (v im-1 ) and (im, im) ∈ W (v im ).

So consider an edge v i -v j of Π with j < i. Let W (v i , v j ) be the set of vertices of the path connecting (i, i) ∈ W (v i ) to (j, j) ∈ W (v j ) inside the grid part of U σ(Π),m , and defined by W (v i , v j ) = {(i, i), (i, i + 1), . . . , (i, j -1) , (i, j), (i -1, j), . . . , (j, j)} .

) is possible. If j = 0, then we are done because b(0) = |σ(Π)|m and thus (i, 0) ∈ R(0). We are also done if |i -j| = 1 since (i, i) ∈ W (v i ) and (j, j) ∈ W (v j ) are neighbors in U σ(Π),m . Since the inner faces of Π are triangulated, the edge v i -v j (that lies inside Π since v i and v j are not consecutive in the cycle outerface) share a triangle v i -v j -v k such that k / ∈ [i, j], the case j = 0 and i = |σ(Π)|m -1 being excluded. If k > i, then a(j) > i and (i, j) ∈ C(v j ). If k < j, then b(i) > j and (i, j) ∈ R(v i ).

We have therefore proved that W (v i , v j ) ⊆ W (v i )∪W (v j ) and that W (v i )∪W (v j ) induces a connected component in U σ(Π),m . It follows that Π is a minor of U σ(Π),m with same outerface,


domain: nlin math

text:
Alexander N Pchelintsev
email: pchelintsev.an@yandex.ru

On the Poisson Stability to Study a Fourth-Order Dynamical System with Quadratic Nonlinearities

Keywords: attractor, Poisson stability, Lyapunov exponents MSC:

This article discusses the search procedure for the Poincaré recurrences to classify solutions on an attractor of a fourth-order nonlinear dynamical system using a previously developed high-precision numerical method. For the resulting limiting solution, the Lyapunov exponents are calculated using the modified Benettin's algorithm to study the stability of the found regime and confirm the type of attractor.

Introduction

As it is known [START_REF] Sprott | Elegant Chaos: Algebraically Simple Chaotic Flows[END_REF], the calculation of the Lyapunov exponents is used for the classification of attractors of dynamical systems. The combinations of signs of such values determine the attractor type: an equilibrium position, limit cycle, multidimensional torus or strange attractor. The dynamics of the corresponding solutions are called static, periodic, quasiperiodic, or chaotic. The computational procedure for determining the Lyapunov exponents is mainly based on the Benettin's algorithm [START_REF] Tancredi | A comparison between methods to compute Lyapunov exponents[END_REF]. However, in [START_REF] Lozi | A new reliable numerical method for computing chaotic solutions of dynamical systems: the Chen attractor case[END_REF][START_REF] Lozi | A new accurate numerical method of approximation of chaotic solutions of dynamical model equations with quadratic nonlinearities[END_REF] the points of limiting solutions were investigated for the Poisson stability, which made it possible to understand whether we have a quasiperiodic or chaotic regime. In [START_REF] Lozi | A new accurate numerical method of approximation of chaotic solutions of dynamical model equations with quadratic nonlinearities[END_REF], it was done for the Jafari-Sprott system [START_REF] Jafari | Recent new examples of hidden attractors[END_REF].

It is noteworthy [START_REF] Nemytskii | Qualitative Theory of Differential Equations[END_REF] that a point y of the phase space is called positively Poisson stable (P + ), if for any neighborhood U of the point y and for any T P > 0 there is a time value t ≥ T P for which the trajectory of the dynamical system enters into the neighborhood U. Similarly, if there is t ≤ -T P such that the trajectory enters into the neighborhood U, then the point y is negatively Poisson stable (P -). A point is called Poisson stable, if it is P + and P --stable.

The boundedness of the limiting solutions of dissipative systems implies [START_REF] Nemytskii | Qualitative Theory of Differential Equations[END_REF][START_REF] Afanas'ev | Method for constructing minimal sets of dynamical systems[END_REF][START_REF] Afanas'ev | Construction of the minimal sets of differential equations with polynomial right-hand side[END_REF][START_REF] Afanas'ev | On recurrent trajectories, minimal sets, and quasiperiodic motions of dynamical systems[END_REF][START_REF] Kuznetsov | Dynamical Chaos[END_REF] that any steady-state oscillation mode is described by Poisson-stable trajectories. This also applies to dynamical chaos. A trajectory different from the equilibrium position is said to be Poisson stable if it verifies the property of returning in an arbitrarily small εneighborhood of each of its points an infinite number of times. Such returns are called the Poincaré recurrences. In [10, p. 146], the author notes that "the study of the statistics of the Poincaré recurrences is a powerful tool for the analysis and classification of dynamic modes. Apparently, the potentialities of this approach have not yet been fully exhausted in modern nonlinear dynamics". For example, the returns follow one another regularly for quasiperiodic regimes. Then [10, p. 145] "the dynamic chaos is a situation when the Poincaré recurrences to the ε-neighborhood of the initial point do not show regularity, the time interval between two successive returns turns out to be different each time, and some statistical distribution of the times of return is arised" [START_REF] Afraimovich | Fractal dimension for Poincaré recurrences as an indicator of synchronized chaotic regimes[END_REF][START_REF] Afraimovich | Fractal and multifractal properties of exit times and Poincaré recurrences[END_REF]. An example of the Poincaré recurrences analysis based on the Kac's theorem [START_REF] Kac | Probability and Related Topics in Physical Sciences[END_REF][START_REF] Penné | Fractal and statistical characteristics of recurrence times[END_REF][START_REF] Afraimovich | Fractal Dimensions for Poincaré Recurrences[END_REF][START_REF] Anishchenko | Statistical characteristics of the Poincaré return times for a onedimensional nonhyperbolic map[END_REF][START_REF] Anishchenko | Poincaré recurrences time and local dimension of chaotic attractors[END_REF][START_REF] Anishchenko | Local and global approaches to the problem of Poincaré recurrences[END_REF], for a discrete dynamical system with a chaotic non-hyperbolic attractor is given in [START_REF] Anishchenko | Statistical characteristics of the Poincaré return times for a onedimensional nonhyperbolic map[END_REF].

As it is known, for the most of dissipative systems, the formulas for a general solution of the system in a class of any known functions with well-studied properties have not yet been found. Therefore, numerical methods are used. The use of classical numerical methods (such as Euler method, Runge-Kutta 4th order method, Adams methods, etc.) for constructing approximate solutions in attractors of dynamical systems leads to significant errors over large sections of time due to the instability of the studied chaotic regimes. In general, the problem of numerical modeling with control of the accuracy of obtained solution and the choice of a platform for computer implementation is relevant today [START_REF] Guedes | Some remarks on the performance of Matlab, Python and Octave in simulating dynamical systems. Anais do 14[END_REF], since small errors introduced at each integration step cause exponential divergence of close trajectories. More recently, the numerical FGBFI method (Firmly Grounded Backward-Forward Integration) has been proposed in [START_REF] Lozi | A new reliable numerical method for computing chaotic solutions of dynamical systems: the Chen attractor case[END_REF][START_REF] Lozi | A new accurate numerical method of approximation of chaotic solutions of dynamical model equations with quadratic nonlinearities[END_REF][START_REF] Pchelintsev | An accurate numerical method and algorithm for constructing solutions of chaotic systems[END_REF], based on the power series method for dynamical systems with quadratic nonlinearities. The main advantages of this method are as follows:

1.

The recurrence relations are obtained for calculating the coefficients of the expansion of solutions in a power series for any dynamical system with quadratic nonlinearities in a general form; 2.

The convergence of the power series is studied. A simple formula for calculations is derived (in comparison with that [START_REF] Babadzhanjanz | Taylor series method for dynamic systems with control: convergence and error estimates[END_REF] obtained in the known literature) for calculating the length of the integration step in a general form; 3.

The criteria for checking the accuracy of the approximate solution are obtained. The control of the accuracy and configuration of the approximate solution of a dynamical system uses forward and backward time, which makes the numerical method reliable (degrees of piecewise polynomials, the value of the maximum integration step, etc.); 4.

The FGBFI method allows to construct high-precision approximations to non-extendable solutions of a system of autonomous differential equations with a quadratic right-hand side. Like for instance the system of the form

ẋ = 1 + x 2 .
In this case, the numerical solution computed with FGBFI will never cross the asymptote and will approach it arbitrarily close.

Thus, it is a high-precision method for constructing the trajectory arc of the system for any time interval. This makes it possible to track the Poincaré recurrences to any neighborhood of the trajectory points, since the resulting computational error can be smaller than the radius of the monitored neighborhoods.

This article considers a fourth-order dissipative system [START_REF] Dong | Ultimate boundary estimation and topological horseshoe analysis on a parallel 4D hyperchaotic system with any number of attractors and its multi-scroll[END_REF], in which, in the opinion of its authors, there is a hyperchaos. The goals of the research are: (1) to find the Poincaré recurrences in the attractor of this system for the approximate solutions obtained by the FGBFI method and to analyze the statistics of return times; (2) to apply the modified Benettin's algorithm [START_REF] Pchelintsev | An accurate numerical method and algorithm for constructing solutions of chaotic systems[END_REF] to refine the computation of the Lyapunov exponents.

Bohr's almost periodic functions

It is noteworthy [23, pp. 368, 418, 419] that the function (trajectory of a dynamical system) F(t) ∈ R n is called Bohr's almost periodic in the sense that if for any ε > 0 there exists a relatively dense set of almost periods τ F = τ F (ε) of the function F(t) with ε-accuracy, i.e. there is the positive number L = L(ε) such that any interval [α; α + L] contains at least one number τ F for which the inequality |F(t + τ F ) -F(t)| < ε holds at t ∈ R. Note that for an almost periodic function (and, in general, for a Poisson stable trajectory), different from a periodic function, with decreasing number ε the number L(ε) must increase indefinitely, otherwise the function would be periodic. Usually (see, for example, [START_REF] Kuznetsov | Dynamical Chaos[END_REF]) in nonlinear dynamics, almost periodic functions are called quasiperiodic.

Note an important theorem in the theory of dynamical systems [10, p. 147]: if a non-periodic trajectory is Poisson and Lyapunov stable, then it is almost periodic.

Since we will be investigating the Poincaré recurrences, it is noteworthy that for an almost periodic trajectory (as follows from the above definition), such returns must form a relatively dense set. An example of a sequence that is not a relatively dense set is: 0, 1 2 , 2 2 , 3 2 , . . . , since sup

q |(q + 1) 2 -q 2 | = ∞,
i.e. in the given example, the distance between neighboring elements grows with the growth of q. A similar situation was observed in a numerical experiment for the Poincaré recurrences in the case of a chaotic solution of the Chen system [START_REF] Lozi | A new reliable numerical method for computing chaotic solutions of dynamical systems: the Chen attractor case[END_REF].

Dissipative system of the fourth order

Let us consider a fourth-order system [22]

       ẋ1 = ax 2 -ax 1 -ex 4 , ẋ2 = bx 1 -x 2 -x 1 x 3 -f x 4 , ẋ3 = x 1 x 2 -cx 3 , ẋ4 = kx 2 x 3 -dx 4 , (1)
where a, b, c, d, e, f and k are positive parameters. The divergence of the vector field G defined by the right-hand side of the system (1) is equal to div G = -(a + 1 + c + d) < 0.

Then the system (1) is dissipative. Hence, there is the ball B a ⊂ R 4 , into which each trajectory of the system (1) is immersed forever, after a while. Therefore, there is a limit set (attractor) such that all trajectories of the dynamical system are attracted when t → ∞ [START_REF] Nemytskii | Qualitative Theory of Differential Equations[END_REF].

Let us investigate the dynamics of the system (1) for the values of the parameters a = 7, b = 50, c = 3, d = 10, e = 5, f = 5 and k = 1.5. In this case, it is possible to determine the position and radius of the ball B a from the data obtained in [START_REF] Dong | Ultimate boundary estimation and topological horseshoe analysis on a parallel 4D hyperchaotic system with any number of attractors and its multi-scroll[END_REF]. Also, in this article the Lyapunov exponents were determined: LE 1 = 2.1040, LE 2 = 0.3563, LE 3 = 0 and LE 4 = -23.4508, indicating a hyperchaos in the system (1). Therefore, its solutions are highly unstable, which requires the use of special numerical procedures for large time intervals.

Let us consider the application of the FGBFI method for constructing approximations to the solutions of the system (1) using high-precision calculations.

The FGBFI method

Following [START_REF] Lozi | A new accurate numerical method of approximation of chaotic solutions of dynamical model equations with quadratic nonlinearities[END_REF], we will rewrite the system (1) as

Ẋ = AX + Φ(X), (2)
where

X(t) = [x 1 (t) x 2 (t) x 3 (t) x 4 (t)] T , Φ(X) = [φ 1 (X) φ 2 (X) φ 3 (X) φ 4 (X)] T , φ p (X) = ⟨Q p X, X⟩, p = 1; 4, A =     -a a 0 -e b -1 0 -f 0 0 -c 0 0 0 0 -d     , Q 1 = 0, Q 2 =     0 0 -1 0 0 0 0 0 0 0 0 0 0 0 0 0     , Q 3 =     0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0     , Q 4 =     0 0 0 0 0 0 k 0 0 0 0 0 0 0 0 0     .
Let us expand the solution

X(t) = ∞ ∑ i=0 Υ i t i , (3)
where X(0) = Υ 0 is an initial condition for the system (2), the vector Υ 0 is given,

Υ i = [υ i,1 υ i,2 υ i,3 υ i,4 ] T , υ i,p ∈ R, p = 1; 4. Let Ψ i = [ψ i,1 ψ i,2 ψ i,3 ψ i,4 ] T and ψ i,p = i ∑ j=0 ⟨Q p Υ j , Υ i-j ⟩, p = 1; 4.
The recurrence relations for calculating the coefficients of the power series have the form [4]

Υ j = AΥ j-1 + Ψ j-1 j (4)
for j ∈ N.

Since the criteria for checking the accuracy of the obtained numerical solution require repeated calculations in forward and backward time, we need to have a guaranteed estimate of the convergence interval of the power series (3) for a given vector Υ 0 . In [START_REF] Lozi | A new accurate numerical method of approximation of chaotic solutions of dynamical model equations with quadratic nonlinearities[END_REF], a theorem is proved on the estimate of this interval for systems with a quadratic right-hand side.

The norms are calculated:

∥A∥ = ∥A∥ 1 = max{a + b, a + 1, c, e + f + d} = 57, ∥Q 1 ∥ = 0, ∥Q 2 ∥ = ∥Q 3 ∥ = 1, ∥Q 4 ∥ = k, µ = 4 • max p=1;4 ∥Q p ∥ = 4 • 1.5 = 6.
Next, two auxiliary numbers are calculated as functions of Υ 0 :

h 1 (Υ 0 ) = ∥Υ 0 ∥ = 4 ∑ p=1 |υ 0,p |, h 2 (Υ 0 ) = µh 2 1 (Υ 0 ) + (∥A∥ + 2µ)h 1 (Υ 0 ), if h 1 > 1, ∥A∥ + µ otherwise, or h 2 (Υ 0 ) = 6 • ∥Υ 0 ∥ 2 + 69 • ∥Υ 0 ∥, if ∥Υ 0 ∥ > 1, 63 otherwise.
Then the series (3) converges for t ∈ [-τ, τ], where

τ = τ(Υ 0 ) = 1 h 2 (Υ 0 ) + δ , (5)
δ is any positive number (can be taken very small).

It can be seen from the formula (5) that the length of the convergence interval depends on the choice of the initial condition for the system (2). This fact is confirmed by computational experiments, for example, for the Lorenz system [START_REF] Pchelintsev | Numerical and physical modeling of the dynamics of the Lorenz system[END_REF].

It is worth noting that a similar approach can be applied to systems with cubic nonlinearities, but when the derivative on the left side of the system is of the second (or higher) order. For example, in [START_REF] Pchelintsev | Solution of the Duffing equation by the power series method[END_REF] this was done for the Duffing equation.

Usually, in calculations, many researchers work with subnormal real numbers of single or double precision, presented in the IEEE 754 format [START_REF] Overton | Numerical Computing with IEEE Floating Point Arithmetic[END_REF]. The main drawback here is the fixed accuracy of the representation of real numbers, which may not allow us numerically constructing approximations to unstable solutions of systems of differential equations on the large time intervals. Therefore, for high-precision calculations, the GNU MPFR library [START_REF]GNU MPFR Library for Multiple-Precision Floating Point Computations with Correct Rounding[END_REF][START_REF] Fousse | MPFR: A multiple-precision binary floating-point library with correct rounding[END_REF] is used.

Next, consider the algorithm [START_REF] Pchelintsev | An accurate numerical method and algorithm for constructing solutions of chaotic systems[END_REF] for the numerical construction of an approximate solution of a dynamical system in forward and backward time:

1.
Set the number b m of bits under the mantissa of a real number and the precision ε pw for estimate of the common term of power series. Note that b m defines the machine epsilon ε m . Let us choose b m so that the precision of representation of the real number is with a margin, i.e.

ε m ≪ ε pw ; 2. t := 0; 3.
Set X(0) ∈ B a for the system (2) and value the number way that determines the direction in time: way = 1 is gone forward in time, way = -1 is gone backward in time; 4.

Set T: the length of the time interval on which the numerical integration will be performed; 5.

ended := false; 6.

Calculate the integration step ∆t := τ(X(0)) by the formula (5); 7.

If ∆t > Tt, then ∆t := Tt, t := T Else t := t + ∆t; 8.

∆t := way • ∆t; 9.

Let Υ 0 := X(0); As can be seen, the algorithm is universal in the direction of numerical integration in time. Next, we consider the criteria for verifying the accuracy of the resulting solution. For this we assume

Ω = [t 0 , t 1 ] ∪ [t 1 , t 2 ] ∪ . . . ∪ [t N-1 , t N ], t 0 = 0, t N = way • T, ∆t {way} l = t l -t l-1 = τ X l-1 ∆t {way} l-1 , X 0 ≡ X(0), l = 1; N {way} ,
where l is an index of the interval [t l-1 , t l ], on which the series (3) is converged, N {way} is a number of such intervals for the direction way on time. Since on each interval [t l-1 , t l ] during calculations we truncate the series (3) to some polynomial X l (t) approximating the solution of the system (2) on it (i.e. locally), then we denote by n {way} l the degree of the polynomial X l (t). Then we introduce the notation:

n {way} max = max l n {way} l , l {way} max = indmax l n {way} l , ∆t {way} max = way • max l ∆t {way} l , d {way} max = indmax l ∆t {way} l . (7)
In [START_REF] Lozi | A new accurate numerical method of approximation of chaotic solutions of dynamical model equations with quadratic nonlinearities[END_REF], the following criteria for checking the accuracy of the obtained solution were proposed:

1.

The accuracy ε a of approximation at way = 1 is a frequently used criterion in applications of numerical methods for solving differential equations. When the inequality ( 6) is true, it is necessary to increase the degrees of all polynomials n {1} l , obtaining the next approximation, and compare the distance δ a between the obtained approximate solutions on the interval [0, T] with the value ε a . If δ a > ε a , then we increase the powers of n l , otherwise we have to use the obtained solution; 2.

The radius ε R of the neighborhood of the initial point, to which the approximate solution should return in backward time is another criterion. In other words, we need to select the accuracy of ε pw so that the following inequality holds

X {-1} N ∆t {-1} N -X(0) < ε R ,
where N = N {-1} . Note here that we do not know the exact solution to the system (2), but we do know its initial point. Then we can set how many digits are in each component of the vector

X {-1} N ∆t {-1} N
after the point must match the digits of the corresponding components X(0). The problem is that a large amount of computation is required (due to repetitions of the forward and backward traverse in time). The solutions of the system (2), as a rule, are highly unstable in the backward time: they immediately leave from the attractor, because in our calculations, we are near it, no matter how accurate ε pw is taken. Therefore, in the above algorithm, we control the finding of the trajectory within the boundaries of the set B a ; 3.

The comparison of configurations of approximate solutions in forward and backward time is necessary to determine the numbers [START_REF] Afanas'ev | Method for constructing minimal sets of dynamical systems[END_REF], describing approximate solutions. Next, check

N ≈ N, n {1} max ≈ n {-1} max , t l {1} max + t l {-1} max ≈ T, d {1} max + d {-1} max ≈ N, t d {1} max + t d {-1} max ≈ T.
Unlike criteria 1 and 2, we control here the arguments of the approximate solutions.

Analysis of the Poincaré recurrences for a fourth-order system

In [22, p. 3221], the authors indicated that there is a hyperchaos for the considered parameters of the system (1). For example, the authors found it for the initial condition We get the next point at the end time: X(15) = [6.2355509634533960831 2.0140572482317481452 35.4929323328531102196 -43.5507482101916799734] T . [START_REF] Afanas'ev | On recurrent trajectories, minimal sets, and quasiperiodic motions of dynamical systems[END_REF] Note that in the numerical experiment by to the first criterion for checking the accuracy, decreasing the value of ε pw , the number of indicated correct signs in coordinates is preserved.

Next, we take the point (9) as the initial point and track the Poincaré recurrences for it. In order to reduce the amount of calculations without moving backward in time (so as not to take too small a value of ε pw ), we look at the times when there is a return to the neighborhood of the initial point for different values ε pw . This is necessary to control the accuracy of the obtained points, choosing such the value ε pw that these moments of time do not differ significantly.

To determine the Poincaré recurrences, we choose the small time step ∆t P > 0. Let us divide the time interval [0, T P ], where T P is the end time at which the search for the Poincaré recurrences is terminated, into intervals of length ∆t P (let N P be the number of such intervals). Let us denote X k the coordinates of the k-th point of the considered trajectory arc, corresponding to the time moment t k = k∆t P , k = 1; N P . We introduce the distance

d k = |X k -X(0)|.
To track the Poincaré recurrences, we fix such values k = k * , moving in the trajectory arc, when there are the local rapprochements with the point X(0), i.e.

d k * -1 > d k * & d k * < d k * +1 & d k * < 1.
For the point (9) for T P = 10 and ∆t P = 10 -4 , the returns are shown in Tabl. 1. As can be seen from this table, the rapprochement with the initial point occurs at approximately the same time (i.e., we have regularity), which suggests the idea of a limit cycle in the system [START_REF] Sprott | Elegant Chaos: Algebraically Simple Chaotic Flows[END_REF]. Decrease in the value of ∆t P does not cause significant changes in Tabl. 1.

In this case to show the absence of chaos in the system (1), we take the time value T = 40 (instead of T = 15, which was used to obtain the coordinates of the point ( 9)) to reach the attractor. We get the following point at the end time: X(40) = [1.6321991613781496393 8.7300523565474285155 39.6961687172415982460 54.8461996449311966025] T . [START_REF] Kuznetsov | Dynamical Chaos[END_REF] Similarly, for the point [START_REF] Kuznetsov | Dynamical Chaos[END_REF], the Poincaré recurrences shown in Tabl. 2 for T P = 10 and ∆t P = 10 -4 . As we can be seen from this table, returns are getting closer. Therefore, we will decrease the value ∆t P , setting it equal to ∆t P = 10 -7 . We get the returns shown in Tabl. 3.

From this iterative procedure, we see that the limiting trajectory in this case is a cycle with a period T c ≈ 0.3655. Its projections are shown in Figs. 1 and2.

To select the value ∆t P , this step must be reduced until the rapprochement distance d k * changes significantly. This approach was also applied for other initial conditions. As a result, we get the same cycle. Next, we investigate the stability of the resulting cycle by calculating the Lyapunov exponents.

Calculation of the Lyapunov exponents

In [START_REF] Pchelintsev | An accurate numerical method and algorithm for constructing solutions of chaotic systems[END_REF], a modification of the algorithm for calculating the Lyapunov exponents for a model describing the growth of a cancerous tumor was considered. Next, we will consider a generalization of this algorithm for dynamical systems with a quadratic right-hand side of the n-th order. Note that an other algorithm for calculating the Lyapunov exponents based on calculating the fundamental matrix of the dynamical system and its QR decomposition is proposed in [START_REF] Kuznetsov | Finite-time Lyapunov dimension and hidden attractor of the Rabinovich system[END_REF]. The procedure described below can also be applied to it. First, we define the form of the linearized system of equations. Let us denote the perturbations by

x n+1 (t), x n+2 (t), . . . , x 2n (t).

We define the form of the following expression using the rule for differentiating the sums and products:

∂⟨Q p X, X⟩ ∂x j • x n+j = ∂ ∂x j n ∑ i=1 n ∑ k=1 q (p) i,k x i x k • x n+j = = n ∑ k=1 q (p) j,k x k + n ∑ i=1 q (p) i,j x i • x n+j = n ∑ i=1 q (p) i,j + q (p) j,i x i x n+j ,
where p = 1; n, j = 1; n. Let us find the product

∂⟨Q p X, X⟩ ∂x 1 ∂⟨Q p X, X⟩ ∂x 2 . . . ∂⟨Q p X, X⟩ ∂x n •      x n+1 x n+2 . . . x 2n      = = n ∑ j=1 n ∑ i=1 q (p) i,j + q (p) j,i x i x n+j = n ∑ i=1 2n ∑ j=n+1 q (p) i,j-n + q (p) j-n,i x i x j . (11)
Next, we extend the system (1) by adding to it a linearized system, which also has a quadratic right-hand side. To do this, we introduce the extended matrix Ã:

à = ãi,j 2n×2n , ãi,j =    a i,j , if i = 1; n & j = 1; n, a i-n,j-n , if i = n + 1; 2n & j = n + 1; 2n, 0 otherwise.
The vector function X(t) is equal to

X(t) = [x 1 (t) x 2 (t) . . . x n (t) x n+1 (t) . . . x 2n (t)] T .
Based on the formula [START_REF] Afraimovich | Fractal dimension for Poincaré recurrences as an indicator of synchronized chaotic regimes[END_REF], we have a quadratic form on the right-hand side of the linearized system for each equation. Then, to reduce the system to general form, we introduce the following extended matrix:

Ql = q(l) i,j 2n×2n , for l = 1; n : q(l) i,j = q (l) i,j , if i = 1; n & j = 1; n, 0 otherwise, for l = n + 1; 2n : q(l) i,j = q (l-n) i,j-n + q (l-n) j-n,i , if i = 1; n & j = n + 1; 2n, 0 otherwise.
We get an extended system similar to the system (2):

Ẋ = Ã X + Φ( X), (12)
where Φ = φ1 ( X) . . . φ2n ( X) T , φl ( X) = ⟨ Ql X, X⟩, l = 1; 2n.

To obtain the domain of convergence of the power series that is a solution to the system [START_REF] Afraimovich | Fractal and multifractal properties of exit times and Poincaré recurrences[END_REF], the norms of the matrices are calculated ∥ Ã∥, ∥ Q1 ∥, ∥ Q2 ∥, . . . ∥ Q2n ∥ and the number μ = 2n max p=1;2n

∥ Qp ∥. Let X(t) = ∞ ∑ i=0 Υi t i , (13)
where Υ0 = X(0) is the initial condition for the system ( 12 By analogy with the system (2), two auxiliary numbers are calculated as functions of Υ0 :

h1 ( Υ0 ) = ∥ Υ0 ∥, h2 ( Υ0 ) = μ h2 1 ( Υ0 ) + (∥ Ã∥ + 2 μ) h1 ( Υ0 ), if h1 > 1, ∥ Ã∥ + μ otherwise.
Then the power series [START_REF] Kac | Probability and Related Topics in Physical Sciences[END_REF] Divide a given time interval [0, T] (usually the value T is large) into short intervals by the length

τ M = T M ,
where M is the number of such intervals.

2.

Let X * is a point of the researched solution; e.g., (10); 3.

k := 0; 4. Y (k) := X * ; 5. Let Z (k) (p)
is the column vectors of initial perturbations from n-components (the real numbers), p = 1; n (p is a number of the perturbation vector);

6.

Assign to each component of the vector Z

(p) random number in the range [0, 1]; 7.

LE p := 0, p = 1; n (the initial values of the sums at calculating the Lyapunov exponents);

8.

Orthogonalize and normalize to unity the system of vectors Z

(1) , Z

(2) , ..., Z

(n) ; 9.

k := k + 1; 10. For p from 1 to n with step 1:

Beginning of cycle Let

X(0) := Y (k-1) Z (k-1) (p)
.

Find X(τ M ) for the system (12) on algorithm, described in Sec. 4. Assign to vectors Y (k) and Z

(p) the matching vector blocks X(τ M ), i.e.

Y (k) Z (k) (p)
:= X(τ M ).

Note that for all values of p the vector Y (k) is the same. End of cycle 11. For p from 1 to n with step 1:

Beginning of cycle For the point [START_REF] Kuznetsov | Dynamical Chaos[END_REF], the Lyapunov exponents were found using the described algorithm for T = 100, M = 20000 and τ M = 0.005. Next, we give them the approximate values:

S := Z (k) (p) - p-1 ∑ i=1 Z (k) (p) , Z (k)
(i) • Z (k) (
LE 1 ≈ 0, LE 2 ≈ -0.498, LE 3 ≈ -0.499, LE 4 ≈ -20.002.
As can be seen, their signature corresponds to a stable limit cycle in the system (1), which confirms the presence of the cycle, found through the study of the Poincaré recurrences. Note that in [START_REF] Pchelintsev | A numerical-analytical method for constructing periodic solutions of the Lorenz system[END_REF] an alternative numerical-analytical scheme for constructing periodic solutions of systems with a quadratic right-hand side was proposed using an example of the Lorenz system. Now this scheme is actively developing in [START_REF] Luo | Approximate solutions of periodic motions in nonlinear systems via a generalized harmonic balance[END_REF][START_REF] Luo | Toward Analytical Chaos in Nonlinear Systems[END_REF][START_REF] Luo | Analytical solutions of period-1 to period-2 motions in a periodically diffused brusselator[END_REF] to find periodic solutions of nonlinear systems of differential equations.

The graphs of changes in the Lyapunov exponents over time are shown in Figs. 3456.

The source codes of all computing programs are located at [START_REF] Pchelintsev | The reliable calculations for the 4-th order system[END_REF].

Conclusion

In this article, the author investigates the limiting regimes corresponding to different initial points of the phase space. As a result, we get the same cycle shown in Fig. 1 and2. Consequently, there is no hyperchaos in system (1). The analysis of the Poincaré recurrences, described in the goals of research, showed that the rapprochement with the initial point occurs approximately at the same time (i.e., we have a regularity). This indicates a limit cycle in the system [START_REF] Sprott | Elegant Chaos: Algebraically Simple Chaotic Flows[END_REF].

Stabilization of the behavior of the values of the Lyapunov exponents in Figs. 3456indicates the correctness of finding their limiting values over a relatively short time interval [0, 100].

Thus, the limiting trajectory investigated using the analysis of the Poincaré recurrences and the high-precision calculation of the Lyapunov exponents (i.e., different methods) is the limit cycle which contradicts the results in [START_REF] Dong | Ultimate boundary estimation and topological horseshoe analysis on a parallel 4D hyperchaotic system with any number of attractors and its multi-scroll[END_REF] for the values of the parameters a = 7, b = 50, c = 3, d = 10, e = 5, f = 5 and k = 1.5 of the system [START_REF] Sprott | Elegant Chaos: Algebraically Simple Chaotic Flows[END_REF]. The obtained data is available at [START_REF] Pchelintsev | The reliable calculations for the 4-th order system[END_REF]. Most likely, the authors of the article [START_REF] Dong | Ultimate boundary estimation and topological horseshoe analysis on a parallel 4D hyperchaotic system with any number of attractors and its multi-scroll[END_REF] made a typo indicating such values of the parameters at which they found a hyperchaos.

A significant novelty of the article is the high-precision method for searching for the Poincaré recurrences on attractors of dynamical systems with a quadratic right-hand side, also a modification of Benettin's algorithm for finding the Lyapunov exponents in order to check the found limiting regimes.

The methods described in the article can be used to study not only periodic regimes, but also quasiperiodic and chaotic regimes for dynamical systems with a quadratic right-hand side. According to the statistics of Poincaré recurrences, we can determine the dimension of the attractor and compare it with the dimension obtained when calculating the Lyapunov exponents (for example, the articles [START_REF] Anishchenko | Poincaré recurrences time and local dimension of chaotic attractors[END_REF][START_REF] Anishchenko | Local and global approaches to the problem of Poincaré recurrences[END_REF]).

The advantage of the considered algorithm of calculating the Lyapunov exponents is combination of the linearized system of equations and the original dynamical system in general form [START_REF] Afraimovich | Fractal and multifractal properties of exit times and Poincaré recurrences[END_REF] for finding the state and perturbation vectors together.

x 1 8 )

 18 (0) = 10, x 2 (0) = 50 • sin 10 ≈ -27.2011, x 3 (0) = 10, x 4 (0) = 10. (Let us verify it. Using the numerical method FGBFI, we construct the trajectory arc from the point (8) on the time interval [0, T], T = 15 getting close to the attractor. Let us indicate the parameters of the numerical method: b m = 128, ε m = 5.87747 • 10 -39 , ε pw = 10 -20 .

Figure 1 .

 1 Figure 1. Projection of the cycle of the system (1) with the period T c ≈ 0.3655 into the space x 1x 2x 3 .

Figure 2 .

 2 Figure 2. Projection of the cycle of the system (1) with the period T c ≈ 0.3655 into the space x 1x 2x 4 .

⟨

  ), Ψi = [ ψi,1 . . . ψi,2n ] Qp Υj , Υi-j ⟩, p = 1; 2n.

  converges for t ∈ [τ, τ], where τ = τ( Υ0 ) = 1 h2 ( Υ0 ) + δ , δ is any positive number. Let us consider a modification of Benettin's algorithm for calculating the Lyapunov exponents using the Gram-Schmidt orthogonalization [10, pp. 163-165]: 1.

12 .

 12 i) ; LE p := LE p + log |S|; If k < M, then Goto step 9; 13. LE p := LE p M • τ M , p = 1; n; 14. Print LE p , p = 1; n.

Figure 3 .Figure 4 .

 34 Figure 3. The graph of the change in time of the Lyapunov exponent LE 1 .

Figure 5 .Figure 6 .

 56 Figure 5. The graph of the change in time of the Lyapunov exponent LE 3 .

10 .

 10 Calculate the approximate value X(∆t) by summing the terms of the series (3) to such a value i, where the following inequality holds∥Υ i ∥ • |∆t| i < ε pw ; (6)11. Print way • t, X(∆t); 12. If X(∆t) / ∈ B a , then we got out the compact B a . Then write "Decrease the value ε pw and/or ε m "; ended := true; 13. If t = T, then ended := true; 14. If ended, then terminate the algorithm; 15. X(0) := X(∆t); 16. Goto step 6.

Table 1 :

 1 The Poincaré recurrences for the point[START_REF] Afanas'ev | On recurrent trajectories, minimal sets, and quasiperiodic motions of dynamical systems[END_REF], ∆t P = 10 -4 , T = 15 for reaching the attractor.

	Time moment t k *	Value d k *
	0.3655	0.0423764
	0.7310	0.0157840
	1.0965	0.0323213
	1.4620	0.0196461
	1.8275	0.0173865
	2.1930	0.0209083
	2.5585	0.0135900
	2.9240	0.0206764
	3.2895	0.0194286
	3.6550	0.0197196
	4.0205	0.0237665
	4.3860	0.0213908
	4.7515	0.0263917
	5.1170	0.0259931
	5.4825	0.0288078
	5.8480	0.0311426
	6.2135	0.0320319
	6.5790	0.0356643
	6.9445	0.0362604
	7.3101	0.0340862
	7.6756	0.0325096
	8.0411	0.0305410
	8.4066	0.0280697
	8.7721	0.0268155
	9.1376	0.0243143
	9.5031	0.0231314
	9.8686	0.0211890

Table 2 :

 2 The Poincaré recurrences for the point[START_REF] Kuznetsov | Dynamical Chaos[END_REF], ∆t P = 10 -4 , T = 40 for reaching the attractor.

	Time moment t k *	Value d k *
	0.3655	0.00133861
	0.7310	0.00267728
	1.0965	0.00401586
	1.4620	0.00535439
	1.8275	0.00669292
	2.1930	0.00803134
	2.5585	0.00936979
	2.9240	0.01070810
	3.2895	0.01204650
	3.6550	0.01338470
	4.0205	0.01472300
	4.3860	0.01606120
	4.7515	0.01739930
	5.1171	0.01774010
	5.4826	0.01640090
	5.8481	0.01506170
	6.2136	0.01372250
	6.5791	0.01238330
	6.9446	0.01104430
	7.3101	0.00970522
	7.6756	0.00836622
	8.0411	0.00702727
	8.4066	0.00568835
	8.7721	0.00434949
	9.1376	0.00301067
	9.5031	0.00167189
	9.8686	0.00033316

Table 3 :

 3 The Poincaré recurrences for the point[START_REF] Kuznetsov | Dynamical Chaos[END_REF], ∆t P = 10 -6 , T = 40 for reaching the attractor.

	Time moment t k *	Value d k * • 10 -5
	0.365504	1.105410
	0.731007	1.096510
	1.096510	0.348186
	1.462010	0.752419
	1.827520	1.799740
	2.193020	0.697982
	2.558530	0.397935
	2.924030	1.498770
	3.289530	1.051780
	3.655040	0.048264
	4.020540	1.146250
	4.386040	1.403530
	4.751550	0.304245
	5.117050	0.793828
	5.482550	1.754720
	5.848060	0.656571
	6.213560	0.442487
	6.579070	1.541010
	6.944570	1.008280
	7.310070	0.091034
	7.675580	1.189060
	8.041080	1.359960
	8.406580	0.261577
	8.772090	0.837369
	9.137590	1.711900
	9.503100	0.613282
	9.868600	0.485547

Acknowledgments

The reported study was funded by RFBR for the research project No. 20-01-00347.