**1. Introduction**

The logistic equation with delay is the name given to the equation

$$
\dot{u} = \lambda [1 - u(t - T)]u. \tag{1}
$$

Being a natural extension of the classical logistic equation, it arises in many applied problems, but most of all it is used in problems of mathematical ecology (see, for example [1–15]). The value *u* describes the dynamics of changes in the number of specimens in the population or density of an isolated biological population living in a homogeneous environment. Therefore, it makes sense to consider only non-negative solutions of Equation (1). We note that the solution of Equation (1) with a non-negative initial function *<sup>ϕ</sup>*(*s*) ∈ *<sup>C</sup>*[−*T*,0] specified at some value *t*<sup>0</sup> (i.e., *u*(*t*<sup>0</sup> + *s*) = *ϕ*(*s*)) remains non-negative for all *t* > *t*0. The coefficient *λ* is called the Malthusian coefficient of linear growth and the parameter *T* > 0 is called the delay time. It is associated with the age of animal units capable of procreation in the corresponding population.

The dynamic properties of the solutions of Equation (1) are well understood. Under the condition *<sup>λ</sup><sup>T</sup>* <sup>≤</sup> *<sup>π</sup>* <sup>2</sup> , the equilibrium state *<sup>u</sup>*<sup>0</sup> <sup>≡</sup> 1 is asymptotically stable, while for *<sup>λ</sup><sup>T</sup>* <sup>&</sup>gt; *<sup>π</sup>* 2 it is unstable and the periodic solution *<sup>U</sup>*0(*t*, *<sup>λ</sup>*) is stable. Under the condition 0 <sup>&</sup>lt; *<sup>λ</sup><sup>T</sup>* <sup>−</sup> *<sup>π</sup>* <sup>2</sup> 1, we obtain its asymptotics on the Andronov–Hopf bifurcation theory application basis (see, for example [16–18]). It is shown in [7] that for sufficiently large *λ* the periodic solution *U*0(*t*, *λ*) has a pronounced relaxation structure. There is a single wavelet of this function on the segment of the period *T*0(*λ*) length. The duration of the wavelet is close to

**Citation:** Kashchenko, A.; Kashchenko, S. Relaxation Oscillations in the Logistic Equation with Delay and Modified Nonlinearity. *Mathematics* **2023**, *11*, 1699. https://doi.org/10.3390/ math11071699

Academic Editor: Carmen Chicone

Received: 28 February 2023 Revised: 27 March 2023 Accepted: 31 March 2023 Published: 2 April 2023

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

the value of *T* and the period *T*0(*λ*) is equal to (*λT*)−<sup>1</sup> exp(*λT*(1 + *o*(1))). It is important to note the asymptotic equalities.

$$\max\_{t} \mathcal{U}\_0(t, \lambda) \quad = \quad \exp(\lambda T(1 + o(1))), \tag{2}$$

$$\min\_{t} \mathcal{U}\_0(t, \lambda) \quad = \ \exp(-\exp(\lambda T)(1 + o(1))).\tag{3}$$

We also note results of [19,20], in which important results on the dynamics of logistic equations with delay and diffusion are obtained, [21], where a randomized nonautonomous logistic equation is discussed, and Ref. [22], where interesting results were obtained for an equation in which the delay itself is a function of *u*.

Due to the applied significance, it is of interest to study equations with various types of nonlinearity. By virtue of biological meaning, solutions must remain positive. Similar to this restriction, it is natural to require solutions to be bounded not only from below, but also from above. In addition, unlimited asymptotic growth, as in the case of Equation (1), is not always convenient, for example, for describing the dynamics of changes in the number of species.

In this paper we study the dynamics of the new model—the logistic equation with delay and modified nonlinearity

$$
\dot{u} = \lambda \left[ 1 - u(t - T) \right] u(A - u), \quad A > 1. \tag{4}
$$

The main difference between the solutions of Equation (4) and the solutions of Equation (1) is that under the condition 0 ≤ *ϕ*(*s*) ≤ *A* the inequalities

$$0 \le \mu(t, \boldsymbol{\varrho}) \le A$$

