© Revista de Matemática: Teoría y Aplicaciones 2021 28(2) : 183–208 cimpa – ucr issn: 1409-2433 (Print), 2215-3373 (Online) doi: 10.15517/rmta.v28i2.36084

bifurcaciones en modelo gause depredador-presa con discontinuidad

bifurcations in model gause predator-prey with discontinuity

Christian Cortés–García^{∗}

Received: 10/Oct/2020; Revised: 19/Jan/2021; Accepted: 11/Mar/2021

∗Universidad Carlos III de Madrid, Departamento de Matemática - Departamento de Biología de Sistemas, Centro de Investigación en Biotecnología, Madrid, España. E-Mail:

chcortes@math.uc3m.es cc.cortes@cnb.csic.es

184 c. cortés-garcía

Resumen

En este trabajo se presentan las condiciones necesarias para garantizar la existencia de un ciclo límite estable en un modelo de Gause depredador

- presa y algunos aspectos geométricos para realizar un análisis cualitativo en sistemas dinámicos de Filippov bidimensional. Con esos lineamientos definidos, se estudia la dinámica de un modelo depredador–presa cuando la explotación en los depredadores es restringida si la cantidad de presas es inferior a un valor critico. El estudio es llevado a cabo por el análisis de bifurcación con relación a dos parámetros: explotación y protección de las poblaciones a interactuar.

Palabras clave: sistemas planares de Filippov; análisis de bifurcación; ciclo límite; modelo depredador–presa; ciclo de Canard.

Abstract

This paper presents the necessary conditions to guarantee the existence of a stable limit cycle in a predator - prey model and some geometrical aspects to perform a qualitative analysis in two - dimensional Filippov dynamic systems. With these defined guidelines, the dynamics of a predator

- prey model are studied when exploitation in predators is restricted if the number of prey is lower than a critical value. The study is carried out by the bifurcation analysis in relation to two parameters: exploitation and protection of the populations to interact.

Keywords: planar systems Filippov; bifurcation analysis; limit cycles predatorprey systems; cicly Canard.

Mathematics Subject Classification: 34A36, 34C23, 34D20, 34D23, 92D25.

1 Introducción

Muchos modelos ecológicos pueden ser modelados como sistemas dinámicos discontinuos o sistemas no suaves. Uno de estos modelos son conocidos como sistemas de Filippov bidimensionales, descritos por ecuaciones diferenciales de la forma:

ẋ = f_{i}(x), x ∈ S_{i} ⊂ ℝ^{2}  (1)

con f_{1}(x) y f_{2}(x) son funciones continuas, S_{i}, i = 1, 2, son regiones abiertas, no sobrepuestas, separadas por una curva diferenciable Σ = f^{−1}(0), y cuya unión cubren el retrato de fase.

Si las componentes normales f_{1}(x) y f_{2}(x) con respecto a f poseen el mismo signo, sus trayectorias cruzan Σ. Por el contrario, cuando los campos vectoriales f_{1}(x) y f_{2}(x) poseen signo opuesto, el estado del sistema es forzado a escapar o a deslizarse sobre Σ. Como las trayectorias sobre Σ no están definidas, se plantean varios métodos. Una de ellas, es el método convexo de Filippov, donde las soluciones en Σ son dadas por ẋ = Z^{s}(x), con Z^{s}(x) una combinación convexa de f_{1}(x) y f_{2}(x) tangente a Σ en x.

Por otro lado, el análisis cualitativo para los sistemas de Filippov, que además de las bifurcaciones genéricas en sistemas continuos, podrían presentar bifurcaciones llamadas deslizantes, donde las variaciones en el parámetro de bifurcación provocan alteraciones sobre la dinámica de Σ. Todas las posibles bifurcaciones en sistemas de Filippov bidimensional fueron listadas por Kuznetsov et al. [9], utilizando el abordaje clásico de equivalencia topologica propuesta por Arnold [1]. Con este abordaje, todas las bifurcaciones, genéricas y deslizantes, pueden ser interpretadas como colisiones entre las trayectorias, puntos con características especiales y/o alteraciones sobre la dinámica de Σ.

Por tal motivo, en este trabajo se estudia la teoría básica de los sistemas planares de Filippov con el fin de realizar un análisis de bifurcación para un modelo discontinuo tipo depredador–presa, equivalente al presentado por Yang [12], cuando los parámetros de protección y explotación de las especies son alteradas.

Para comprender este proceso, en los preliminares se muestran las condiciones necesarias para que un modelo de Gause depredador-presa presenta un único ciclo límite globalmente asintóticamente estable, probado por Kuang & Freedman [7]. De igual forma se realiza una descripción, como lo plantea Buzzi et al. [2] y Guardia et al. [6], de los elementos básicos tales como puntos de equilibrio, puntos tangentes, trayectorias y órbitas periódicas en sistemas planares de Filippov. Por último, se analizan las bifurcaciones locales y globales de un modelo depredador–presa con discontinuidad, equivalente al presentado por Yang et al. [12], e interacción de dos parámetros: explotación y protección de las poblaciones a interactuar. Cabe resaltar que a diferencia del trabajo de Yang et al. [12], se han abarcado todos los posibles casos de bifurcaciones de codimensión dos junto con las curvas cuya dinámica del modelo planteado es alterada.

2 Preliminares

2.1 Sistemas de Filippov bidimensional

Consideren X e Y campos vectoriales de clase C^{r}, con r > 1, definidos en un conjunto abierto U ⊂ ℝ^{2} tal que (0, 0) ∈ U. Sea f una función de clase C^{r}, r > 1, tal que gradf(x, y) ≠ 0 para todo (x, y) ∈ U y

Σ = f^{−1}(0) ∩ U = {(x, y) ∈ U : f(x, y) = 0}

una curva abierta, diferenciable y que divide al conjunto U en dos regiones abiertas

Σ^{+} = {(x, y) ∈ U : f(x, y) > 0} y Σ^{−} = {(x, y) ∈ U : f(x, y) < 0}

con \overline{Σ^{+}} y \overline{Σ^{−}} sus adherencias.

Según Buzzi et al. [2] y Guardia et al. [6], un sistema planar de Filippov Z = (X, Y) es un campo vectorial definido por