hold for the solution *u*(*t*, *ϕ*) of Equation (4) with the initial condition *u*(*t*<sup>0</sup> + *s*) = *ϕ*(*s*) as *t* > *t*0. Thus, the use of the factor (*A* − *u*) on the right-hand side of Equation (4) leads to the restriction *u*(*t*, *ϕ*) ≤ *A* for the solutions of Equation (4). Therefore, the solutions of Equation (4) (in contrast to the solutions of Equation (1)) do not take asymptotically large values for large values of *λ*.

Model (4) is quite complex. At the moment, there are no analytical methods that would allow for absolutely any values of the parameters *λ*, *T*, and *A* to obtain results on the qualitative and quantitative behavior of all solutions of Equation (4) on the entire positive semiaxis *t* ∈ [0, +∞). In addition, numerical methods are also not very effective in the study of such models, since, firstly, it is impossible to enumerate all possible initial conditions in order to draw a qualitative conclusion about the dynamics of the model, and secondly, for some values of the model parameters, the solutions are very close will approach zero (see Sections 3 and 4), so a relatively small error in the calculations can lead to a completely wrong result.

Therefore, in the paper, the research is carried out by analytical (asymptotic) methods. At present, the most effective among them are the methods of the theory of bifurcations in the study in the neighbourhood of the equilibrium state and the methods of a large parameter (for *λ* 1) in the study of nonlocal processes. Large parameter methods are based on the use of singular perturbation theory, but even in cases where the *λ* parameter (Malthusian coefficient in model Equation (4)) is not large enough, nevertheless, conclusions can be drawn about the trends in fluctuations with an increase in this parameter. In addition, we point out that the relaxation nature of the oscillations obtained for model Equation (4) in the case of sufficiently large values of parameter *λ* corresponds to ideas about the dynamic behavior of population density dynamics in mathematical ecology [4,11].

The paper is organized as follows. In Section 2, the conditions for the parameter *λ* are found using the methods of bifurcation analysis. Under them, the Andronov–Hopf bifurcation occurs in Equation (4) and a cycle arises from the equilibrium state.

The main content of Section 3 and all subsequent sections is devoted to the study of the solutions of Equation (4) on the interval *t* ∈ [0, +∞) provided that the parameter *λ* is sufficiently large:

$$
\lambda \gg 1.\tag{5}
$$

In this case, the equation under consideration is singularly perturbed. After dividing by the left-hand and right-hand parts of Equation (4) by *λ*, we obtain an equation with a small parameter at the derivative:

$$
\lambda^{-1}\dot{u} = [1 - u(t - T)]u(A - u). \tag{6}
$$

However, the reduced equation

$$0 = [1 - \mu(t - T)]\mu(A - \mu)$$

does not describe the behavior of the Equation (6) solutions (except the simplest equilibrium states) under the condition (5).

It is important to mention that from a computational point of view, this problem is difficult, since the relaxation solution approaches *A* and 0 very closely, so even a small error in the calculations will take us out of the class of solutions under consideration.

Therefore, a special analytical research method was developed. First, we briefly describe it for the simplest case, for studying slowly oscillating solutions. Recall that slowly oscillating solutions are those solutions for which the time distance between adjacent roots of the equation *u*(*t*) = 1 is greater than the delay time. Consider Equation (6). Denote by *<sup>S</sup>* ⊂ *<sup>C</sup>*[−*T*,0] the set of all such functions *<sup>ϕ</sup>*(*s*) ∈ *<sup>S</sup>* (*<sup>s</sup>* ∈ [−*T*, 0]) that satisfy the conditions (see Figure 1).

$$
\varphi(0) = 1; \quad 0 \le \varphi(s) \le \exp\left(\frac{\lambda}{2}s\right).
$$

$$
\sum\_{s=0}^{\infty} \frac{1}{s} \le \frac{1}{2}
$$

**Figure 1.** The set *S*.

Denote by *u*(*t*, *ϕ*) a solution of (6) with initial condition *ϕ*(*s*):*u*(*s*, *ϕ*) = *ϕ*(*s*) (*s* ∈ [−*T*, 0]). The method of steps will be used: for *t* ∈ [0, *T*], taking into account the initial function *ϕ*(*s*), to determine *u*(*t*, *ϕ*) we obtain an ordinary differential equation of the first order. It is easy to find the asymptotics of this solution as *λ* → ∞ on the indicated interval. Then we consider *u*(*t*, *ϕ*) on the interval [*T*, 2*T*], on which we also obtain a first-order ordinary differential equation, and find the asymptotics of the solutions. After that, we will carry out the same actions on the segment [2*T*, 3*T*], then on the segment [3*T*, 4*T*], and so on.

In Section 3.1 we study the asymptotic behavior at *λ* → ∞ of all solutions *u*(*t*, *ϕ*) with initial conditions from *S*. In particular, the asymptotics of the first two positive roots *t*1(*ϕ*) and *t*2(*ϕ*) of the equation *u*(*t*, *ϕ*) = 1 will be found. The main conclusion is that after the time interval *t*2(*ϕ*) the solution *u*(*t*, *ϕ*) will again fall on the set *S*:*u*(*s* + *t*2(*ϕ*), *ϕ*) ∈ *S*. This gives reason to introduce the Π operator:

$$
\Pi(\boldsymbol{\varrho}(\mathbf{s})) = \boldsymbol{\iota}(t\_2(\boldsymbol{\varrho}) + \boldsymbol{s}, \boldsymbol{\varrho}).
$$

and justify the inclusion Π*S* ⊂ *S*. From this and the well-known results of functional analysis [23], we conclude that *S* has a fixed point *ϕ*0(*s*) of the operator Π: Π(*ϕ*0(*s*)) = *ϕ*0(*s*). Then the function *u*0(*t*) = *u*(*t*, *ϕ*0) is periodic with period *t*2(*ϕ*0). In Section 3.1, asymptotic formulas for *u*0(*t*) will be given.

In Section 3.2 we investigate the asymptotic behavior of rapidly oscillating solutions. As a set of initial functions, we consider the set *S*(*τ*1, *τ*2) depending on two parameters *τ*<sup>1</sup> and *τ*2. It is shown schematically in Figure 3. In Section 3.2, a rigorous description of this set is given and the asymptotic behavior of all solutions (6) with initial conditions from *S*(*τ*1, *τ*2) is studied. Again, by *t*1(*ϕ*) and *t*2(*ϕ*) we denote the first and second positive roots of the equation *u*(*t*, *ϕ*) = 1. We again introduce the operator Π: Π(*ϕ*(*s*)) = *u*(*t*2(*ϕ*) + *s*, *ϕ*) and it will be shown that up to *o*(1) the inclusion Π(*ϕ*(*s*)) ∈ *S*(*τ*¯1, *τ*¯2), where *τ*¯1,2 are explicitly expressed in terms of *τ*1,2. This gives grounds to take the function *u*(*s* + *t*2(*ϕ*), *ϕ*) as the initial one and pass to the iteration process. It is possible to find a fixed point (*τ*<sup>0</sup> <sup>1</sup> , *<sup>τ</sup>*<sup>0</sup> <sup>2</sup> ) and a function *<sup>ϕ</sup>*0(*s*) ∈ *<sup>S</sup>*(*τ*<sup>0</sup> <sup>1</sup> , *<sup>τ</sup>*<sup>0</sup> <sup>2</sup> ) such that the solution *u*0(*t*, *ϕ*0) will be periodic with period *t*2(*ϕ*0). It is important to note that the rapidly oscillating periodic solutions found in this way are unstable.

The results of Section 4 are more interesting. We consider equations with two delays

$$
\dot{u} = \lambda \left[ 1 - \alpha u(t - T) - (1 - \alpha) u(t - h) \right] u(A - u), \tag{8}
$$

where *α* ∈ (0, 1) and *h* < *T*.

It is commonly supposed that the use of two (or more) delays allows one to take into account the influence of the population age structure on the dynamics of population level changes [24–27].

One succeeded to clearly describe the nonlocal dynamic properties of Equation (8) under the condition (5). In contrast to the results of the Section 3, the relaxation step solutions for (8) can be stable.

We note that for the simpler logistic equation with two delays

$$\dot{u} = \lambda [1 - \alpha u(t - T) - (1 - \alpha)u(t - h)]u \quad (u \in (0, 1), \ h < T) \tag{9}$$

the results were obtained in Refs. [7] as *λ* 1. It follows from them that for sufficiently large *λ* in (9) there is an orbitally stable periodic solution which is similar in appearance to the solution of Equation (1) with the same initial conditions, and the asymptotic representation (2) holds for it.

## **2. Andronov–Hopf Bifurcations in Equations with One Delay**

We shall fix the parameters *λ*0, *A*<sup>0</sup> and *T*<sup>0</sup> so that the following equality is fulfilled:

$$
\lambda\_0 (A\_0 - 1) T\_0 = \frac{\pi}{2}.\tag{10}
$$

**Lemma 1.** *Suppose that in Equation* (4)

$$
\lambda (A - 1)T < \frac{\pi}{2}.\tag{11}
$$

*Then the equilibrium state u*<sup>0</sup> ≡ 1 *is asymptotically stable.*

**Proof.** To prove this statement we linearize (4) on *u*<sup>0</sup> and consider the characteristic quasipolynomial

$$
\mu = -\lambda (A - 1) \exp(-\mu T). \tag{12}
$$

Under the condition (11) all its roots have negative real parts. From this, and from known results (see Ref. [5]) on the stability in the first approximation, the proof of Lemma 1 follows.

Next we introduce a small positive parameter *ε*:

$$0 < \varepsilon \ll 1,\tag{13}$$

and assume that in Equation (4), for certain constant *λ*1, *A*<sup>1</sup> and *T*1, the following relations are fulfilled:

$$
\lambda = \lambda\_0 + \varepsilon \lambda\_1, \quad A = A\_0 + \varepsilon A\_1, \quad T = T\_0 + \varepsilon T\_1. \tag{14}
$$

We set

$$b = \left(1 + \frac{\pi^2}{4}\right)^{-1} \left[\left(\frac{\pi}{2} + i\right)(\lambda\_1(A\_0 - 1) + \lambda\_0 A\_1) + \lambda\_0^2 (A\_0 - 1)^2 T\_1 \left(1 - i\frac{\pi}{2}\right)\right],\tag{15}$$

$$d = -\lambda\_0 \left( 1 + \frac{\pi^2}{4} \right)^{-1} \left[ \frac{\pi}{2} + i + \frac{3(A\_0 - 2)^2}{5(A\_0 - 1)} \left( \frac{\pi}{2} - 1 + i(\frac{\pi}{2} + 1) \right) \right],\tag{16}$$

$$
\omega = \pi (2T\_0)^{-1} \,\text{\textquotedblleft}\tag{17}
$$

$$\begin{aligned} \xi\_0(\tau) &= \xi\_0 \exp(i\psi\tau), \quad \psi = \text{Im } b + \xi\_0^2 \text{Im } d, \\ \xi\_0 &= \left[ \left( \frac{\pi}{2} (\lambda\_1 (A\_0 - 1) + \lambda\_0 A\_1) + \lambda\_0^2 (A\_0 - 1)^2 T\_1 \right) \lambda\_0^{-1} \left( \frac{\pi}{2} + \frac{3(A\_0 - 2)^2}{5(A\_0 - 1)} \left( \frac{\pi}{2} - 1 \right) \right)^{-1} \right]^{\frac{1}{2}}. \end{aligned} \tag{18}$$

**Theorem 1.** *Let conditions* (13)*–*(15)*,* (17)*, and* (18) *hold.*

*1. Let* Re *b* < 0*. Then all solutions of Equation* (4) *with initial conditions from a sufficiently small (and independent of ε) neighbourhood of the equilibrium state u* ≡ 1 *tend to* 1 *as t* → ∞*. 2. Let*

$$\text{Re}\,b > 0.\tag{19}$$

*Then for all sufficiently small ε Equation* (4) *has in the neighbourhood of the unit equilibrium state a stable periodic solution u*0(*t*,*ε*) *for which the following asymptotic equality is fulfilled:*

$$u\_0(t, \varepsilon) = 1 + \varepsilon^{1/2} [\xi\_0(\tau) \exp(i\omega t) + \xi\_0(\tau) \exp(-i\omega t)] + O(\varepsilon),$$

*where τ* = *εt.*

**Proof.** Under the condition (10) the quasi-polynomial (12) has a pair of purely imaginary roots *μ*1,2 = ±*iω*, where the value of *ω* is given in (17), and all its other roots have negative real parts. Then, under the condition (14), and for all sufficiently small *ε*, Equation (4) in the neighbourhood of *u* ≡ 1 has (see Refs. [16–18,28]) a two-dimensional stable local invariant integral manifold, on which this equation can be written, to within higher-order terms, in the form of a complex scalar equation—a normal form:

$$\frac{d\tilde{\xi}}{d\tau} = B\tilde{\xi}(\tau) + D\tilde{\xi}(\tau)|\tilde{\xi}(\tau)|^2,\tag{20}$$