Z(x, y) = \{ \begin{matrix}X(x, y), & (x, y) ∈ Σ^{+} \\ Y (x, y), & (x, y) ∈ Σ^{−} \end{matrix}

donde X e Y son de clase C^{r}, r > 1, en \overline{Σ^{+}} y \overline{Σ^{−}}, respectivamente.

(a) Condición Xf(p), Y f(p) > 0.  (b) Condición Xf(p), Y f(p) < 0.

Figura 1: Segmento de costura Σ^{c}.

(a) Segmento deslizante.  (b) Segmento de Escape.

Figura 2: Segmento Σ^{s} e Σ^{e}.

Para aquellos puntos p ∈ Σ^{+} o p ∈ Σ^{−}, la trayectoria en Z = (X, Y) está definida por los campos vectoriales X o Y, respectivamente, sin embargo, se debe establecer una trayectoria para los puntos p ∈ Σ. Para ello, se divide Σ en tres regiones disyuntas dados por:

• segmento de costura: Σ^{c} = {p ∈ Σ : Xf(p) · Y f(p) > 0} como se observa en la Figura 1,

• segmento deslizante: Σ^{s} = {p ∈ Σ : Xf(p) < 0, Y f(p) > 0} dada por la Figura 2(a),

• segmento de escape: Σ^{e} = {p ∈ Σ : Xf(p) > 0, Y f(p) < 0} dada por la Figura 2(b),

con Xf(p) = X(p) · gradf(p) la derivada de Lie de f con respecto al campo vectorial X en p.

Se denomina p ∈ Σ como punto tangente si Xf(p) = 0 o Y f(p) = 0, los cuales aparecen en la frontera de Σ^{c}, Σ^{s} o Σ^{e}, denotadas por ∂Σ^{c}, ∂Σ^{s} o ∂Σ^{e}, respectivamente. También p ∈ Σ es punto tangente si X(p) = 0 o Y (p) = 0, esto es, cuando uno de los dos campos vectoriales tiene un punto crítico en Σ.

Los puntos tangentes p ∈ Σ se pueden clasificar en cuadrático si Xf(p) = 0 y X^{2}f(p) = X(p) · gradXf(p) ≠ 0, o cúbico si Xf(p) = X^{2}f(p) = 0 y X^{3}f(p) = X(p) · gradX^{2}f(p) ≠ 0.

Se dice que un punto cuadrático p ∈ Σ es regular si Xf(p) = 0, X^{2}f(p) ≠ 0 y Y f(p) ≠ 0; o Y f(p) = 0, Y^{2}f(p) ≠ 0 y Xf(p) ≠ 0. Para el primer caso, un punto cuadrático regular es visible si X^{2}f(p) > 0 e invisible si X^{2}f(p) < 0. Para el segundo caso, p ∈ Σ es visible si Y^{2}f(p) < 0 e invisible si Y^{2}f(p) > 0.

Como lo observado en la Figura 1, la trayectoria en p ∈ Σ^{c} apuntan en la misma dirección de las trayectorias en X e Y. Sin embargo, en Σ^{s} o Σ^{e}, su trayectoria es dada por una combinación convexa de los campos vectoriales X e Y, esto es,

Z^{s}(p) = λ(p)X(p) + (1 − λ(p))Y (p).

En vista de la Figura 3,

Z^{s}(p) · gradf(p) = 0,

entonces

λ(p) = \frac{Y f(p)}{Y f(p) − Xf(p)},

Rev.Mate.Teor.Aplic. (ISSN print: 1409-2433; online: 2215-3373) Vol. 28(2): 183–208, Jul–Dec 2021

Figura 3: Construcción de trayectorias Z^{s}(p).

y así el campo vectorial deslizante Z^{s} es dado por

Z^{s}(p) = \frac{1}{Y f(p) − Xf(p)}(Y f(p)X(p) − Xf(p)Y (p)),

definido en Σ^{e} ∪ Σ^{s} y la trayectoria para p ∈ Σ^{s} ∪ Σ^{e} es dada por la solución de ese campo vectorial. Por lo tanto, las trayectorias en el sistema Z = (X, Y) son definidas por:

Definición 1 Sean φ_{X} y φ_{Y} las trayectorias de los campos vectoriales X y Y de clase C^{r}, r > 1, en U, respectivamente. La trayectoria φ_{Z}(t, p) de un sistema de Filippov Z = (X, Y) para un punto p ∈ U es definida por:

• Para p ∈ Σ^{+} con X(p) ≠ 0 es dada por φ_{Z}(t, p) = φ_{X}(t, p). De igual manera, para p ∈ Σ^{−} con Y (p) ≠ 0 la trayectoria es φ_{Z}(t, p) = φ_{Y} (t, p), respectivamente, donde t ∈ I ⊂ ℝ. Este tipo de trayectorias se denominan regulares.

• Para p ∈ Σ^{c} tal que Xf(p) > 0, Y f(p) > 0 y al tomar t = 0 en p, la trayectoria es dada por φ_{Z}(t, p) = φ_{Y} (t, p) para t ∈ I ∩ {t ≤ 0} e φ_{Z}(t, p) = φ_{X}(t, p) para t ∈ I ∩ {t ≥ 0}. Para el caso Xf(p) < 0, Y f(p) < 0, la trayectoria es dada por φ_{Z}(t, p) = φ_{Y} (t, p) para t ∈ I ∩ {t ≥ 0} y φ_{Z}(t, p) = φ_{X}(t, p) para t ∈ I ∩ {t ≤ 0},

• Para p ∈ Σ^{e} ∪ Σ^{s} tal que Z^{s}(p) ≠ 0, φ_{Z}(t, p) = φ_{Z^{s}}(t, p) para t ∈ I ⊂ ℝ,

• Para p ∈ ∂Σ^{c} ∪ ∂Σ^{s} ∪ ∂Σ^{e}, la trayectoria por p es φ_{Z}(t, p) = {p}.

• Para cualquier otro punto de tangencia en Σ o puntos de equilibrio en X sobre Σ^{+}, Y en Σ^{−} o Z^{s} en Σ^{s} ∪ Σ^{e}, φ_{Z}(t, p) = {p} para todo t ∈ ℝ.

Se denomina pseudo-equilíbrio en Z = (X, Y) al punto p ∈ Σ^{s} ∪ Σ^{e} tal que Z^{s}(p) = 0, los cuales pueden ser pseudo-nodo estable si p ∈ Σ^{s} y (Z^{s})′(p) < 0, como lo observado en la Figura 4(a), pseudo-nodo inestable si p ∈ Σ^{e} y (Z^{s})′(p) > 0, como lo observado en la Figura 4(b), y una pseudo-silla si p ∈ Σ^{s} y (Z^{s})′(p) > 0 o p ∈ Σ^{e} y (Z^{s})′(p) < 0, como lo observado en la Figura 4(c).

(a) Pseudo-nodo estable.

(b) Pseudo-nodo inestable.

(c) Pseudo-silla

Figura 4: Ejemplos de Pseudo-equilibrios del Sistema Z = (X, Y).

Como lo observado en la Figura 5(a), una órbita Γ = {φ_{Z}(t, p) : t ∈ ℝ ⊂ Σ^{+} ∪ Σ^{−} ∪ \overline{Σ^{c}} es periódica si φ_{Z}(t + T, p) = φ_{Z}(t, p) para algún T > 0. De igual forma, se denomina ciclo límite en Z = (X, Y) a una órbita periódica Γ que esta contenida en Σ^{+} o Σ^{−}.

(a) Órbita periódica  (b) Ciclo  (c) Pseudo-ciclo

Figura 5: Ejemplos de una órbita periódica, ciclo y un pseudo-ciclo.

A diferencia de los sistemas dinámicos continuos, las órbitas periódicas en los sistemas de Filippov pueden presentar los siguientes casos:

Definición 2 Un ciclo es la unión de un conjunto finito de trayectorias γ_{1}, …, γ_{n} tal que γ_{2k} es un trayectoria contenida en el segmento deslizante Σ^{s} y γ_{2k+1} una trayectoria contenida en Σ^{+} o Σ^{−} el cual sus puntos de llegada y partida pertenecen a \overline{γ_{2k}} y \overline{γ_{2k+2}}, respectivamente. En la Figura 5(a) se observa un ciclo para n = 2.

Definición 3 Un pseudo-ciclo es la unión de un conjunto de trayectorias γ_{1}, …, γ_{n}, contenidas Σ^{+} o Σ^{−}, tal que una de las extremidades de γ_{i} coincide con una extremidad de γ_{i−1} y la otra con una extremidad de γ_{i+1}. En la Figura 5(b) se observa un pseudo-ciclo para n = 2.

El siguiente resultado, probado por Buzzi et al. [2], muestran condiciones necesarias que garantiza la existencia de un ciclo para un sistema de Filippov Z = (X, Y).

Teorema 4 Considere Z = (X, Y) un sistema de Filippov con un único punto cuadrático regular visible A ∈ Σ. Sea γ_{1} la trayectoria en Σ^{+} o Σ^{−} que pasa a través del punto A y B como punto transversal de γ_{1} con Σ. Entonces, Z posee un ciclo Γ, como lo observado en la Figura 5(b), si y solo si satisfacen las siguientes condiciones:

(i) La componente γ_{1} no posee puntos cuadrático regulares entre A y B,

(ii) Xf(p) · Y f(p) < 0 para todo p ∈ (A, B],

(iii) {X, Y } es un conjunto linealmente independiente em [A, B].

2.2 Modelo Gause depredador–presa

Un modelo de Gause depredador–presa es descrito como

\{ \begin{matrix}ẋ = xg(x) − yp(x) \\ ẏ = y[−γ + q(x)], \end{matrix}  (2)

donde x(t) ≥ 0 e y(t) ≥ 0 representan la cantidad de presas y depredadores, respectivamente para un t ≥ 0, de tal forma que satisfacen las siguientes condiciones:

• En ausencia de los depredadores, el aumento de las presas converge para alguna constante K > 0. Esto es, para ẋ = xg(x) satisface g(x) > 0 para 0 ≤ x < K, g(K) = 0 y g(x) < 0 para x > K.

• Los depredadores reducen la tasa de crecimiento al reducir el número de presas. Esto es, yp(x) donde p(x) es la cantidad de presas muertas por un depredador tal que p(0) = 0 y p(x) > 0 para x > 0.

• La tasa de crecimiento en los depredadores es dada por −γ + q(x). El parámetro γ > 0 corresponde a la mortalidad de los depredadores en ausencia de presas y q(x) es una función positiva creciente que satisface q(0) = 0 y q′(x) > 0 para x > 0.

Para garantizar la existencia de un único ciclo límite globalmente asintóticamente estable en el modelo (2) se usa el siguiente resultado, cuya prueba fue realizada por Kuang & Freedman [7].

Teorema 5 Un modelo de Gause

\{ \begin{matrix}ẋ = xg(x) − yp(x) \\ ẏ = y[−γ + q(x)], \end{matrix}

con γ > 0 y g, p, q funciones de clase C^{1}, que satisfacen

(H1) Existe un K > 0 tal que g(x) > 0 si 0 ≤ x < K, g(K) = 0 y g(x) < 0 si x > K,

(H2) p(0) = 0 y p(x) > 0 para x > 0,

(H3) q(0) = 0, q′(x) > 0,

(H4) Existe x^{∗}, 0 < x^{∗} < K tal que q(x^{∗}) = γ,

(H5) \frac{d}{dx} (\frac{xg(x)}{p(x)}) |_{x=x^{∗}} > 0,

(H6) \frac{d}{dx} (^{\frac{xg′(x)+g(x)−xg(x) \frac{p′(x)}{p(x)}}{−γ+q(x)}} ) ≤ 0 para 0 ≤ x < x^{∗} y x^{∗} < x ≤ K,

posee un único punto de equilibrio (x^{∗}, y^{∗}) localmente inestable dentro del primer cuadrante y un único ciclo límite globalmente asintóticamente estable en {(x, y) : x > 0, y > 0} \ {(x^{∗}, y^{∗})}.

3 Modelo discontinuo depredador–presa

Sean x ≥ 0 e y ≥ 0 la cantidad de presas y depredadores, respectivamente, donde solamente los depredadores pueden ser capturados. Si la población de presas cae por debajo de un valor critico P, esto es, x < P, entonces se activa un mecanismo de protección que impide que estas sean consumidas por los depredadores. Este modelo es descrito de la siguiente forma:

\begin{cases}ẋ = x [ρ (1 − \frac{x}{K} )] \\ ẏ = y [−d − qE] , \end{cases}  (3)

donde • f(x) = x [ρ (1 − \frac{x}{K} )] describe el crecimiento logístico de las presas, con

ρ > 0 la tasa de crecimiento sin consumo por los depredadores y K > P la capacidad de soporte medio para las presas en ausencia de depredadores.

• h(y) = y[−d − qE] describe el decrecimiento de los depredadores, con d > 0 la tasa de mortalidad, q > 0 el coeficiente de capturabilidad y E ≥ 0 el esfuerzo de capturabilidad.

Por otro lado, cuando el número de presas es superior a P, esto es, x > P, su mecanismo de protección se desactiva y pasan a ser consumidas por los depredadores. Este modelo es descrito de como:

\{ \begin{matrix}ẋ & = & x [ρ (1 − \frac{x}{K} ) − \frac{ay}{b + x}] \\ ẏ & = & y [\frac{cax}{b + x} − d − qE] , \end{matrix}  (4)

donde

• p(x) = \frac{ax}{b+x} es el aumento desacelerado en el consumo una vez que aumentan las presas consumidas hasta aproximarse asintóticamente a una constante a, donde a > 0 es la taza de consumo de los depredadores y b > 0 la constante media de captura,

• −yp(x) es el consumo de las presas por unidad de tiempo para cada depredador,

• cyp(x) es la conversión de presas consumidas en nuevos depredadores, donde c > 0 es la taza de eficiencia con que el depredador capturó la presa.

Finalmente, este modelo, equivalente a lo mostrado por Yang et al. [12], es descrito como un sistema planar de Filippov:

Z(x, y) = \{ \begin{matrix}X(x, y) = (\begin{matrix}x [ρ (1 − \frac{x}{K} ) − \frac{ay}{b + x}] \\ y [\frac{cax}{b + x} − d − qE] \end{matrix}) & si & x > P, \\ Y (x, y) = (\begin{matrix}x [ρ (1 − \frac{x}{K} )] \\ y [−d − qE] \end{matrix}) & si & 0 ≤ x < P, \end{matrix}

(5)

donde P < K y

\begin{aligned}Σ^{+} & = & \{(x, y) ∈ ℝ^{2} : y ≥ 0, f(x, y) = x − P > 0\}, \\ Σ & = & \{(x, y) ∈ ℝ^{2} : y ≥ 0, f(x, y) = x − P = 0\}, \\ Σ^{−} & = & \{(x, y) ∈ ℝ^{2} : x ≥ 0, y ≥ 0, f(x, y) = x − P < 0\}. \end{aligned}  (6)

Para analizar la dinámica sobre aquellas alteración en los parámetros de protección P y capturabilidad E de las especies en el modelo (5), es necesario encontrar y determinar la estabilidad local de los puntos de equilibrio, y posibles ciclos límites, en los campos vectoriales X e Y, representados por los modelos 3 y 4, respectivamente.

3.1 Estudio del campo vectorial Y

Los puntos de equilibrio de Y son dados por (0, 0) y (K, 0), cuya matriz Jacobiana DY del campo vectorial Y es dada por

DY (x, y) = [\begin{matrix}\frac{ρ(K−2x)}{K} & 0 \\ 0 & −d − qE \end{matrix}] .  (7)

Como

detDY (0, 0) = −ρ(d + qE) < 0,  (8)

entonces, por el Teorema de Hartman Grobman dado en [11], el punto de equilibrio (0, 0) es localmente un punto silla.

Para el punto de equilibrio (K, 0) se tiene

detDY (K, 0) = ρ(d + qE) > 0,
trDY (K, 0) = −(ρ + d + qE) < 0,
△ = [trDY (K, 0)]^{2} − 4detDY (K, 0) = 2d^{2} + 2qdE + (ρ − dE)^{2} > 0,

del cual (K, 0) es localmente un nodo estable. Se observa que el campo Y no posee ciclos límites en el primer cuadrante dado que ẏ = −(d + qE)y ≤ 0 para todo y ≥ 0.

3.2 Estudio del campo vectorial X

Los puntos de equilibrio de X son (0, 0), (K, 0) y

(x^{∗}, y^{∗}) = (\frac{b(d+qE)}{ac−(d+qE)}, \frac{bcρ[acK−(b+K)(d+qE)]}{K[ac−(d+qE)]^{2}}) .

Por la descripción del modelo, (x^{∗}, y^{∗}) debe ser un punto interior en el primer cuadrante.

Como d + qE > 0 y a > 0, se necesita que ac − (d + qE) > 0 para que x^{∗} > 0. Además y^{∗} > 0 si E < \frac{acK−d(b+K)}{q(b+K)}. Por tanto, (x^{∗}, y^{∗}) se encuentra

en el primer cuadrante si E < \frac{acK−d(b+K)}{q(b+K)} y ac > d + qE.

Por otro lado, como los valores propios de la matriz Jacobiana de X en (0, 0) son dados por λ_{1} = ρ > 0 y λ_{2} = −(d − qE) < 0, entonces el punto de equilibrio (0, 0) es localmente un punto silla. Ademas, de los valores propios de la matriz Jacobiana de X evaluada en (K, 0),

\begin{aligned}λ_{1} & = & −ρ < 0, \\ λ_{2} & = & \frac{acK−(b+K)(d+qE)}{b+K}, \end{aligned}

si λ_{2} > 0, equivalente a E < \frac{acK−d(b+K)}{q(b+K)}, el punto (K, 0) es localmente un

punto silla. Sin embargo, si λ_{2} < 0, esto es, \frac{acK−d(b+K)}{q(b+K)} < E, el punto (K, 0) es localmente estable y el punto (x^{∗}, y^{∗}) no coexistiría.

El determinante y la traza de la matriz Jacobiana de X en (x^{∗}, y^{∗}) son

\begin{aligned}detDX(x^{∗}, y^{∗}) & = & \frac{ρ(d + qE)[acK − (b + K)(d + qE)]}{acK}, \\ trDX(x^{∗}, y^{∗}) & = & \frac{ρ(d + qE)[ac(b − K) + (b + K)(d + qE)]}{acK[(d + qE) − ac]}. \end{aligned}

Como acK > (b + K)(d + qE) tenemos que detDX(x^{∗}, y^{∗}) > 0 y, si x^{∗} < \frac{K−b}{2} entonces ac(b−K)+(b+K)(d+qE) < 0 y así trDX(x^{∗}, y^{∗}) > 0. Análogamente, si x^{∗} > \frac{K−b}{2} entonces ac(b − K) + (b + K)(d + qE) > 0 y trDX(x^{∗}, y^{∗}) < 0.

Si △ = [trDX(x^{∗}, y^{∗})]^{2} − 4detDX(x^{∗}, y^{∗}) se concluye que la estabilidad local de (x^{∗}, y^{∗}) es

• Un nodo estable si x^{∗} > \frac{K−b}{2} y △ ≥ 0,
• Un nodo inestable si x^{∗} < \frac{K−b}{2} y △ ≥ 0,

• Un foco estable si x^{∗} > \frac{K−b}{2} y △ < 0,

• Un foco inestable si x^{∗} < \frac{K−b}{2} y △ < 0.

El siguiente resultado muestra la existencia del ciclo límite del campo vectorial X.

Lema 1 Si x^{∗} ≤ \frac{K−b}{2} y △ < 0, entonces (x^{∗}, y^{∗}) es un foco inestable y existe un único ciclo límite globalmente asintóticamente estable Γ_{X} de X en el primer cuadrante.

Prueba. Por hipótesis, (x^{∗}, y^{∗}) es un foco inestable.

Dado que campo X es equivalente al modelo 2, con

g(x) = ρ (1 − \frac{x}{K} ) , p(x) = \frac{ax}{b + x}, q(x) = \frac{cax}{b + x}, γ = d + qE,

para garantizar la existencia de un único ciclo límite globalmente asintóticamente estable en el campo vectorial X, se usa el Teorema 5. Nótese que,

i) g(x) > 0 si 0 ≤ x < K, g(K) = 0 y g(x) < 0 si x > K,

ii) p(0) = 0, p′(x) = \frac{ab}{(x+b)^{2}} > 0,

iii) q(0) = 0, q′(x) = \frac{abc}{(x+b)^{2}} > 0,

iv) Como

\frac{xg′(x) + g(x) − xg(x)\frac{p′(x)}{p(x)}}{−γ + q(x)} = \frac{−\frac{ρx}{K} + ρ (1 − \frac{x}{K} ) − ρ (1 − \frac{x}{K} ) \frac{b}{x+b}}{\frac{cax}{b+x} − d − qE}
= \frac{ρx(2x + b − K)}{K \{b(d + qE) − x[ac − (d + qE)]\}}
= \frac{ρx(2x + b − K)}{K[ac − (d + qE)](x^{∗} − x)},

al calcular su derivada se tiene que:

\frac{1}{ac − (d + qE)(x^{∗} − x)^{2}} (−\frac{2ρ}{K}x^{2} + \frac{4ρx^{∗}}{K} x + ρx^{∗} + \frac{ρbx^{∗}}{K} ) .

Ya que

\frac{1}{ac − (d + qE)(x^{∗} − x)^{2}} ≥ 0,

Rev.Mate.Teor.Aplic. (ISSN print: 1409-2433; online: 2215-3373) Vol. 28(2): 183–208, Jul–Dec 2021

considere

F (x) = −\frac{2ρ}{K}x^{2} + \frac{4ρx^{∗}}{K} x + ρx^{∗} + \frac{ρbx^{∗}}{K} ,

esto es, una parábola que abre hacia abajo con discriminante

△_{F} = \frac{8ρ^{2}x^{∗}}{K^{2}} (b − K + 2x^{∗}) .

Si x^{∗} ≤ \frac{K−b}{2} se tiene que b − K + 2x^{∗} ≤ 0, esto es, △_{F} ≤ 0. De igual forma F (x) posee a lo mucho una raíz real que garantiza F (x) ≤ 0 y

\frac{d}{dx} (\frac{xg′(x) + g(x) − xg(x)\frac{p′(x)}{p(x)}}{−γ + q(x)} ) ≤ 0.

Por lo tanto, el campo X posee un único ciclo límite globalmente asintóticamente estable en el primer cuadrante.

3.3 Estudio del sistema de Filippov Z

Para todo p = (P, y) ∈ Σ con y ≥ 0, se tiene que

Xf(p) = X(p) · gradf(p) = P [ρ (1 − \frac{P}{K} ) − \frac{ay}{b + P} ] ,
Y f(p) = Y (p) · gradf(p) = P ρ (1 − \frac{P}{K} ) .

Como P < K, entonces Y f(p) = P ρ (1 − \frac{P}{K} ) > 0 y si

Xf(p) = P [ρ (1 − \frac{P}{K} ) − \frac{ay}{b + P} ] > 0,

entonces

y_{P} = \frac{ρ (K − P) (b + P)}{aK} > y.

Por lo tanto,

\begin{aligned}Σ^{s} & = & \{(x, y) ∈ ℝ^{2} : x = P, y_{P} < y\} , \\ Σ^{c} & = & \{(x, y) ∈ ℝ^{2} : x = P, y_{P} > y\} , \\ Σ^{e} & = & ∅. \end{aligned}

Rev.Mate.Teor.Aplic. (ISSN print: 1409-2433; online: 2215-3373) Vol. 28(2): 183–208, Jul–Dec 2021

El campo vectorial del segmento deslizante Z^{s}(p), con p ∈ Σ^{s}, es dado por

Z^{s}(p) = ( 0 ; \frac{P cρ(K − P) − Ky(d + qE)}{K}) .

El pseudo-equilibrio del sistema de Filippov (5)

P N = (P, y_{1}) = (P, \frac{ρcP (K − P)}{K(d + qE)}) ,

existe si y_{1} > y_{p} pues P N ∈ Σ^{s}.

Por otro lado, si x^{∗} > P, esto es, E > \frac{acP −d(b+P)}{q(b+P)}, entonces (x^{∗}, y^{∗}) es un equilibrio del sistema de Filippov (5). Sin embargo, si x^{∗} < P, esto es E < \frac{acP −d(b+P)}{q(b+P)}, entonces (x^{∗}, y^{∗}) no pertenece al sistema de Filippov (5).

Se observa que para el segmento deslizante Z^{s}(p),

\frac{d}{dy} (\frac{P cρ(K − P) − Ky(d + qE)}{K}) = −(d + qE) < 0,

luego P N ∈ Σ^{s} es un pseudo-nodo estable. Como Y f(p) > 0 para todo p ∈ Σ, existe un punto tangente cuadrático cuando Xf(p) = 0, esto es,

T = (P, \frac{ρ(b + P)(K − P)}{aK}) ,

con

X^{2}f(T) = X(T) · gradXf(T) = \frac{P ρ(P − K)(acP − (b + P)(d + qE))}{K(b + P)}.

El punto T es cuadrático visible si X^{2}f(T) > 0, esto es, P < x^{∗}, e invisible si X^{2}f(T) < 0, esto es, P > x^{∗}. Por otro lado, si (x^{∗}, y^{∗}) es un equilibrio del sistema de Filippov (5), esto es, x^{∗} > P, entonces

y_{p} − y_{1} = \frac{ρ (K − P) (b + P)}{aK} − \frac{ρcP (K − P)}{K(d + qE)} = \frac{X^{2}f(T)}{aP} > 0,

así y_{p} > y_{1} del cual P N no existe en el sistema de Filippov (5). Análogamente, si P N existe, esto es, y_{1} > y_{P}, se tiene que acP − (b + P)(d + qE) > 0 y

P − x^{∗} = P − \frac{b(d + qE)}{ac − (d + qE)} = \frac{acP − (b + P)(d + qE)}{ac − (d + qE)} > 0,

así (x^{∗}, y^{∗}) no ésta definido en el sistema de Filippov (5). De estas deducciones se ha demostrado el siguiente resultado.

(a) First  (b) Second

Figura 6: Retrato de fase del sistema de Filippov (5) cuando (x^{∗}, y^{∗}) es un foco inestable, E = 0.2 y (a) P = 3.3; (b) P = 10.3.

Lema 2 El equilibrio (x^{∗}, y^{∗}) y el pseudo-equilibrio P N no pueden existir al mismo tiempo. Además, P N es estable si Σ^{s} existe.

Al realizar un análisis global al sistema de Filippov (5), equivalente a lo mostrado por Yang et al. [12], se procede a estudiar dos casos:

Caso uno: El punto (x^{∗}, y^{∗}) es un foco inestable, esto es, x^{∗} ≤ \frac{K−b}{2} y △ < 0, y Γ_{X} es un ciclo límite en X.

Sea Γ un ciclo límite del sistema de Filippov (5) en Σ^{+}, L = (x_{L_{X}} , y_{L_{X}}) y R = (x_{R_{X}} , y_{R_{X}}) puntos que se encuentran mas a la izquierda y derecha del ciclo límite Γ_{X} en X, respectivamente. Se cumplen los siguientes resultados.

Teorema 6 El sistema de Filippov (5) posee un único ciclo límite Γ globalmente asintóticamente estable si, y solo si 0 ≤ P < x_{L_{X}}.

Prueba. Si 0 ≤ P < x_{L_{2}}, entonces P < x^{∗} y (x^{∗}, y^{∗}) es un equilibrio en Σ^{+}, el punto tangente T es visible y el pseudo-equilibrio P N no existe. Además, el ciclo límite Γ_{X} en X se encuentra completamente al lado derecho de la línea x = P, esto es, en Σ^{+}. Luego Γ_{X} = Γ como se observa en la Figura 6(a).

Por otro lado, cualquier trayectoria que inicia en p ∈ Σ^{s} se desplaza por el segmento deslizante hasta T, y debido a la estabilidad de Γ_{X}, converge a Γ_{X}. Además, si la trayectoria parte en Σ^{+} y llega a Σ^{s}, se desliza hasta T a lo largo del segmento deslizante y converge a Γ_{X}. De igual forma, si la trayectoria parte en Σ^{+} mas no llega a Σ^{s}, debe converger a Γ_{X}. Por lo tanto, cualquier trayectoria que inicia en Σ^{+} converge a Γ_{X} = Γ.

Ahora, cualquier solución que inicia en Σ^{−} llega a x = P para un tiempo finito. Si las trayectorias llegan a Σ^{s}, se desliza hasta T y converge para Γ_{X}. Análogamente, si tal trayectoria llega a la línea x = P por debajo del punto tangente T, esto es, llega a Σ^{c}, entran a la región Σ^{+} y convergen a Γ_{X}. Por lo tanto, todas las soluciones del sistema de Filippov (5) finalmente se aproximan a Γ_{X} = Γ.

Recíprocamente, es claro que P < 0 no es posible. Si P = x_{L_{X}}, se tiene un ciclo tangente con Σ, mas no un ciclo límite Γ para el Sistema de Filippov (5). Se P > x_{L_{X}}, se obtiene un ciclo mas no un ciclo límite Γ en X. Por tanto 0 ≤ P < x_{L_{X}}.

Los siguientes resultados muestran las condiciones para garantizar la existencia de un ciclo.

Teorema 7 Existe un ciclo tangente a Σ globalmente asintóticamente estable en el sistema de Filippov (5) si, y solo si x_{L_{X}} = P.

Prueba. Si x_{L_{X}} = P, el ciclo límite Γ_{X} intercepta a Σ en T como se observa en la Figura 6(b). Luego existe un ciclo tangente Γ a Σ en el sistema de Filippov (5). Como Γ_{X} es globalmente estable, cualquier trayectoria que inicia en Σ^{+}, Σ o Σ^{−} debe converger a Γ. Por otro lado, P = x_{L_{X}} es una condición necesaria para garantizar la existencia y estabilidad global del ciclo tangente a Σ en el sistema de Filippov (5).

Teorema 8 El sistema de Filippov (5) posee un ciclo Γ globalmente asintóticamente estable si, y solo si x_{L_{X}} < P < x^{∗}.

Prueba. Si x_{L_{X}} < P < x^{∗}, entonces (x^{∗}, y^{∗}) es un equilibrio del sistema de Filippov (5) y el punto tangente T es visible. Además, el ciclo límite Γ_{X} no se encuentra completamente contenido en Σ^{+} y el pseudo-equilibrio P N no existe conforme al Lema 2. Para la existencia del ciclo, es necesario probar que la solución del sistema de Filippov (5) que comienza a partir del punto tangente T llega a un punto B ∈ Σ^{s}.

En efecto, la trayectoria γ_{1} que inicia en T rodea una parte del ciclo contenido en Σ^{+}, debido a la inestabilidad de (x^{∗}, y^{∗}) como se muestra en la Figura 7(a), y se encuentra con un punto B ∈ Σ^{s}. Se observa que B ≠ T, caso contrario, Γ sería un ciclo tangente a Σ. Además, γ_{1} es un arco tipo focal, pues no posee puntos tangentes en [B, T), y Xf(p) · Y f(p) < 0 en p ∈ (T, A] ⊂ Σ^{s}. Por el Lema 2, no existe el pseudo-equilibrio P N en Σ^{s}, indicando que {X(p), Y (p)} es un conjunto linealmente independiente para p ∈ [B, T].

(a)  (b)

(c)

Figura 7: Retrato de fase del sistema de Filippov (5) cuando (x^{∗}, y^{∗}) es un foco inestable, E = 0.2 y (a) P = 10 y (x^{∗}, y^{∗}) estable, (b) P = 20 y (x^{∗}, y^{∗}) estable, (c) P = 30 y P N estable.

De esa forma, por el Teorema 4, el sistema de Filippov (5) posee un ciclo dado por Γ = γ_{1} ∪ \overline{BT} como se observa en la Figura 7(a). Cualquier solución en el interior de Γ converge al segmento \overline{BT} ⊂ Σ^{s} pues (x^{∗}, y^{∗}) es un foco inestable. Además, cualquier solución que inicia en un punto por fuera de Γ debe converger para el ciclo Γ. Por consiguiente, Γ es un ciclo globalmente asintóticamente estable en el sistema Filippov (5).

Recíprocamente, se observa que x_{L_{X}} < P, caso contrario, la condición P ≤ x_{L_{X}} garantiza la existencia y estabilidad del ciclo tangente con Σ o ciclo límite Γ_{X} del sistema de Filippov (5), sin embargo contradice el hecho que Γ es un ciclo. Por otro lado P < x^{∗}, caso contrario, si P ≥ x^{∗}, entonces (x^{∗}, y^{∗}) no pertenece al sistema de Filippov (5), y por el Lema 2, el pseudo-equilibrio P N ∈ Σ^{s} existe y es estable. Esto contradice la estabilidad global del ciclo y la prueba ésta completa.

Teorema 9 El pseudo-equilibrio P N del sistema de Filippov (5) es globalmente asintóticamente estable si, y solo si, x^{∗} ≤ P < K.

Prueba. Si P N es globalmente asintóticamente estable, entonces (x^{∗}, y^{∗}) no ésta definido en el sistema de Filippov (5), esto es, x^{∗} ≤ P. Además, P < K. Por tanto x^{∗} ≤ P < K.

Recíprocamente, si x^{∗} ≤ P < K, entonces (x^{∗}, y^{∗}) no existe en el sistema de Filippov (5), el punto tangente T es invisible por Lema 1 y el pseudoequilibrio P N existe y es estable en Σ^{s} como se observa en la Figura 7(c). Por lo tanto, es suficiente mostrar que cualquier trayectoria que llega a x = P para y_{1} < y_{p} debe converger a Σ^{s}. En efecto, cualquier trayectoria que parte en Σ^{−} llega a Σ^{s} o Σ^{c}. Si entra en Σ^{c}, la trayectoria recorre Σ^{+} hasta llegar a Σ^{s}.

Por lo tanto, tal trayectoria se desliza sobre Σ^{s} y converge al pseudo-equilibrio P N. Luego P N es globalmente asintóticamente estable.

Caso dos: El punto (x^{∗}, y^{∗}) es un foco estable en Σ^{+}, esto es, x^{∗} > \frac{K−b}{2} y △ < 0.

Al observar la Figura 8(a) y la prueba del Teorema 6 se obtiene lo siguiente.

Teorema 10 Suponga que 0 ≤ P < x^{∗}, entonces (x^{∗}, y^{∗}) del sistema de Filippov (5) es globalmente asintóticamente estable.

Análogamente, la prueba del siguiente resultado es semejante al Teorema 9 y de observar la Figura 8(c).

Teorema 11 El pseudo-equilibrio P N en el sistema de Filippov (5) es globalmente asintóticamente estable si, y solo si, x^{∗} ≤ P < K.

3.4 Análisis de bifurcación

Dado la clasificación de las bifurcaciones en sistemas de Filippov presentados por Kuznetsov et al. [9], el sistema de Filippov (5) presenta una bifurcación Σ-foco cuando los puntos (x^{∗}, y^{∗}) y T colisionan simultáneamente siempre que P cruza a través de x^{∗}. Para P = x^{∗}, el ciclo estable se encoge y el equilibrio (x^{∗}, y^{∗}) con el punto tangente visible T colisionan y forman un punto estable

Parámetro Unidad  Valor
a  Por día  1
b  Por presa  5
c  Por día  0.45
d  Por día  0.2
ρ  Por día  1
K  Número por unidad de área 50
q Por día 0.8
P  Número por unidad de área 0 - 50

Tabla 1: Valores de los parámetros para el sistema de Filippov (5).

como se observa en la Figura 7(b). Cuando x_{L_{2}} < P < x^{∗}, el ciclo Γ pasa a través del punto tangente visible T y rodea el foco inestable (x^{∗}, y^{∗}).

Además, el foco inestable (x^{∗}, y^{∗}) y el punto tangente visible T coexisten como lo observado en la Figura 7(a). Para P > x^{∗}, el pseudo-equilibrio estable P N y el punto tangente invisible T aparecen en el sistema de Filippov (5) como lo observado en la Figura 7(c).

Existe otra bifurcación Σ-foco para el sistema de Filippov (5) cuando (x^{∗}, y^{∗}) es un foco estable y P = x^{∗} como se observa en la Figura 8(b). El foco estable (x^{∗}, y^{∗}) y el punto tangente visible T coexisten cuando 0 < P < x^{∗} como se observa en la Figura 8(a), colisionan si P = x^{∗} y forma un pseudo-equilibrio estable P N y el punto tangente invisible T como lo mostrado en la Figura 8(c). Esa bifurcación muestra como un foco estable se transforma en un pseudo-equilibrio estable.

De igual forma, el sistema de Filippov (5) presenta una bifurcación Σ-nodo cuando el nodo estable (x^{∗}, y^{∗}) colisiona con la línea x = P. En la Figura 9(a), cuando P = 20, la trayectoria que llegan a Σ^{s} se deslizan hasta el punto tangente visible T y tiende para el nodo estable (x^{∗}, y^{∗}), en cuanto las trayectorias que no llegan al segmento deslizante también tiende para el nodo estable. Por tanto, existe un nodo globalmente asintóticamente estable (x^{∗}, y^{∗}). Cuando P = x^{∗}, (x^{∗}, y^{∗}) colisiona con T formando un punto globalmente asintóticamente estable como lo observado en la Figura 8(b). Cuando P > x^{∗}, el nodo estable (x^{∗}, y^{∗}) no esta definido en el Sistema de Filippov (5) y existe un pseudo-equilibrio P N globalmente asintóticamente estable cercano al punto tangente invisible T como lo observado en la Figura 9(c).

Por otro lado, el sistema de Filippov (5) presenta una bifurcación Grazing cuando el ciclo límite Γ_{X} es tangente a T ∈ Σ^{s} en P = x_{L_{X}}. Por ejemplo, al considerar los valores para cada parámetro dados en la Tabla 1, si E = 0.2 y P ≈ 3.3, el ciclo límite estable Γ_{X} se convierte en un ciclo tangente a Σ como se observa en la Figura 6(b). En este caso, el ciclo tangente es atractor con una trayectoria de entrada deslizante estable.

Cuando P < x_{L_{X}}, el ciclo límite estable Γ_{X} coexiste con el punto tangente visible T como se observa en la Figura 6 (a). El ciclo tangente estable a Σ se convierte en un ciclo para x_{L_{X}} < P < x^{∗} como se observa en la Figura 7(a). Como el ciclo límite Γ_{X} no es conocido de forma explícita, tampoco es posible determinar una función P (E) que explique la curva P = x_{L_{X}}.

Si (x^{∗}, y^{∗}) es un foco inestable, entonces una bifurcación de Hopf se presenta al colisionar el ciclo límite Γ_{X} con el foco inestable, donde el ciclo se encoge y se transforma en un foco estable cuando los parámetros varían. La bifurcación de Hopf [8], ocurre cuando los valores propios de la matriz DX(x^{∗}, y^{∗}) tiene parte real nula, esto es, trDX(x^{∗}, y^{∗}) = 0, equivalente

E = E_{H}(P) = \frac{ac(b − K) + d(b + K)}{q(b + K)}

y la curva P = x_{L_{2}} separa las regiones 1 y 2 dadas en la Figura 10. Si (x^{∗}, y^{∗}) es un foco inestable, esto es, x^{∗} < \frac{K−b}{2} y △ < 0, entonces E_{H} < E y las regiones 1 y 2 presentan un foco inestable y un único ciclo, ciclo tangente, o ciclo límite globalmente estable cuando P < x_{L_{2}}, P = x_{L_{2}} o x_{L_{2}} < P, respectivamente. Si E_{H} > E y △ < 0, entonces x^{∗} > \frac{K−b}{2} y (x^{∗}, y^{∗}) es un foco estable.

A diferencia de los resultados de Yang et al. [12], cuando △ = 0, el foco estable (x^{∗}, y^{∗}) se transforma en un nodo estable que al colisionar con el punto silla (K, 0), del campo vectorial X, sobre el eje x es una bifurcación transcrítica determinada por x^{∗} = K, esto es,

E = E_{T C}(P) = \frac{acK − d(b + K)}{q(b + K)}.

De igual forma, si (x^{∗}, y^{∗}) colisiona con Σ, esto es, P = x^{∗}, entonces, 5

E = E_{F −N}(P) = \frac{acP − d(b + P)}{q(b + P)}

5 donde E_{F −N}(K) = E_{T C}(K) y se intercepta con la curva E_{H} cuando P = \frac{b(ac(b−K)+2d(b+K))}{2(acK−d(b+K))}, dividiendo la región 3 en las regiones 4 y 5.

Rev.Mate.Teor.Aplic. (ISSN print: 1409-2433; online: 2215-3373) Vol. 28(2): 183–208, Jul–Dec 2021

(a)  (b)

(c)

Figura 8: Retrato de fase do Sistema Planar de Filippov (5) cuando (x^{∗}, y^{∗}) es un foco estable, E = 0.22 y (a) P = 10 y (x^{∗}, y^{∗}) estable, (b) P ≈ 25.405 y (x^{∗}, y^{∗}) estable, (c) P = 30 y P N estable.

bifurcaciones en modelo gause depredador-presa con ... 205

(a)  (b)

(c)

Figura 9: Retrato de fase del Sistema de Filippov (5) cuando (x^{∗}, y^{∗}) es un nodo estable, E = 0.24 y (a) P = 20 y (x^{∗}, y^{∗}) estable, (b) P ≈ 33.793 y (x^{∗}, y^{∗}) estable, (c) P = 45 y P N estable.

(a) Plano (E, P)

(b) Región 1  (c) Región 2  (d) Región 3

(e) Región 4  (f) Región 5  (g) Región 6

Figura 10: Curvas de bifurcaciones para el Sistema de Filippov (5) en el plano (E, P) con parámetros mostrados en la Tabla 1. En la región 1 existe un ciclo; la región 2 un ciclo límite en X; la región 3 un pseudo-equilibrio estable; región 4 un foco estable; región 5 un nodo estable; región 6 una silla.

Al verificar si se han detectado todas las bifurcaciones para el sistema de Filippov (5), resta analizar si existe una colisión de T con (x^{∗}, y^{∗}) sobre el eje x y una colisión de un nodo inestable con Σ. Para el primer caso implica que P = K, lo cual, bajo las hipótesis del sistema de Filippov (5), no es probable. Para el segundo caso, si x^{∗} < \frac{K−b}{2}, esto es ac(b − K) + (b + K)(d + qE) < 0, y △ > 0, el sistema de Filippov (5) no podría presentar un nodo inestable.

En la Figura 10, y a diferencia del trabajo de Yang et al. [12], se muestra el diagrama de bifurcaciones en el sistema de Filippov (5) en el plano (E, P) y con valores en los parámetros mostrados en la Tabla 1, donde E = E_{H} ≈ 0.2102, E = E_{F −N} = \frac{5(P −4)}{16(P +5)} y E = E_{T C} ≈ 0.2613. De igual forma, los retratos de fase para el sistema de Filippov (5) son mostrados en diferentes regiones de la Figura 10.

Agradecimientos

El autor agradece tanto a los revisores como a los editores por su esfuerzos para la publicación del presente articulo.

Financiación

El autor agradece a la CAPES, a la Universidad Carlos III de Madrid y al CSIC por las becas de posgrado asignada.

Referencias

[1] V.I. Arnold, Ordinary Differential Equations, The MIT Press, translated from Russian and edited by Richard A. Silverman, Cambridge MA, London, 1973. https://mitpress.mit.edu/books/ordinary-differential-equations

[2] C.A. Buzzi , T. de Carvalho, P.R. da Silva, Closed poly-trajectories and Poincaré index of non-smooth vector fields on the plane, Journal of Dynamical and Control Systems, 19(2013), no. 2, 173–193. Doi: 10.1007/s10883-013-9169-4

[3] L. Edelstein-Keshet, Mathematical Models in Biology, Society for Industrial and Applied Mathematics, Philadelphia PA, 2005. Doi: 10.1137/1.9780898719147 [4] A.F. Filippov, Differential Equations with Discontinuous Righthand Sides,

Mathematics and Its Applications 18 (Soviet Series), Kluwer Academic Publishers Group, Dordrecht, 1988. Doi: 10.1007/978-94-015-7793-9

[5] H.I. Freedman, Deterministic Mathematical Models in Population Ecology. Marcel Dekker, New York, 1980. Doi: 10.2307/3556198

[6] M. Guardia, T.M. Seara, M.A. Teixeira, Generic bifurcations of low codimension of planar Filippov Systems, Journal of Differential Equations, 250(2010), no. 4, 1967–2023. Doi: 10.1016/j.jde.2010.11.016

[7] Y. Kuang, H.I. Freedman, Uniqueness of limit cycles in Gause-type models of predator-prey systems, Mathematical Biosciences, 88(1988), no. 1, 67–84. Doi: 10.1016/0025-5564(88)90049-1

[8] Y. Kuznetsov, Elements of Applied Bifurcation Theory, Applied Mathematical Sciences 112, Springer, New York, 1995. Doi: 10.1007/978-1-4757-3978-7

[9] Y. Kuznetsov, S. Rinaldi, A. Gragnani, One-parameter bifurcations in planar Filippov systems, International Journal of Bifurcation and Chaos, 13(2003), no. 8, 2157–2188. Doi: 10.1142/S0218127403007874

[10] J.D. Murray, Mathematical Biology: I. An Introduction, 3rd Edition, Springer, New York, 2002. Doi: 10.1007/b98868

[11] J. Sotomayor, Lições de equações diferenciais ordinárias, Projeto Euclides, Instituto de Matemática Pura e Aplicada, Rio de Janeiro, 11(1979).

[12] J. Yang, S. Tang, R.A. Cheke, Global stability and sliding bifurcations of a non-smooth Gause predator-prey system, Applied Mathematics and Computation, 224(2013), no. 1, 9–20. Doi: 10.1016/j.amc.2013.08.024