where values *B* and *D* are to be determined.

In the case Re *B* = 0 and Re *D* = 0 normal form completely determines the behavior of all solutions in a neighbourhood of *u* ≡ 1 (see Refs. [16–18,28]). The solutions of (20) are connected with the solutions of (4) by the relation

$$\begin{split} u(t, \varepsilon) &= 1 + \varepsilon^{1/2} [\tilde{\xi}(\tau) \exp(i\omega t) + \tilde{\xi}(\tau) \exp(-i\omega t)] + \\ &+ \varepsilon u\_2(t, \tau) + \varepsilon^{3/2} u\_3(t, \tau) + \dots \end{split} \tag{21}$$

Here, *τ* = *εt*, the *uj*(*t*, *τ*) are periodic in *t* with period 2*π*/*ω*. Substituting (21) into (4) and collecting coefficients of equal powers of *ε*, we successively find all the elements appearing there. At *ε*1/2 we obtain true identity, and collecting coefficients of the first power of *ε* and taking into account (10), (12), and (17), we find equation for determining the function *u*2(*t*, *τ*):

$$\frac{\partial \mu\_2}{\partial t} = -\lambda\_0 (A\_0 - 1)\mu\_2 (t - T\_0, \tau) + i\lambda\_0 (A\_0 - 2)\tilde{\xi}^2(\tau) \exp(2i\omega t) - i\lambda\_0 (A\_0 - 2)\tilde{\xi}^2(\tau) \exp(-2i\omega t).$$

Thus, taking into account that *u*2(*t*, *τ*) is periodic in *t* with period 2*π*/*ω*, we find that

$$u\_2(t, \tau) = \frac{(A\_0 - 2)(2 - i)}{\mathfrak{F}(A\_0 - 1)} \tilde{\xi}^2(\tau) \exp(2i\omega t) + \frac{(A\_0 - 2)(2 + i)}{\mathfrak{F}(A\_0 - 1)} \tilde{\xi}^2(\tau) \exp(-2i\omega t).$$

Collecting coefficients of *ε*3/2 we obtain an equation for *u*3(*t*, *τ*):

$$\frac{\partial u\_3}{\partial t} = -\lambda \mathbf{o}(A\mathbf{o} - \mathbf{1}) u\_3 \mathbf{t} (t - T\_0, \tau) + B\_1 \exp(i\omega t) + \bar{B}\_1 \exp(-i\omega t) + B\_3 \exp(3i\omega t) + \bar{B}\_3 \exp(-3i\omega t), \tag{22}$$

where *B*<sup>1</sup> and *B*<sup>3</sup> are complex values.

From the condition that Equation (22) be solvable in the class of functions periodic in *t* with period 2*π*/*ω* (this condition is *B*<sup>1</sup> = 0) we arrive at Equation (20) for the unknown amplitude *ξ*(*τ*), in which *B* = *b* and *D* = *d* (coefficient b is defined in (15) and the coefficient *d* is defined in (16)).

Multiplying both sides of Equation (20) by ¯ *ξ*(*τ*) and taking into account equalities *B* = *b* and *D* = *d* we obtain a scalar real ordinary differential equation

$$\frac{1}{2}\dot{\rho}(\tau) = \text{Re }b\rho(\tau) + \text{Re }d\rho^2(\tau),\tag{23}$$

where *ρ*(*τ*) = |*ξ*(*τ*)| 2.

We note that Re *d* < 0.

That is why, if Re *b* < 0, all solutions of (23) tend to zero as *t* → +∞. Thus, all solutions of (4) from the neighbourhood of *u* ≡ 1 tend to 1 as *t* → +∞. This completes the proof of the first part of the Theorem.

Under the condition (19) Equation (20) with *B* = *b* and *D* = *d* has a stable periodic solution *ξ*0(*τ*)=*ξ*<sup>0</sup> exp(*iψτ*), where the values *ξ*<sup>0</sup> and *ψ* are given in (18). Taking into account here the asymptotic equality (21) we have completed the proof of the Theorem.

It is interesting to note that the bifurcation effect considered here can be realized for fixed *λ* = *λ*<sup>0</sup> and *T* = *T*0, i.e., for *λ*<sup>1</sup> = *T*<sup>1</sup> = 0, and only upon variation of the parameter *A*.
