EUROPEAN JOURNAL OF PURE AND APPLIED MATHEMATICS Vol. 17, No. 2, 2024, 870-904 ISSN 1307-5543 – ejpam.com Published by New York Business Global Modeling and multi-objective optimal control of the dynamics of counterterrorism in the Sahel region in Africa Mathieu Romaric Pooda1,∗, Yacouba Simpore2,3, Oumar Traore1,2 1 Laboratoire de Sciences et Technologie, Université Thomas SANKARA, 12 BP 417 Ouagadougou 12, Burkina Faso. 2 Laboratoire d’Analyse Mathématiques et d’Informatique, Université Joseph KI-ZERBO, 03 BP 7021 Ouagadougou 03, Burkina Faso. 3 DeustoTech Fundación Deusto Avda Universidades, 24, 48007, Bilbao, Basque Country, Spain. Abstract. Terrorist activity in the Sahel region has been on the increase for almost a decade. Groups advocating extremist ideologies with a political or religious base are carrying out attacks against states with the aim of imposing a totalitarian ideology, sometimes against a backdrop of political crisis and famine. In this state of crisis and sometimes communal conflict, it is questionable whether ideological terrorism can really be eradicated in the Sahel. In this study, we develop a model of the dynamics of ideological terrorism based on the terrorism situation in the Sahel. In particular, this model incorporates popular resistance to terrorism through the class of volunteers for the defense of the homeland, but we also take into account the indoctrination of certain sections of the population vulnerable to fanatical ideology. We estimate R0 , the number of elementary replications of extremist behavior, which enables us to predict the evolution of extremism. We also identify four thresholds R1, R2, R3 and R4 of sufficient conditions for the eradication of ideological terrorism and brigandage. A multi-objective optimal control of a counter-terrorism strategy is also presented. Finally, we perform a numerical simulation of the analysis and control results to test our hypotheses. 2020 Mathematics Subject Classifications: 49K15, 93B05, 93C15, 93D23 Key Words and Phrases: Ideological terrorism, fanatical behavior, local and global asymptotic stability, global threshold, optimal control 1. Introduction Africa’s Sahel region faces major security challenges, not least terrorism, which poses a considerable threat to regional stability. Several terrorist groups operate in the region, ∗Corresponding author. DOI: https://doi.org/10.29020/nybg.ejpam.v17i2.5122 Email addresses: math7roma8@gmail.com (M. R. Pooda), simplesaint@gmail.com (Y. Simpore), oumar.traore@uts.bf (O. Traore) https://www.ejpam.com 870 © 2024 EJPAM All rights reserved. M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 871 including Al-Qaeda in the Islamic Maghreb (AQIM), the Group for Support of Islam and Muslims (GSIM), Boko Haram and the Islamic State in the Greater Sahara (ISGS). Their motivations vary, ranging from radical ideological claims to local rivalries for control of resources and territory, to funding opportunities through illicit trafficking. These groups employ a variety of tactics, including bomb attacks, ambushes and kidnappings, often tar- geting security forces, critical infrastructure and civilian populations. Fragile institutions, poverty, high unemployment and ethnic tensions contribute to the region’s vulnerability. International responses include efforts to support Sahel governments through training, ca- pacity building and security cooperation, but a holistic, long-term approach is needed to effectively tackle the root causes of instability and insecurity. Let us note that, evolving in time and space, terrorism is not an easy concept to de- fine. International law is unable to give it a clear definition. When we speak of terrorism, we generally mean terrorist methods. Oscillating and unpredictable, the means used by terrorists have varied greatly over time, as have their objectives. The threat seems to come from nowhere to strike anywhere. This raises the question of the purpose of terrorist organizations, which operate in the shadows, making them difficult to catch. There are four main types of terrorism: individual terrorism, caused by rebels, anarchists or nihilists (who accept moral freedom); organized terrorism, advocated by groups defending different ideologies; state terrorism and cyber terrorism. These various forms of terrorism are moti- vated by vengeful hatred (hatred based on a person’s determination to avenge the abuses for which their enemies are responsible), deterrence (so that the terrorized population can exert pressure on their government), propaganda (to strike an emotional chord), and provocation (to push a government to react excessively). Thus, terrorist methods have evolved, and although no definitive definition is currently established, we can attempt to understand what terrorist methods and motivations entail. Additionally, terrorist groups are generally driven by significant ideologies that form the basis of their recruitment. Terrorism represents a significant threat in the contemporary world, garnering seri- ous attention from numerous countries. Among those affected are Burkina Faso, Mali, Niger, Nigeria, and the United States of America, among others. Consequently, these nations spare no effort in implementing measures to mitigate the risk of attacks. They tailor local and global security initiatives to their specific circumstances in order to protect themselves effectively. Despite serious attempts to establish frameworks for analyzing the dynamics of human behavior [14], [6], [26], [18], there is still no comprehensive mathemat- ical framework or approach for systematically studying human behavior. Similarly, there is no comprehensive mathematical framework or approach for the systematic study of the transmission of ideas [9]. Nevertheless, a wide range of mathematical modeling problems related to evolution, contagion or propagation phenomena have been explored. These ef- forts encompass programs using paradigms rooted in evolutionary biology, with notable contributions from sources such as [6], [26], [18]. However, there is no longer any doubt that mathematical modeling can play a very important role in understanding and solving these critical problems. We can cite work on migration and crowd behavior [7],[28],[20], M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 872 work on crime [12], [13], [24], [17], [15], [16], [27], [22], gang membership [30], [2], [8]. We can also cite works on the dynamics of war [3], [11], [10] and many others. In this work, we study the dynamics of the transmission of extreme behaviors as a kind of contact process in epidemiology combined with local and global security initiatives, as well as the mode of recruitment of terrorists. In particular, we have designed a model of the dynamics of ideological terrorism in a population developing anti-terrorism and anti- brigandage initiatives. This model, which describes the dynamics of ideological terrorism or fanatical insurgency, is inspired on the one hand by the terrorism situation in the Sahel countries, and is characterized by a sequence of thresholds that provide quantitative and qualitative information on the evolution of the permanent insecurity situation. A first local threshold makes it possible to evaluate the number of elementary reproductions of the behavior extremismR0, well known in epidemiology, and to predict the evolution of violent extremism. In addition to this first local threshold, we identify four global thresholds R1, R2, R3, and R4 that essentially give us sufficient conditions for the eradication or extinction of the radical core subpopulation of terrorism, fanatical ideology and brigandage respectively. We conclude this study with a numerical simulation and an optimal control analysis. The specifics of the model are described in the next section. 2. Model formulation To aid comprehension of this model, we introduce the following definitions and assump- tions: the host population is categorized into two subpopulations. One is a non-radical subpopulation denoted as G(t), while the other is a radical subpopulation characterized by violent extremism and fanaticism, defined as D(t). The non-radical subpopulation comprises individuals not indoctrinated into extreme ideologies and not involved in criminal activities, enjoying their freedom. This subset G(t) is further divided into four classes: civilians not indoctrinated in extreme ideology and non-combatants denoted as C(t), the defense and security forces referred to as A(t), a highly organized group formed as part of a state’s global security initiative, composed of individuals trained in warfare techniques, combat, and counter-terrorism activities, along with elites. Additionally, there are volunteers for homeland defense denoted as V (t), part of a local initiative for vigilance and defense against terrorism and criminal activities, and those who have been removed from the ranks of the defense and security forces designated as R(t). On the other hand, the radical subpopulation, mainly comprising individuals indoc- trinated with extremist ideologies, is divided into six classes, with four being hierarchical. This hierarchy is determined by the level of individual commitment to the ideology, with fanatical individuals assumed to be the most effective in propagating the ideology to vul- nerable members of the population, making recruitment a crucial factor. Recruitment in this context is modeled based on previous works [9], [21], [23], [31]. The classes in the M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 873 radical subpopulation include: the naive or vulnerable class S(t), consisting of individuals frequently in contact with members of the radical core but not yet converted; the semi- fanatical class E(t), comprising individuals newly converted to the ideology or not fully committed; the fanatical class F (t), consisting of individuals fully embracing the extreme ideology; and the terrorist class T (t), encompassing all armed individuals participating in acts of terrorism targeting symbols of the state. In addition to these hierarchical classes within the radical subpopulation, there is the brigand class B(t), comprising individuals not indoctrinated with radical ideologies but engaged in violent criminal activities such as robbery and looting, often operating in gangs. Finally, there is the prisoner class P (t). The dynamics of ideological terrorism or fanatical insurgency are described in the diagram below. Figure 1: Diagram of the dynamics of ideological terrorism or fanatic insurgency Where Λ represents the turnover constant of class C, D = S + E + F + T designates the radical core strongly dominated by indoctrination, I = A + V + B + T the combatant class, and N = C + R + A + V + S + E + F + B + T + P designates the total pop- ulation. The conversion rate from class C, V,A,B,R, and P to class S is β1 D N , where β1 = π1C + π2V + π3A+ π4B+ π5R+ π6P and πi, i ∈ {1, 2, 3, 4, 5, 6} represent the inten- sity of the recruitment force associated with D on C for i = 1, V for i = 2, A for i = 3, B for i = 4, R for i = 5, P for i = 6. This conversion or recruitment occurs during individual contacts between members of D and individuals not in class D. A contact has a broad meaning here as it includes distant contacts (phone calls, emails, etc.). This type of con- tact is distinguished from contacts in epidemiology. We define β2 E + F + T N , β3 F + T N , M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 874 β4 T N as the conversion rates from class S to class E, from class E to class F , and from class F to class T respectively, where β2, β3, and β4 designate the contact frequency or recruitment force intensity associated with D for classes S,E, F , and T respectively. Similarly, ω1 T I , ω2 T I , ω3 T R+ I , θ1 T P + I represent the recruitment rates into class T of individuals from class B,A,R, and P respectively. Additionally, ω4 B R+ I denotes the conversion rate from class R to class B, with ω4 representing the intensity of the conversion force in class B acting on individuals in class R. We designate by θ2 B P + I the conversion rate from class P to class B, ν3 B I the conversion rate from class A to class B, where ν3 indicates the intensity of the conversion force in class B acting on individuals in class A. The transfer rates in class P for individuals of classes F,B, and T under the influence of individuals in class A or V are τ1 A+ V F + I , τ2 A+ V I , τ3 A+ V I respectively. Furthermore, α2 B C + I is the rate of conversion from class C to class B, and α1 T +B C + I is the rate of conversion from class C to class V . These rates are proportional to the number of contacts per unit of time and the probability of success. Moreover, ν1, σ1, σ2 represent the recruitment rates of class B, V , and C respectively into class A, while ν2 denotes the attrition rate in class A. Additionally, γi, i ∈ {1, 2, 3, 4, 5, 6, 7, 8, 9}, represents the recovery rate or return to normal civilian life for individuals in classes S,E, F,A, P,R,B, T, V respectively. The assumption about the hierarchy of certain classes in the radical subpopulation implies that γ3 is extremely small, indicating that the av- erage durations of individuals in the fanatical and terrorist classes are very long ( 1 γ3 >> 0). We denote by δi, i ∈ {1, 2, 3, 4, 5}, the probability of dying following a fight of indi- viduals in class A for i = 1, of class V for i = 2, of class F for i = 3, of class B for i = 4, and of class T for i = 5. These probabilities are assumed to be proportional to the numbers of contacts between combatants and the intensity of the nuisance force that each combatant class exerts on its opponents. Hence, we have δ1 = ζ1 T +B I , δ2 = ζ2 T +B I , δ3 = ζ3 A+ V I , δ4 = ζ4 A+ V I , δ5 = ζ5 A+ V F + I , where ζi, i ∈ {1, 2, 3, 4, 5}, represents the intensity of the nuisance force of classes T and B on class A and V respectively for i = 1 and i = 2, and of A and V on classes F,B, and T respectively for i = 3, i = 4, and i = 5. Additionally, η represents the probability of dying in prison as a result of torture or any other form of maltreatment. Finally, we assume that all individuals in the total population have the same natural mortality rate µ. Therefore, the model is expressed as follows: M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 875 dC dt = Λ + γ1S + γ2E + γ3F + γ4A + γ5P + γ6R + γ7B + γ8T + γ9V − ( π1 D N + α1 T + B C + I + α2 B C + I + σ2 + µ ) C (1) dR dt = ν2A − ( π5 D N + ω3 T R + I + ω4 B R + I + γ6 + µ ) R (2) dA dt = σ1V + σ2C + ν1B − ( π3 D N + ν3 B I + ω2 T I + γ4 + ν2 + µ + ζ1 T + B I ) A (3) dV dt = α1C T + B C + I − ( π2 D N + γ9 + σ1 + µ + ζ2 T + B I ) V (4) dS dt = ( π1C + π2V + π3A + π4B + π5R + π6P ) D N − ( β2 E + F + T N + γ1 + µ ) S (5) dE dt = β2S E + F + T N − ( β3 F + T N + γ2 + µ ) E (6) dF dt = β3E F + T N − ( γ3 + τ1 A + V F + I + β4 T N + µ + ζ3 A + V F + I ) F (7) dB dt = α2C B C + I + ω4R B R + I + ν3A B I + θ2P B P + I − ( π4 D N + ω1 T I + τ2 A + V I + γ7 + ν1 + µ + ζ4 A + V I ) B (8) dT dt = β4F T N + ω1B T I + ω2A T I + ω3R T R + I + θ1P T P + I − ( τ3 A + V I + γ8 + µ + ζ5 A + V I ) T (9) dP dt = τ1F A + V F + I + τ2B A + V I + τ3T A + V I − ( π6 D N + θ1 T P + I + θ2 B P + I + γ5 + µ + η ) P (10) With non-negative initial conditions given by : C(0) > 0;S(0) ≥ 0;E(0) ≥ 0;F (0) ≥ 0;V (0) ≥ 0;A(0) > 0;R(0) ≥ 0;B(0) ≥ 0;P (0) ≥ 0;T (0) ≥ 0, N(0) ⩽ Λ µ (11) The parameters of the system (1)− (10) are assumed to be all non-negative. 3. Mathematical analysis of the model 3.1. Existence and uniqueness of solution The (1)−(10) model is described by a system of first order nonlinear differential equations. It is rewritten as follows: X ′(t) = f(X(t)) (12) where X(t) is a column vector of the number of individuals by class, and f : R10 → R10 is a fonction. X(t) =  C(t) R(t) A(t) V (t) S(t) E(t) F (t) B(t) T (t) P (t)  (13) M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 876 and f(x) =  Λ + γ1x5 + γ2x6 + γ3x7 + γ4x3 + γ5x10 + γ6x2 + γ7x8 + γ8x9 + γ9x4 − ( π1 x12 x13 + α1 x9 + x10 x1 + x11 + α2 x8 x1 + x11 + σ2 + µ ) x1 ν2x3 − ( π5 x12 x13 + ω3 x9 x2 + x11 + ω4 x8 x2 + x11 + γ6 + µ ) x2 σ1x4 + σ2x1 + ν1x8 − ( π3 x12 x13 + ν3 x8 x11 + ω2 x9 x11 + γ4 + ν2 + µ + ζ1 x9 + x8 x11 ) x3 α1x1 x9 + x8 x1 + x11 − ( π2 x12 x13 + γ9 + σ1 + µ + ζ2 x9 + x8 x11 ) x4 ( π1x1 + π2x4 + π3x3 + π4x8 + π5x2 + π6x10 ) x12 x13 − ( β2 x6 + x7 + x9 x13 + γ1 + µ ) x5 β2x5 x6 + x7 + x9 x13 − ( β3 x7 + x9 x13 + γ2 + µ ) x6 β3x6 x7 + x9 x13 − ( γ3 + τ1 x3 + x4 x7 + x11 + β4 x9 x13 + µ + ζ3 x3 + x4 x7 + x11 ) x7 α2x1 x8 x1 + x11 + ω4x2 x8 x2 + x11 + ν3x3 x8 x11 + θ2x10 x8 x10 + x11 − ( π4 x12 x13 + ω1 x9 x11 + τ2 x3 + x4 x11 + γ7 + ν1 + µ + ζ4 x3 + x4 x11 ) x8 β4x7 x9 x13 + ω1x8 x9 x11 + ω2x3 x9 x11 + ω3x2 x9 x2 + x11 + θ1x10 x9 x10 + x11 − ( τ3 x3 + x4 x11 + γ8 + µ + ζ5 x3 + x4 x11 ) x9 τ1x7 x3 + x4 x7 + x11 + τ2x8 x3 + x4 x11 + τ3x9 x3 + x4 x11 − ( π6 x12 x13 + θ1 x9 x10 + x11 + θ2 x8 x10 + x11 + γ5 + µ + η ) x10  (14) with x = (x1, x2, x3, x4, x5, x6, x7, x8, x9, x10) ∈ R10 and  x11 = x3 + x4 + x8 + x9 x12 = x5 + x6 + x7 + x9 x13 = x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 + x10 The function f is clearly locally lipschitzian with respect to x. We then deduce the existence and the uniqueness of the maximal solution to the Cauchy problem associated to the differential equation (1)− (10) related to the initial condition (11). 3.2. Positivity of the solutions For this model of the dynamics of ideological terrorism to be realistic, it is necessary to show that all state variables remain positive at all times. Proposition 1. (Positivity) The positive orthan R10 ≥0 is positively invariant for the system (1) − (10), and the initial condition (11) ensures the positivity of the solutions of the system (1)− (10) for any time t > 0. Proof: We use the barrier theorem [5]. Let us show that the set { C ≥ 0 } is positively invariant. Let x = ( C,R,A, V, S,E, F,B, T, P ) and consider L an application defined by L(x) = −C (15) M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 877 The application L thus defined is differentiable and we have: ∇L(x) = (−1, 0, 0, 0, 0, 0, 0, 0, 0, 0) ̸= 0R10 . (16) The vector field for { C = 0 } is given by X(x) =  Λ + γ1S + γ2E + γ3F + γ4A+ γ5P + γ6R+ γ7B + γ8T + γ9V ν2A− ( π5 D N + ω3 T R+ I + ω4 B R+ I + γ6 + µ ) R σ1V + ν1B − ( π3 D N + ν3 B I + ω2 T I + γ4 + ν2 + µ+ ζ1 T +B I ) A − ( π2 D N + γ9 + σ1 + µ+ ζ2 T +B I ) V ( π2V + π3A+ π4B + π5R+ π6P ) D N − ( β2 E + F + T N + γ1 + µ ) S β2S E + F + T N − ( β3 F + T N + γ2 + µ ) E β3E F + T N − ( γ3 + τ1 A+ V F + I + β4 T N + µ+ ζ3 A+ V F + I ) F ω4R B R+ I + ν3A B I + θ2P B P + I − ( π4 D N + ω1 T I + τ2 A+V I + γ7 + ν1 + µ+ ζ4 A+ V I ) B β4F T N + ω1B T I + ω2A T I + ω3R T R+ I + θ1P T P + I − ( τ3 A+ V I + γ8 + µ+ ζ5 A+ V I ) T τ1F A+ V I + τ2B A+ V I + τ3T A+ V I − ( π6 D N + θ1 T P + I + θ2 B P + I + γ5 + µ+ η ) P  (17) From (16) and (17), we have ⟨X(x),∇L(x)⟩ = − ( Λ + γ1S + γ2E + γ3F + γ4A+ γ5P + γ6R+ γ7B + γ8T + γ9V ) ≤ 0 (18) From (16) and (18), we infer that the sets { C ≥ 0 } , { R ≥ 0 } , { A ≥ 0 } , { V ≥ 0 } ,{ S ≥ 0 } , { E ≥ 0 } , { F ≥ 0 } , { B ≥ 0 } , { T ≥ 0 } , and { P ≥ 0 } are positively invariant, as established by the application of the barrier theorem. Hence, R10 ≥0 is positively invariant. Additionally, by the initial condition (11), we have x(0) ∈ R10 ≥0. Since R10 ≥0 is positively invariant, this ensures that all solutions of the system (1)−(10) remain positive for all time t > 0. □ 3.3. Invariant region Theorem 1. For initial conditions (11), the solutions of the system (1)−(10) are contained in the positively invariant, compact and attractive region Ψ = {( C, S,E, F, V,A,R,B, P, T ) ∈ R10 ≥0 : N(t) ≤ Λ µ } (19) Proof: Summing the equations (1) to (10), we find : M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 878 dN dt = Λ− µN − δ1V − δ2A− δ3B − δ4T − δ5F − ηP Since A, V,B, T, F, P are positive functions and using the positivity of the functions δ1, δ2, δ3, δ4, δ5, given that the constants ζ1, ζ2, ζ3, ζ4, ζ5 and η are strictly positive as well, we get: dN dt ≤ Λ− µN Then d dt ( N − Λ µ ) ≤ −µ ( N − Λ µ ) So the Gromwall inequality give N(t)− Λ µ ≤ ( N(0)− Λ µ ) e−µt Thus N(t) ≤ Λ µ + ( N(0)− Λ µ ) e−µt Since N(0) ≤ Λ µ , then 0 ≤ N(t) ≤ Λ µ . Therefore, all feasible solutions of the model (1)−(10) converge in the region Ψ. □ 4. Analysis of the equilibrium without terrorist, nor brigand nor fanatic x∗, and basic reproduction rumber R0 4.1. Equilibrium without terrorist, nor brigand, nor fanatic x∗ The uninfected compartments are C, R, A, V and the infected compartments are S, E, F, B, T, P. Given that we are at equilibrium without terrorists, fanatics, or robbers then we can discard the P compartment and the infected compartments being S, E, F, B, T, then an equilibrium solution with S = E = F = B = T=0 has the form: x∗ = ( C∗, R∗, A∗, 0, 0, 0, 0, 0, 0 ) (20) with C∗ = Λ(γ6 + µ)(γ4 + ν2 + µ) µ [ (γ6 + µ)(γ4 + µ+ ν2 + σ2) + σ2ν2 ] R∗ = Λν2σ2 µ [ (γ6 + µ)(γ4 + µ+ ν2 + σ2) + σ2ν2 ] A∗ = Λσ2(γ6 + µ) µ [ (γ6 + µ)(γ4 + µ+ ν2 + σ2) + σ2ν2 ] M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 879 4.2. Matrix of next generation K, and basic reproduction number R0 The Jacobian matrix of the system (1) − (10) is decomposed into Jx(x ∗) = DF(x∗) + DV(x∗) with F =  0 0 0 0( π1C + π2V + π3A+ π4B + π5R ) D N 0 0 α2C B C + I + ω4R B R+ I + ν3A B I β4F T N + ω1B T I + ω2A T I + ω3R T R+ I  and V =  Λ + γ1S + γ2E + γ3F + γ4A+ γ6R+ γ7B + γ8T + γ9V − ( π1 D N + α1 T +B C + I + α2 B C + I + σ2 + µ ) C ν2A− ( π5 D N + ω3 T R+ I + ω4 B R+ I + γ6 + µ ) R σ1V + σ2C + ν1B − ( π3 D N + ν3 B I + ω2 T I + γ4 + ν2 + µ+ ζ1 T +B I ) A α1C T +B C + I − ( π2 D N + γ9 + σ1 + µ+ ζ2 T +B I ) V − ( β2 E + F + T N + γ1 + µ ) S β2S E + F + T N − ( β3 F + T N + γ2 + µ ) E β3E F + T N − ( γ3 + τ1 A+ V F + I + β4 T N + µ+ ζ3 A+ V F + I ) F − ( π4 D N + ω1 T I + τ2 A+V I + γ7 + ν1 + µ+ ζ4 A+ V I ) B − ( τ3 A+ V I + γ8 + µ+ ζ5 A+ V I ) T  DF(x∗) = [ 0 0 0 F ] ; DV(x∗) = [ J1 J2 0 V ] with F = [ ∂Fi(x ∗) ∂xj ] 5≤i,j≤9 J1 = [ ∂Vi(x ∗) ∂xj ] 1≤i,j≤4 and J2 = [ ∂Vi(x ∗) ∂xj ] 1 ≤ i ≤ 4; 5 ≤ j ≤ 9 ; V = [ ∂Vi(x ∗) ∂xj ] 5≤i,j≤9 M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 880 Let: f = π1 C∗ C∗ +R∗ +A∗ + π5 R∗ C∗ +R∗ +A∗ + π3 A∗ C∗ +R∗ +A∗ g = α2 C∗ C∗ +A∗ + ω4 R∗ R∗ +A∗ + ν3 h = ω2 + ω3 R∗ R∗ +A∗ We get: F =  f f f 0 f 0 0 0 0 0 0 0 0 0 0 0 0 0 g 0 0 0 0 0 h  J1 =  −(σ2 + µ) γ6 γ4 γ9 0 −(γ6 + µ) ν2 0 σ2 0 −(γ4 + µ+ ν2) σ1 0 0 0 −(γ9 + µ+ σ1)  J2 =  γ1 γ2 − π1C ∗ γ3 − π3C ∗ ϖ1 ϖ2 −π5R ∗ −π5R ∗ −π5R ∗ −ω4 R∗ R∗ +A∗ ϖ3 −π3A ∗ −π3A ∗ −π3A ∗ ν1 − ν3 −π3A ∗ − ω2 0 0 0 α1 C∗ C∗ +A∗ α1 C∗ C∗ +A∗  with ϖ1 = γ7 − α1 C∗ C∗ +A∗ − α2 C∗ C∗ +A∗ ϖ2 = γ8 − π1 C∗ C∗ +R∗ +A∗ − α1 C∗ C∗ +A∗ ϖ3 = −π5 C∗ C∗ +R∗ +A∗ − ω3 R∗ R∗ +A∗ − ω4 R∗ R∗ +A∗ Note that J1 is a non-singular Metzler matrix (see [4]). M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 881 V =  −a 0 0 0 0 0 −b 0 0 0 0 0 −c 0 0 0 0 0 −d 0 0 0 0 0 −e  with a = γ1 + µ b = γ2 + µ c = γ3 + µ+ τ1 + ζ3 d = γ7 + µ+ τ2 + ν1 + ζ4 e = γ8 + µ+ τ3 + ζ5 We also note that V is a Metzler-Hurwitz matrix and V−1 =  − 1 a 0 0 0 0 0 −1 b 0 0 0 0 0 −1 c 0 0 0 0 0 −1 d 0 0 0 0 0 −1 e  V−1 =  − 1 a 0 0 0 0 0 −1 b 0 0 0 0 0 −1 c 0 0 0 0 0 −1 d 0 0 0 0 0 −1 e  ⇒ K = −FV−1 =  f a 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 g d 0 0 0 0 0 h e  where f a = (γ6 + µ) [ π1(γ4 + ν2 + µ) + π3σ2 ] + π5σ2ν2[ (γ6 + µ)(γ4 + µ+ ν2) + σ2ν2 ] (γ1 + µ) g d = ( 1 γ7 + τ2 + ν1 + µ+ ζ4 )( α2 γ4 + ν2 + µ γ4 + ν2 + σ2 + µ + ω4 ν2 γ6 + ν2 + µ + ν3 ) h e = ω2(γ6 + ν2 + µ) + ω3ν2 (γ6 + ν2 + µ)(γ8 + τ3 + µ+ ζ5) and R0 = ρ(K) = max { f a ; g d ; h e } (21) M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 882 Theorem 2. The equilibrium without terrorist, nor brigand nor fanatic x∗, is locally asymptotically stable if R0 < 1 and is unstable if R0 > 1. Proof : See [32, 33]. Thus R0 < 1, then an average radical indoctrinates or recruits less than one, which means that radical fanaticism as well as terrorism and robbery will disappear from this population over time. Conversely, if R0 > 1, then radical fanaticism and insecurity can spread through the population. Theorem 3. The equilibrium without terrorist, brigand or fanatic x∗, is globally asymp- totically stable if R0 < 1 and is unstable if R0 > 1. Proof : From Theorem 2, when R0 < 1, the compartments S,E, F,B, T tend to 0 as t → ∞. By setting S,E, F,B, and T to zero, it follows that (C,R,A, V, S,E, F,B, T, P ) → x∗ as t → ∞, where x∗ is the unique point in the positively invariant, compact, and attractive solution region Ψ, such that S = E = F = B = T = 0. □ 5. Global thresholds The (1) − (10) model is characterized qualitatively by a sequence of thresholds that are based on global dynamics. The first global threshold identifies the conditions for the extinction of the radical population. The second threshold identifies the conditions of the eradication of ideological terrorism. The third gives us the conditions of the extinction of fanatical ideology, and the fourth gives us the conditions of the eradication of insecurity without terrorism. Remark 1. In all this part, without loss of generality we note γi = γi + µ, ∀i ∈ {1, 2, 3, 4, 5, 6, 7, 8, 9} 5.1. A sufficient condition of the extinction of radical indoctrinated sub- populations Now we give a necessary and sufficient condition for extinction or stabilization of the radical core subpopulation. Theorem 4. Let λ1 = π1 + π2 + π3 + π4 + π5 + π6 + ω1 + ω2 + ω3 + θ1 and γ = min { γ1, γ2, γ3, γ8 } . So for all R1 = λ1 γ < 1, we have lim t→∞ D(t) = 0. Proof: We have D = S + E + F + T . Then: dD dt = dS dt + dE dt + dF dt + dT dt = ( π1C + π2V + π3A+ π4B + π5R+ π6P )D N + ( ω1B T I + ω2A T I + ω3R T R+ I + θ1P T P + I ) M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 883 − ( γ1S + γ2E + γ3F + γ8T ) − τ1 A+ V F + I F − ζ3 A+ V F + I F − τ3 A+ V I T − ζ5 A+ V I T ≤ ( π1 + π2 + π3 + π4 + π5 + π6 ) D + ( ω1 + ω2 + ω3 + θ1 ) D − ( γ1S + γ2E + γ3F + γ8T ) − τ1 A+ V F + I F − ζ3F − τ3 A+ V I T − ζ5T ≤ ( π1 + π2 + π3 + π4 + π5 + π6 + ω1 + ω2 + ω3 + θ1 ) D − γD ≤ ( λ1 − γ ) D The inequality λ1 < γ implies that the rate of decay γ exceeds the association strength λ1. Consequently, D, representing the radical core subpopulation, decreases exponentially to- wards zero. This decline is swift and inevitable, highlighting the vulnerability of the radical core subpopulation under these conditions. □ This result underscores a crucial insight: when the combined influence of association strength and recruitment capacity within the radical core subpopulation is outweighed by the rate of recovery or reintegration into civilian life, the radical core subpopulation is destined for extinction. In essence, if the mechanisms driving radicalization and recruit- ment cannot outpace the natural tendency of individuals to return to normal societal roles, the radical core subpopulation will inevitably perish. This insight underscores the signif- icance of addressing factors that promote radicalization and recruitment, as well as the importance of rehabilitation and reintegration efforts in countering extremist movements. 5.2. A sufficient condition of the eradication of ideological terrorism Considering the capacities of indoctrination to the fanatic ideology and the capacity of re- cruitment of the terrorist, we establish in this part a condition of eradication of ideological terrorism. Theorem 5. Let λ2 = β4 + ω1 + ω2 + ω3 + θ1 , λ3 = (τ3 + ζ5)κ+ γ8 with κ the infimum of A+ V I . So for all R2 = λ2 λ3 < 1, we have lim t→∞ T (t) = 0. Proof: From the equation (9) we have: dT dt = β4F T N + ω1B T I + ω2A T I + ω3R T R+ I + θ1P T P + I − ( τ3 A+ V I + γ8 + ζ5 A+ V I ) T = ( β4 F N + ω1 B I + ω2 A I + ω3 R R+ I + θ1 P P + I ) T − ( τ3 A+ V I + γ8 + ζ5 A+ V I ) T ≤ ( β4 F N + ω1 B I + ω2 A I + ω3 R R+ I + θ1 P P + I ) T − ( τ3κ+ γ8 + ζ5κ ) T ≤ ( β4 + ω1 + ω2 + ω3 + θ1 ) T − ( (τ3 + ζ5)κ+ γ8 ) T = (λ2 − λ3)T The last inequality implies that T decreases exponentially towards zero as soon as λ2 < λ3. Thus, when R2 = λ2 λ3 < 1, it provides a sufficient condition for the stabilization or eradica- tion of ideological terrorism. □ M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 884 This observation is crucial as it highlights the pivotal role of recruitment and conver- sion rates in the dynamics of ideological terrorism propagation. When the recruitment rate is lower than the conversion rate, the number of individuals joining the terrorist cause diminishes over time, eventually leading to the extinction of this form of terrorism. This underscores the importance of implementing policies and strategies aimed at reduc- ing recruitment opportunities and effectively countering extremist ideological propaganda. Consequently, controlling these parameters offers a potentially effective pathway to miti- gate and prevent the threats posed by ideological terrorism. The fact that the invariant superplane T = 0 is globally attractive for R2 < 1, we can then reduce the dimension of the model (1)− (10). In fact the model (1) − (10) is reduced to the equivalent nine dimensional system using the limit equation approach as follows: dC dt = Λ + γ1S + γ2E + γ3F + γ4A + γ5P + γ6R + γ7B + γ9V − ( π1 D N + α1 B C + I + α2 B C + I + σ2 ) C (22) dR dt = ν2A − ( π5 D N + ω4 B R + I + γ6 ) R (23) dA dt = σ1V + σ2C + ν1B − ( π3 D N + ν3 B I + γ4 + ν2 + ζ1 B I ) A (24) dV dt = α1C B C + I − ( π2 D N + γ9 + σ1 + ζ2 B I ) V (25) dS dt = ( π1C + π2V + π3A + π4B + π5R + π6P ) D N − ( β2 E + F N + γ1 ) S (26) dE dt = β2S E + F N − ( β3 F N + γ2 ) E (27) dF dt = β3E F N − ( γ3 + τ1 A + V F + I + ζ3 A + V F + I ) F (28) dB dt = α2C B C + I + ω4R B R + I + ν3A B I + θ2P B P + I − ( π4 D N + τ2 A + V I + γ7 + ν1 + ζ4 A + V I ) B (29) dP dt = τ1F A + V F + I + τ2B A + V I − ( π6 D N + θ2 B P + I + γ5 + η ) P (30) where D = S + E + F ; I = B +A+ V 5.3. A sufficient condition of the eradication of fanatical ideology without terrorism In this part we give a condition of eradication of fanatic ideology without terrorism. We first note the complexity of conducting an armed struggle against a fanatical ideology since it is difficult to identify a fanatical individual in a given society. We also note that some- times the need to identify fanatical individuals can lead to stigmatization and exactions in certain cases. Through the result below, we observe that inter-religious dialogue as well as values of good governance and good communication could have a very important weight in the factors that can lead to the extinction of a fanatical ideology by giving more weight to γ3 compared to β3 by also taking into account the effect of a citizen’s vigilance and police measures of proximity through τ1 and ζ3. M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 885 Theorem 6. Let λ4 = γ3 + (τ1 + ζ3)κ ′ with κ′ infimum of A+ V F + I . So for all R3 = β3 λ4 < 1, we have lim t→∞ F (t) = 0. Proof: From the equation (28) we have: dF dt = β3E F N − ( γ3 + τ1 A+ V F + I + ζ3 A+ V F + I ) F ≤ β3E F N − ( γ3 + τ1κ ′ + ζ3κ ′ ) F ≤ β3F − ( γ3 + (τ1 + ζ3)κ ′ ) F ≤ (β3 − λ4)F Hence F decreases exponentially to zero as β3 < λ4. □ In the same way R3 = β3 λ4 < 1, gives a sufficient condition of the stabilization or the erad- ication of the fanatical ideology. The fact also that the invariant superplane F = 0 is also a global attractor for R3 < 1, one can then reduce the dimension of the model (22)− (30). In fact the model (22)− (30) is reduced to the equivalent eight dimensional system using the limit equation approach as follows: dC dt = Λ+ γ1S + γ2E + γ4A+ γ5P + γ6R+ γ7B + γ9V − ( π1 D N + α1 B C + I + α2 B C + I + σ2 ) C (31) dR dt = ν2A− ( π5 D N + ω4 B R+ I + γ6 ) R (32) dA dt = σ1V + σ2C + ν1B − ( π3 D N + ν3 B I + γ4 + ν2 + ζ1 B I ) A (33) dV dt = α1C B C + I − ( π2 D N + γ9 + σ1 + ζ2 B I ) V (34) dS dt = ( π1C + π2V + π3A+ π4B + π5R+ π6P ) D N − ( β2 E N + γ1 ) S (35) dE dt = ( β2 S N − γ2 ) E (36) dB dt = α2C B C + I + ω4R B R+ I + ν3A B I + θ2P B P + I − ( π4 D N + τ2 A+ V I + γ7 + ν1 + ζ4 A+ V I ) B (37) dP dt = τ2B A+ V I − ( π6 D N + θ2 B P + I + γ5 + η ) P (38) where D = S + E; I = B +A+ V M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 886 5.4. A sufficient condition of the eradication of brigandage without ter- rorism nor ideological fanaticism. First of all, let us note that the fight against brigandage is one of the essential steps in the construction of a state. One of the indicators of a state’s strength lies in its ability to enforce its laws throughout its territory. In this section we give a condition for eradicating brigandage without terrorism nor ideological fanaticism. Theorem 7. Let λ5 = α2 +ω4 + ν3 + θ2 and λ6 = π4κ ′′ + τ2κ+ γ7 + ν1 + ζ4κ with κ, κ′′ respective infimum of A+ V I and of D N . So for all R4 = λ5 λ6 < 1, we have lim t→∞ B(t) = 0. Proof: From the equation (37) we have: dB dt = α2C B C + I + ω4R B R+ I + ν3A B I + θ2P B P + I − ( π4 D N + τ2 A+ V I + γ7 + ν1 + ζ4 A+ V I ) B ≤ ( α2 + ω4 + ν3 + θ2 ) B − ( π4κ ′′ + τ2κ+ γ7 + ν1 + ζ4κ ) B ≤ ( λ5 − λ6 ) B From this last inequality, B decreases exponentially if λ5 < λ6. □ This result underscores the critical notion that the effectiveness of countermeasures against brig- andage hinges on the balance between conversion forces within society and the state’s capacity to combat insecurity through citizen vigilance and effective governance. Specifically, if the cumulative force of conversion into brigandage, stemming from various sources such as non-indoctrinated civil- ians, discharged personnel from defense and security forces, individuals released from the justice system, and even the defense and security forces themselves, falls short of the state’s capability to address all forms of insecurity, including through citizen engagement and robust governance structures, then we are likely to witness either stabilization or complete eradication of brigandage. This observation underscores the interconnectedness between societal dynamics and state capacity in combating criminal activities like brigandage. It suggests that a holistic approach, encom- passing not only law enforcement and judicial efforts but also social, economic, and governance reforms, is essential for effectively tackling such security challenges. By addressing root causes, fostering community resilience, and ensuring robust governance frameworks, states can create an environment where brigandage finds little fertile ground for recruitment and operation. Therefore, achieving stability and eventual eradication of brigandage requires a concerted effort that goes beyond traditional security measures and encompasses broader societal and governance reforms. 6. Numerical simulation In this section, we present numerical results to provide further insights into the dynamics of the propagation of fanatic ideology. We conducted numerical simulations using MATLAB with the finite difference method, employing an explicit numerical scheme with a discretization step of 0.02. Figure 2 illustrates scenarios where the number of basic reproductions of extreme behavior R0 is less than 1, indicating a situation where the spread of fanatic ideology cannot sustain itself over time. Here, we observe the stabilization or extinction of various classes affected by the propagation of fanatic ideology, consistent with theoretical expectations. This indicates that under conditions M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 887 where the rate of propagation is lower than the rate of recovery, the spread of fanatic ideology tends to diminish over time, ultimately leading to its eradication. Conversely, Figure 3 depicts scenarios where R0 is greater than 1, indicating situations where the spread of fanatic ideology can persist and propagate within the population. In such cases, we observe the persistence of the affected classes, as the rate of propagation outweighs the rate of re- covery. This highlights the potential for sustained propagation of fanatic ideology under conditions where recruitment and indoctrination efforts outpace efforts to counter and rehabilitate individuals. These numerical simulations provide valuable insights into the dynamics of the propagation of fanatic ideology and underscore the importance of understanding and controlling key parameters such as recruitment, indoctrination, and recovery rates in mitigating its spread. Additionally, they emphasize the critical role of effective counterterrorism strategies and governance structures in ad- dressing and preventing the proliferation of radical movements. To initialize our simulations, we set the following initial conditions: C(0) = 150000, R(0) = 8, A(0) = 150, V (0) = 150, S(0) = 25000, E(0) = 1500, F (0) = 400, B(0) = 100, T (0) = 150, P (0) = 20. These initial conditions reflect the initial composition of the population in the different classes affected by the spread of fanatical ideology. In addition, the parameter values used in our simula- tions have been carefully estimated and are detailed in Table 1. Table 1: Parameter values estimated Parameters value for extinction value for persistence Λ 22500 22500 γ1 0.46 0.46 γ2 0.28 0.28 γ3 0.000111 0.00111 γ4 0.12 0.137 γ5 0.0000016 0.0016 γ6 0.026 0.09 γ7 0.002 0.06 γ8 0.000011 0.00011 γ9 0.011 0.011 π1 0.000007 0.0000006 π2 0.0000534 0.000534 π3 0.0000002 0.0000002 π4 0.000074 0.0074 π5 0.0000004 0.0000004 π6 0.44 0.44 θ1 0.000032 0.000032 θ2 0.000032 0.000032 η 0.15 0.15 ζ1 0.27 0.27 ζ2 0.27 0.27 ζ3 0.17 0.17 ζ4 0.17 0.47 ζ5 0.17 0.37 µ 0.08 0.08 ν1 0.002 0.02 ν2 0.0002 0.0002 ν3 0.2 0.2 τ1 0.2 0.2 τ2 0.45 0.45 τ3 0.45 0.45 β2 0.95 0.95 β3 0.72 0.72 β4 0.0078 0.0078 σ1 0.011 0.011 σ2 0.0001 0.0001 α1 0.0002 0.0002 α2 0.000534 4 ω1 0.25 1.2 ω2 0.8 0.2 ω3 0.38 1.2 ω4 0.38 2 M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 888 0 10 20 30 40 50 60 70 time(year) 1.5 2 2.5 3 ef fe ct if 105 Evolution of the subpopulation of the class C C 0 10 20 30 40 50 60 70 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 10 20 30 40 50 60 70 time(year) 0 1 2 3 4 5 6 7 8 ef fe ct if Evolution of the subpopulation of the class R R 0 10 20 30 40 50 60 70 time (year) 0 50 100 150 200 250 300 350 400 ef fe ct if Evolution of the subpopulation of the class F F 0 10 20 30 40 50 60 70 time(year) 80 90 100 110 120 130 140 150 e ff e c ti f Evolution of the subpopulation of the class A A 0 10 20 30 40 50 60 70 time (year) 0 10 20 30 40 50 60 70 80 90 100 ef fe ct if Evolution of the subpopulation of the class B B 0 10 20 30 40 50 60 70 time(year) 0 50 100 150 ef fe ct if Evolution of the subpopulation of the class V V 0 10 20 30 40 50 60 70 time (year) 0 50 100 150 e ff e c ti f Evolution of the subpopulation of the class T T 0 10 20 30 40 50 60 70 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 10 20 30 40 50 60 70 time(year) 0 20 40 60 80 100 120 140 ef fe ct if Evolution of the subpopulation of the class P P Figure 2: Evolution of the different classes of the model (1) − (10) with the extinction values. We get R0 = 0.3666, which is less than unity. M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 889 0 20 40 60 80 100 time(year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class C C 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time(year) 0 1 2 3 4 5 6 7 8 ef fe ct if Evolution of the subpopulation of the class R R 0 20 40 60 80 100 time (year) 200 400 600 800 1000 1200 1400 1600 1800 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time(year) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000 ef fe ct if Evolution of the subpopulation of the class A A 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time(year) 0 50 100 150 ef fe ct if Evolution of the subpopulation of the class V V 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class T T 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time(year) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 ef fe ct if Evolution of the subpopulation of the class P P Figure 3: Evolution of the different classes of the model (1) − (10) with the persistence values. We get R0 = 11.5538, which is greater than unity. M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 890 7. Mathematical analysis of a strategy to fight against terrorism, fanaticism and brigandage 7.1. Strategy to fight against terrorism, fanaticism and brigandage Based on the findings from our analysis, we apply optimal control theory to the model (1)− (10) in order to combat fanatical insurgency and brigandage. This involves introducing five time- dependent control variables, denoted as u1(t), u2(t), u3(t), u4(t), and u5(t), each representing specific strategies aimed at addressing radicalization, violent extremism, and insecurity. Let’s delve into the details of these strategies: (i) u1(t) represents a preventive strategy designed to hinder the indoctrination into radical ideology. This strategy may involve advocacy efforts to bolster social cohesion, socio-economic integration of vulnerable communities, educational initiatives, and facilitating employment op- portunities. Notably, u1(t) = 1 indicates the strategy’s effectiveness against radicalization, while u1(t) = 0 signifies its failure. The focus here is on significantly reducing the vulnerability of certain populations to radicalization. (ii) Similarly, u2(t) serves as another preventive strategy, which, in addition to the measures outlined in (i), emphasizes the state’s presence in vulnerable communities and endeavors to instill hope among disenchanted youth. This involves offering more attractive alternatives than those provided by radical groups, thereby fostering a sense of optimism. A value of u2(t) = 1 denotes the strategy’s effectiveness, while u2(t) = 0 indicates failure. (iii) u3(t) represents deradicalization strategies, which may entail intensified efforts towards inter-religious dialogue, community reconciliation, and engagement with radical groups by the state. A value of u3(t) = 1 signifies the effectiveness of the deradicalization strategy, while u3(t) = 0 indicates failure. (iv) The control variable u4(t) is a strategy aimed at combating organized crime, brigandage, and corruption. It encompasses police actions, such as community policing, investigations, and protection measures. This strategy also involves improving the capacity and equipment of defense and security forces and facilitating the rehabilitation and reintegration of offenders into society. A value of u4(t) = 1 denotes the effectiveness of the strategy against brigandage, while u4(t) = 0 signifies failure. (v) Finally, u5(t) represents the counterterrorism strategy, which, in addition to preventive measures and efforts against organized crime, focuses on disrupting terrorist financing and strength- ening the capabilities of defense and security forces. This entails enhancing firepower, intelligence systems, and executing coordinated actions to mitigate terrorist threats. A value of u5(t) = 1 indicates the strategy’s effectiveness against terrorism, while u5(t) = 0 denotes failure. 7.2. Mathematical analysis of strategy optimality Let ci(t) = 1− ui(t), ∀i ∈ {1, 2, 3, 4, 5}. (39) Consequently, the optimal control model with the five aforementioned time-dependent variables is given by the following differential equations M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 891  dC dt = Λ − γ1S + γ2E + γ3F + γ4A + γ5P + γ6R + γ7B + γ8T + γ9V − ( c1π1 D N + α1 T + B C + I + c4α2 B C + I + σ2 + µ ) C dR dt = ν2A − ( c1π5 D N + c5ω3 T R + I + c4ω4 B R + I + γ6 + µ ) R dA dt = σ1V + σ2C + ν1B − ( c1π3 D N + c4ν3 B I + c5ω2 T I + γ4 + ν2 + µ + ζ1 T + B I ) A dV dt = α1C T + B C + I − ( c1π2 D N + γ9 + σ1 + µ + ζ2 T + B I ) V dS dt = c1 ( π1C + π2V + π3A + π4B + π5R + π6P ) D N − ( c2β2 E + F + T N + γ1 + µ ) S dE dt = c2β2S E + F + T N − ( c3β3 F + T N + γ2 + µ ) E dF dt = c3β3E F + T N − ( γ3 + τ1 A + V F + I + c5β4 T N + µ + ζ3 A + V F + I ) F dB dt = c4 ( α2 C C + I + ω4 R R + I + ν3 A I + θ2 P P + I ) B − ( c1π4 D N + c5ω1 T I + τ2 A+V I + γ7 + ν1 + µ + ζ4 A + V I ) B dT dt = c5 ( β4 F N + ω1 B I + ω2 A I + ω3 R R + I + θ1 P P + I ) T − ( τ3 A + V I + γ8 + µ + ζ5 A + V I ) T dP dt = τ1F A + V F + I + τ2B A + V I + τ3T A + V I − ( c1π6 D N + c5θ1 T P + I + c4θ2 B P + I + γ5 + µ + η ) P (40) With non-negative initial conditions given by (11). This system can be written in matrix form as follows: X ′(t) = g(t,X, c) (41) Where X is defined in (13), c = (c1(t), c2(t), c3(t), c4(t), c5(t)) ∈ R5 verifies (39), and g : R×R10 × R5 → R10 is a nonlinear function which is written as in (14) but introducing the command c(t) so as to verify (40). The purpose of introducing the five control variables is to seek the optimal solution required to minimize the number of individuals in the radical subpopulation or core of violent extremism and fanatical behavior as well as brigands. Therefore, the objective function for this control problem is given by J (u1, u2, u3, u4, u5) = min 0⩽u1,u2,u3,u4,u5⩽1 ∫ Tf 0 ( j(t) + 1 2 k(t) ) dt (42) where j(t) = w1S(t) + w2E(t) + w3F (t) + w4B(t) + w5T (t) + w6P (t) k(t) = [ w7u 2 1(t) + w8u 2 2(t) + w9u 2 3(t) + w10u 2 4(t) + w11u 2 5(t) ] where the constants wi, i = 1, 2, ..., 11 are positive weights needed to balance the correspond- ing terms of the objective function. We choose quadratic costs on the orders, where 1 2 w7u 2 1(t), 1 2 w8u 2 2(t), 1 2 w9u 2 3(t), 1 2 w10u 2 4(t), 1 2 w11u 2 5(t) are the total cost of implementing the preventive mea- sure and the police-military response to manage active cases of armed insurgency over the time M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 892 interval [0, Tf ]. Specifically, we look for the optimal quintuple control u∗ = ( u∗ 1, u ∗ 2, u ∗ 3, u ∗ 4, u ∗ 5 ) is sought such that J ( u∗ 1, u ∗ 2, u ∗ 3, u ∗ 4, u ∗ 5 ) = min { J ( u1, u2, u3, u4, u5 ) : u1, u2, u3, u4, u5 ∈ U } , (43) where, U is the non-empty control set defined by U = {( u1, u2, u3, u4, u5 ) ∣∣∣∣ ui(t) is a piecewise continuous function on [0, Tf ] and 0 ⩽ ui ⩽ 1, ∀ ∈ t ∈ [0, Tf ], i = 1, 2, 3, 4, 5 } (44) Thus, to determine the necessary conditions that the optimal control quintuple must satisfy, we use Pontryagin’s maximum principle [29], which transforms the control problem (43) subject to model (40) into a pointwise minimization problem of a Hamiltonian H. This Hamiltonian is given by H = w1S + w2E + w3F + w4B + w5T + w6P + 1 2 [ w7u2 1(t) + w8u2 2(t) + w9u2 3(t) + w10u2 4(t) + w11u2 5(t) ] + λ1 [ Λ + γ1S + γ2E + γ3F + γ4A+ γ5P + γ6R+ γ7B + γ8T + γ9V − ( c1π1 D N + α1 T +B C + I + c4α2 B C + I + σ2 + µ ) C ] + λ2 [ ν2A− ( c1π5 D N + c5ω3 T R+ I + c4ω4 B R+ I + γ6 + µ ) R ] + λ3 [ σ1V + σ2C + ν1B − ( c1π3 D N + c4ν3 B I + c5ω2 T I + γ4 + ν2 + µ+ ζ1 T +B I ) A ] + λ4 [ α1C T +B C + I − ( c1π2 D N + γ9 + σ1 + µ+ ζ2 T +B I ) V ] + λ5 [ c1 ( π1C + π2V + π3A+ π4B + π5R+ π6P ) D N − ( c2β2 E + F + T N + γ1 + µ ) S ] + λ6 [ c2β2S E + F + T N − ( c3β3 F + T N + γ2 + µ ) E ] + λ7 [ c3β3E F + T N − ( γ3 + τ1 A+ V F + I + c5β4 T N + µ+ ζ3 A+ V F + I ) F ] + λ8 [ c4 ( α2 C C + I + ω4 R R+ I + ν3 A I + θ2 P P + I ) B − ( c1π4 D N + c5ω1 T I + τ2 A+ V I + γ7 + ν1 + µ+ ζ4 A+ V I ) B ] + λ9 [ c5 ( β4 F N + ω1 B I + ω2 A I + ω3 R R+ I + θ1 P P + I ) T − ( τ3 A+ V I + γ8 + µ+ ζ5 A+ V I ) T ] + λ10 [ τ1F A+ V F + I + τ2B A+ V I + τ3T A+ V I − ( c1π6 D N + c5θ1 T P + I + c4θ2 B P + I + γ5 + µ+ η ) P ] (45) where λi, i = 1, 2, ..., 10, represent the adjoint variables associated with the state variables of the model (40). The standard existence result for minimizing control problem as appeared in [19] is adapted as follows. M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 893 Theorem 8. (Existence and well-posedness of control problem) There exits an optimal quintuple control ( u∗ 1, u ∗ 2, u ∗ 3, u ∗ 4, u ∗ 5 ) ∈ U satysfing (42) subject to the control system (40) with non-negative initial conditions given by (11). Proof: The existence of an optimal control is obtained thanks to the Fleming and Rishel [19]. Thanks to a result of Lukes [25] which ensures the existence of solutions for the state system (40) with constant coefficients, the set of controls and corresponding solutions is non empty. In addition the set of controls U is a closed convex by definition and the vector field of the system (40) is bounded. Also the integrand of the objective function is clearly convex and g(t,X, c) in (41) is convex with respect to c. On the other hand there exist a1, a2 > 0 and β > 1 such that w1S + w2E + w3F + w4B + w5T + w6P + 1 2 [ w7u 2 1(t) + w8u 2 2(t) + w9u 2 3(t) + w10u 2 4(t) + w11u 2 5(t) ] ≥ a1 ( |u1|2 + |u2|2 + |u3|2 + |u4|2 + |u5|2 )β 2 − a2 since the state variables are bounded. Then, we deduce the existence of an optimal control (u∗ 1, u ∗ 2, u ∗ 3, u ∗ 4, u ∗ 5) that minimizes the objective function J(u1, u2, u3, u4, u5). □ Theorem 9. Given that ( u∗ 1, u ∗ 2, u ∗ 3, u ∗ 4, u ∗ 5 ) minimizes the objective functional (42) subject to the corresponding state system (40), then the adjoint variables λi, i = 1, 2, ..., 10, satisfy the following system: dλ1 dt = (λ1 − λ5)c1π1 D(N − C) N2 + (λ1 − λ4)α1 (T +B)I (C + I)2 + (λ1 − λ8)c4α2 BI (C + I)2 + (λ1 − λ3)σ2 + λ1µ +(λ5 − λ2)c1π5 DR N2 + (λ5 − λ3)c1π3 DA N2 + (λ5 − λ4)c1π2 DV N2 + (λ8 − λ5)c1π4 DB N2 + (λ5 − λ10)c1π6 DP N2 +(λ6 − λ5)c2β2 S(E + F + T ) N2 + (λ7 − λ6)c3β3 E(F + T ) N2 + (λ9 − λ7)c5β4 TF N2 dλ2 dt = (λ2 − λ1)γ6 + (λ2 − λ1)c1π1 DC N2 + (λ2 − λ1)c1π5 D(N −R) N2 + (λ2 − λ6)c5ω3 TI (R+ I)2 + λ2µ +(λ2 − λ8)c4ω4 BI (R+ I)2 + (λ5 − λ3)c1π3 DA N2 + (λ5 − λ4)c1π2 DV N2 + (λ5 − λ8)c1π4 DB N2 +(λ5 − λ10)c1π6 DP N2 + (λ6 − λ5)c2β2 S(E + F + T ) N2 + (λ7 − λ6)c3β3 E(F + T ) N2 + (λ9 − λ7)c5β4 TF N2 (46) M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 894 dλ3 dt = (λ3 − λ1)γ4 + (λ5 − λ1)c1π1 DC N2 + (λ4 − λ1)α1 (T +B)C (C + I)2 + (λ8 − λ1)c4α2 BC (C + I)2 + (λ3 − λ2)ν2 +(λ5 − λ2)c1π5 DR N2 + (λ9 − λ2)c5ω3 TR (R+ I)2 + (λ8 − λ2)c4ω4 BR (R+ I)2 + (λ3 − λ5)c1π3 D(N −A) N2 +(λ3 − λ8)c4ν3 B(V + T +B) I2 + (λ3 − λ9)c5ω2 T (V + T +B) I2 + λ3ζ1 (T +B)(V + T +B) I2 + λ3µ (λ5 − λ4)c1π2 DV N2 − λ4ζ2 (T +B)V I2 + (λ5 − λ8)c1π1 DB N2 + (λ5 − λ10)c1π6 DP N2 + (λ6 − λ5)c2β2 S(E + F + T ) N2 +(λ7 − λ6)c3β3 E(F + T ) N2 + (λ9 − λ7)c5β4 TF N2 + (λ7 − λ10)τ1 F (F + T +B) (F + I)2 + (λ8 − λ10)c4θ2 PB (P + I)2 +λ7ζ3 F (F + T +B) (F + I)2 + (λ8 − λ10)τ2 B(T +B) I2 + (λ9 − λ8)c5ω1 TB I2 + λ8ζ4 B(T +B) I2 +(λ9 − λ10)θ1 TP (P + I)2 + (λ9 − λ10)τ3 T (T +B) I2 + λ9ζ5 T (T +B) I2 dλ4 dt = (λ3 − λ1)γ9 + (λ5 − λ1)c1π1 DC N2 + (λ4 − λ1)α1 (T +B)C (C + I)2 + (λ8 − λ1)c4α2 BC (C + I)2 + (λ4 − λ3)σ1 +(λ5 − λ2)c1π5 DR N2 + (λ9 − λ2)c5ω3 TR (R+ I)2 + (λ8 − λ2)c4ω4 BR (R+ I)2 + (λ5 − λ3)c1π3 DA N2 +(λ8 − λ3)c4ν3 BA I2 + (λ9 − λ3)c5ω2 TA I2 − λ3ζ1 (T +B)A I2 + λ4µ+ (λ4 − λ5)c1π2 D(N − V ) N2 +λ4ζ2 (T +B)(A+ T +B) I2 + (λ5 − λ8)c1π4 DB N2 + (λ5 − λ10)c1π6 DP N2 + (λ6 − λ5)c2β2 S(E + F + T ) N2 +(λ7 − λ6)c3β3 E(F + T ) N2 + (λ9 − λ7)c5β4 TF N2 + (λ7 − λ10)τ1 F (F + T +B) (F + I)2 + (λ8 − λ10)c4θ2 PB (P + I)2 +λ7ζ3 F (F + T +B) (F + I)2 + (λ8 − λ10)τ2 B(T +B) I2 + (λ9 − λ8)c5ω1 TB I2 + λ8ζ4 B(T +B) I2 +(λ9 − λ10)θ1 TP (P + I)2 + (λ9 − λ10)τ3 T (T +B) I2 + λ9ζ5 T (T +B) I2 dλ5 dt = −w1 + (λ5 − λ1)γ1 + (λ1 − λ5)c1π1 C(N −D) N2 + (λ2 − λ5)c1π5 R(N −D) N2 + (λ3 − λ5)c1π3 A(N −D) N2 +(λ4 − λ5)c1π2 V (N −D) N2 + (λ8 − λ5)c1π4 B(N −D) N2 + (λ10 − λ5)c1π6 P (N −D) N2 + λ5µ +(λ5 − λ6)c2β2 (E + F + T )(N − S) N2 + (λ7 − λ6)c3β3 E(F + T ) N2 + (λ9 − λ7)c5π4 TF N2 M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 895 dλ6 dt = −w2 + (λ6 − λ1)γ2 + (λ1 − λ5)c1π1 C(N −D) N2 + (λ2 − λ5)c1π5 R(N −D) N2 + (λ3 − λ5)c1π3 A(N −D) N2 +(λ4 − λ5)c1π2 V (N −D) N2 + (λ8 − λ5)c1π4 B(N −D) N2 + (λ10 − λ5)c1π6 P (N −D) N2 + λ6µ +(λ5 − λ6)c2β2 S ( N − (E + F + T ) ) N2 + (λ6 − λ7)c3β3 (F + T )(N − E) N2 + (λ9 − λ7)c5π4 TF N2 dλ7 dt = −w3 + (λ7 − λ1)γ2 + (λ1 − λ5)c1π1 C(N −D) N2 + (λ2 − λ5)c1π5 R(N −D) N2 + (λ3 − λ5)c1π3 A(N −D) N2 +(λ4 − λ5)c1π2 V (N −D) N2 + (λ8 − λ5)c1π4 B(N −D) N2 + (λ10 − λ5)c1π6 P (N −D) N2 + λ7µ +(λ5 − λ6)c2β2 S ( N − (E + F + T ) ) N2 + (λ6 − λ7)c3β3 (F + T )(N − E) N2 + (λ7 − λ9)c5β4 T (N − F ) N2 +(λ7 − λ10)τ1 (A+ V )I (F + I)2 + λ7ζ3 (A+ V )I (F + I)2 dλ8 dt = −w4 + (λ8 − λ1)γ7 + (λ5 − λ1)c1π1 DC N2 + (λ5 − λ2)c1π5 DR N2 + (λ9 − λ2)c5ω3 TR (R+ I)2 + (λ8 − λ3)ν1 +(λ5 − λ10)c1π6 DP N2 + (λ5 − λ3)c1π3 DA N2 + (λ5 − λ4)c1π2 DV N2 + (λ1 − λ4)c4α1 C(C +A+ V ) (C + I)2 +(λ1 − λ8)c4α2 C(C +A+ V + T ) (C + I)2 + (λ2 − λ8)c4ω4 R(R+A+ V + T ) (R+ I)2 + λ8µ+ (λ3 − λ9)c5ω2 TA I2 +(λ3 − λ8)c4ν3 A(A+ V + T ) I2 + λ3ζ1 A(A+ V ) I2 + λ4ζ2 V (A+ V ) I2 + (λ6 − λ5)c2β2 S(E + F + T ) N2 +(λ7 − λ6)c3β3 E(F + T ) N2 + (λ9 − λ7)c5β4 TF N2 + (λ10 − λ7)τ1 F (A+ V ) (F + I)2 − λ7ζ3 F (A+ V ) (F + I)2 +(λ8 − λ10)τ2 (A+ V )(A+ V + T ) I2 + λ8ζ4 (A+ V )(A+ V + T ) I2 + (λ8 − λ10)c4θ2 P (P +A+ V + T ) (P + I)2 +(λ9 − λ8)c5ω1 TB I2 + (λ9 − λ10)θ1 TP (P + I)2 + (λ10 − λ9)τ3 T (A+ V ) I2 − λ9ζ5 T (A+ V ) I2 +(λ8 − λ5)c1π4 D(N −B) N2 + (λ8 − λ9)c5ω1 T (A+ V + T ) I2 M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 896 dλ9 dt = −w5 + (λ9 − λ1)γ8 + (λ1 − λ5)c1π1 C(N −D) N2 + (λ1 − λ4)c1α1 C(C +A+ V ) (C + I)2 + (λ8 − λ1)c4α2 BC (C + I)2 +(λ2 − λ5)c1π5 R(N −D) N2 + (λ2 − λ9)c5ω3 R(R+A+ V +B) (R+ I)2 + (λ8 − λ2)c4ω4 BR (R+ I)2 +(λ3 − λ5)c1π3 A(N −D) N2 + (λ8 − λ3)c4ν3 BA I2 + (λ3 − λ9)c5ω2 A(A+ V +B) I2 + λ3ζ1 A(A+ V ) I2 +(λ4 − λ5)c1π2 V (N −D) N2 + λ4ζ2 V (A+ V ) I2 + (λ5 − λ6)c2β2 S ( N − (E + F + T ) ) N2 +(λ6 − λ7)c3β3 E ( N − (F + T ) ) N2 + (λ10 − λ7)τ1 F (A+ V ) (F + I)2 + (λ7 − λ10)c5β4 F (N − T ) N2 − λ7ζ3 F (A+ V ) (F + I)2 +(λ8 − λ10)c4θ2 PB (P + I)2 + (λ8 − λ5)c1π4 B(N −D) N2 + (λ8 − λ9)c5ω1 B(A+ V +B) I2 − λ8ζ4 B(A+ V ) I2 +(λ10 − λ8)τ2 B(A+ V ) I2 + (λ10 − λ9)c5θ1 P (P +A+ V +B) (P + I)2 + (λ9 − λ10)τ3 (A+ V )(A+ V +B) I2 +λ9 + (λ10 − λ5)c1π6 P (N −B) N2 + λ9ζ5 (A+ V )(A+ V +B) I2 dλ10 dt = −w6 + (λ10 − λ1)γ5 + (λ5 − λ1)c1π1 DC N2 + (λ5 − λ2)c1π5 DR N2 + (λ5 − λ3)c1π3 DA N2 + (λ5 − λ4)c1π2 DV N2 +(λ5 − λ8)c1π4 DB N2 + (λ10 − λ5)c1π6 D(N − P ) N2 + (λ6 − λ5)c2β2 S(E + F + T ) N2 + (λ7 − λ6)c3β3 E(F + T ) N2 +(λ9 − λ7)c5β4 TF N2 + (λ10 − λ8)c4θ2 BI (P + I)2 + (λ10 − λ9)c5θ1 TI (P + I)2 + λ10µ+ λ10η with transversality conditions λi(Tf ) = 0, i = 1, 2, ..., 10. (47) Further, the optimal control quintuple ( u∗ 1, u ∗ 2, u ∗ 3, u ∗ 4, u ∗ 5 ) is given as follows M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 897  u∗ 1 = max { 0,min { 1, ( (λ5 − λ1)π1C + (λ5 − λ2)π5R + (λ5 − λ3)π3A + (λ5 − λ4)π2V + (λ5 − λ8)π4B + (λ5 − λ10)π6P ) D w7N }} u∗ 2 = max { 0,min { 1, (λ6 − λ5)β2S(E + F + T ) w8N }} u∗ 3 = max { 0,min { 1, (λ7 − λ6)β3E(F + T ) w9N }} u∗ 4 = max { 0,min { 1, (λ8 − λ1)α2 BC C + I + (λ8 − λ2)ω4 BR R + I + (λ8 − λ3)ν3 BA I + (λ8 − λ10)θ2 BP P + I w10 }} u∗ 5 = max { 0,min { 1, (λ9 − λ2)ω3 TR R + I + (λ9 − λ3)ω2 TA I + (λ9 − λ7)β4 TF N + (λ9 − λ8)ω1 TB I + (λ9 − λ10)θ1 TP P + I w11 }} (48) Proof: As mentioned earlier, the characterization of the optimal solution is obtained by ap- plying the Pontryagin’s maximum principle to the Hamiltonian of the system H. The system of ordinary differential equations (46) governing the adjoint variables is derived by differentiating the Hamiltonian. Further, the control characterizations in (48) are derived by solving, on the interior of the control set U , the partial differentials of the Hamiltonian H with respect to each of the controls u1, u2, u3, u4 and u5. Hence, by standard arguments involving control bounds, it follows that: u∗ 1 =  0 if r∗1 ≤ 0 r∗1 if 0 < r∗1 < 1 1 if r∗1 ≥ 1 u∗ 2 =  0 if r∗2 ≤ 0 r∗2 if 0 < r∗2 < 1 1 if r∗2 ≥ 1 u∗ 3 =  0 if r∗3 ≤ 0 r∗3 if 0 < r∗3 < 1 1 if r∗3 ≥ 1 u∗ 4 =  0 if r∗4 ≤ 0 r∗4 if 0 < r∗1 < 1 1 if r∗4 ≥ 1 u∗ 5 =  0 if r∗5 ≤ 0 r∗5 if 0 < r∗5 < 1 1 if r∗5 ≥ 1 M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 898 where, r∗1 = ( (λ5 − λ1)π1C + (λ5 − λ2)π5R+ (λ5 − λ3)π3A+ (λ5 − λ4)π2V + (λ5 − λ8)π4B + (λ5 − λ10)π6P ) D w7N r∗2 = (λ6 − λ5)β2S(E + F + T ) w8N r∗3 = (λ7 − λ6)β3E(F + T ) w9N r∗4 = (λ8 − λ1)α2 BC C + I + (λ8 − λ2)ω4 BR R+ I + (λ8 − λ3)ν3 BA I + (λ8 − λ10)θ2 BP P + I w10 r∗5 = (λ9 − λ2)ω3 TR R+ I + (λ9 − λ3)ω2 TA I + (λ9 − λ7)β4 TF N + (λ9 − λ8)ω1 TB I + (λ9 − λ10)θ1 TP P + I w11 This ends the proof. □ 7.3. Numerical simulations of the control system In this section we perform a numerical simulation to illustrate the impacts of counter-radicalization and counter-terrorism strategies on the evolutionary dynamics of the core sub-population of violent extremism. The results we present were obtained by an implementation on MATLAB using an explicit Euler scheme with a discretization step of 0.02. For the readers convenience, they should remembers that each of individuals without control is marked by red line. The individuals with control are marked by blue line. We use the persistence values for the parameters obtained in Table 1. In Figure 4 we illustrate the dynamics of the individuals of the classes S,E, F,B and T in the case where the preventive measures and the deradicalization strategies u1, u2 and u3 have a very low effectiveness while that of the fight against terrorism and brigandage have a high effectiveness. In Figure 5 we illustrate the dynamics of individuals in classes S,E, F,B, and T in the case where preventive measures and deradicalization strategies u1, u2, and u3 have very high effectiveness while those of counterterrorism and brigandage u4 and u5 have very low effectiveness. In Figure 6 we illustrate the dynamics of individuals in the S,E, F,B, and T classes in the case where preven- tive measures and deradicalization strategies u1, u2 and u3 have a very high effectiveness as well as those of counter-terrorism and brigandage u4 and u5. Figures 4 and 5 show that preventive measures and deradicalization strategies u1, u2, and u3 have a very large effect in combating radicalization and indoctrination to a fanatical ideology, but even with maximum effectiveness of u1, u2, and u3, terrorism and brigandage will persist if counterterrorism measures u4 and u5 are weak . While a high effectiveness of counter-terrorism measures, combining the military-police proximity measures u4 and u5, ensures stabilization of all subclasses of the radical core of the subpopulation. By contrast, from Figure 6, we see a faster regression in the S,E, F,B and T compartments compared to Figure 5, when preventive measures and deradicalization strategies (u1, u2, u3), as well as counterterrorism and brigandage combining proximity military-police measures (u4, u5) are very high. M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 899 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time (year) 200 400 600 800 1000 1200 1400 1600 1800 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time (year) 0 50 100 150 200 250 300 350 400 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time (year) 0 10 20 30 40 50 60 70 80 90 100 ef fe ct if Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class T T 0 20 40 60 80 100 time (year) 0 50 100 150 ef fe ct if Evolution of the subpopulation of the class T T Figure 4: Dynamics of individuals in the S,E, F,B, and T classes with and without control when u1 = u2 = u3 = 0 and u4 = u5 = 0.87 . M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 900 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time (year) 200 400 600 800 1000 1200 1400 1600 1800 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time (year) 0 50 100 150 200 250 300 350 400 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class T T 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 e ff e c ti f 10 4 Evolution of the subpopulation of the class T T Figure 5: Dynamics of individuals in the S,E, F,B, and T classes with and without control when u1 = u2 = u3 = 0.87 and u4 = u5 = 0.13 . M. R. Pooda, Y. Simpore, O. Traore / Eur. J. Pure Appl. Math, 17 (2) (2024), 870-904 901 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time (year) 0 0.5 1 1.5 2 2.5 ef fe ct if 104 Evolution of the subpopulation of the class S S 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time (year) 0 500 1000 1500 ef fe ct if Evolution of the subpopulation of the class E E 0 20 40 60 80 100 time (year) 200 400 600 800 1000 1200 1400 1600 1800 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time (year) 0 50 100 150 200 250 300 350 400 ef fe ct if Evolution of the subpopulation of the class F F 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time (year) 0 10 20 30 40 50 60 70 80 90 100 ef fe ct if Evolution of the subpopulation of the class B B 0 20 40 60 80 100 time (year) 0 2 4 6 8 10 12 14 16 18 ef fe ct if 104 Evolution of the subpopulation of the class T T 0 20 40 60 80 100 time (year) 0 50 100 150 ef fe ct if Evolution of the subpopulation of the class T T Figure 6: Dynamics of individuals in the S,E, F,B, and T classes with and without control when u1 = u2 = u3 = 0.87 and u4 = u5 = 0.87. REFERENCES 902 8. Conclusion In this study, we initially introduced the (1) − (10) model to characterize the dynamics of ideological terrorism, drawing inspiration from the terrorism situation in the Sahel region. This approach enabled us to highlight the influence of local and global security initiatives on the dynam- ics of ideological terrorism or fanatical insurgency. Within the framework of this model, we defined the threshold R0, representing the number of fundamental reproductions of extreme behavior, and established the conditions for local and global asymptotic stabilization of the equilibrium devoid of terrorism, brigandage or fanaticism. In addition, we identified four primordial thresholds R1, R2, R3, and R4, from which we de- rived conditions leading to the extinction of the radical core subpopulation with respect to violent extremism, ideological terrorism, fanatical ideology, and brigandage without terrorism or fanati- cism, respectively. Based on the results of the mathematical analysis of the model, we proposed a control strategy to combat ideological terrorism in the Sahel. We then studied the optimality of this control strategy using mathematical tools such as Pontryagin’s maximum principle. Finally, we presented a numerical simulation to complement our theoretical results. It is pertinent to rec- ognize that validation of the model’s predictions relies on available data, a challenge frequently encountered. The model predictions are straightforward and based on a theory of the propagation of extreme ideologies [9], as a function of initial conditions. Overall, these results improve our understanding of class evolution by elucidating how specific parameter values can lead either to the stabilization or elimination of extremism, or to its persistence. To enhance this research, it is worth exploring new mathematical perspectives. This includes examining the influence of social media on the spread of ideological extremism, analyzing the multi- level dynamics of extremism, and conducting social network analysis. By delving deeper into these areas, we can gain a better understanding of the complexity of ideological extremism and develop more effective control strategies. As part of our ongoing project, we are studying the spatio- temporal diffusion of ideological terrorism in the Sahel. We will explore new methods of numerical resolution, drawing on the work [1], [34]. Our aim is to analyze how extremist ideologies spread both spatially and temporally, taking into account factors such as social networks, geographical features, and political climates. Acknowledgements The authors would like to thank the The World Bank and the Higher Education Support Project in Burkina Faso for funding this study. References [1] Emad Az-Zo’bi. A reliable analytic study for higher-dimensional telegraph equation. J. Math. Comput. Sci, 18(4):423–429, 2018. [2] Alethea BT Barbaro, Lincoln Chayes, and Maria R D’Orsogna. Territorial developments based on graffiti: A statistical mechanics approach. Physica A: Statistical Mechanics and its Applications, 392(1):252–270, 2013. [3] Bijan Berenji, Tom Chou, and Maria R D’Orsogna. Recidivism and rehabilitation of criminal offenders: A carrot and stick evolutionary game. PloS one, 9(1):e85531, 2014. REFERENCES 903 [4] A Berman and RJ Plemmons. Nonnegative matrices in the mathematical sciences, ser. Classics in applied mathmatics. Philadelphia: SIAM, 1994. [5] Jean-Michel Bony. Principe du maximum, inégalité de harnack et unicité du probleme de cauchy pour les opérateurs elliptiques dégénérés. In Annales de l’institut Fourier, pages 277– 304, 1969. [6] Robert Boyd and Peter J Richerson. Culture and the evolutionary process. University of Chicago press, 1988. [7] Mario Bunge. Four models of human migration: An exercise in mathematical sociology. GEN. SYSTEMS, 16, 1971. [8] Jose Cadena, Gizem Korkmaz, Chris J Kuhlman, Achla Marathe, Naren Ramakrishnan, and Anil Vullikanti. Forecasting social unrest using activity cascades. PloS one, 10(6):e0128879, 2015. [9] Carlos Castillo-Chavez and Baojun Song. Models for the transmission dynamics of fanatic behaviors. In Bioterrorism: mathematical modeling applications in homeland security, pages 155–172. SIAM, 2003. [10] L Cavalli-Sforza, M Feldman, S Dornbusch, and K-H Chen. Cultural evolution: Anthropology and cultural transmission. Nature, 304(5922):124–124, 1983. [11] Luigi Luca Cavalli-Sforza and Marcus W Feldman. Cultural transmission and evolution: A quantitative approach. Princeton University Press, 1981. [12] Yao-Li Chuang, Tom Chou, and Maria R D’Orsogna. A network model of immigration: Enclave formation vs. cultural integration. arXiv preprint arXiv:1901.09396, 2019. [13] Lawrence E Cohen. Modeling crime trends: A criminal opportunity perspective. Journal of Research in Crime and Delinquency, 18(1):138–164, 1981. [14] DC Dennett. Darwin’s dangerous idea simon & schuster. New York, 1995. [15] Maria R D’Orsogna, Ryan Kendall, Michael McBride, and Martin B Short. Criminal defec- tors lead to the emergence of cooperation in an experimental, adversarial game. PloS one, 8(4):e61458, 2013. [16] Maria R D’Orsogna and Matjaž Perc. Statistical physics of crime: A review. Physics of life reviews, 12:1–21, 2015. [17] Joshua M Epstein. Modeling civil violence: An agent-based computational approach. Pro- ceedings of the National Academy of Sciences, 99(suppl 3):7243–7250, 2002. [18] Marcus W Feldman. Mathematical evolutionary theory. Princeton University Press, 1989. [19] Wendell H Fleming and Raymond W Rishel. Deterministic and stochastic optimal control, volume 1. Springer Science & Business Media, 2012. [20] Babak Fotouhi and Michael G Rabbat. Migration in a small world: A network approach to modeling immigration processes. In 2012 50th Annual Allerton Conference on Communica- tion, Control, and Computing (Allerton), pages 136–143. IEEE, 2012. REFERENCES 904 [21] Karl P Hadeler and Carlos Castillo-Chávez. A core group model for disease transmission. Mathematical biosciences, 128(1-2):41–55, 1995. [22] Rachel A Hegemann, Laura M Smith, Alethea BT Barbaro, Andrea L Bertozzi, Shannon E Reid, and George E Tita. Geographical influences of an emerging network of gang rivalries. Physica A: Statistical Mechanics and its Applications, 390(21-22):3894–3914, 2011. [23] Herbert W Hethcote and James A Yorke. Gonorrhea transmission dynamics and control, volume 56. Springer, 2014. [24] Mark I Lichbach. Nobody cites nobody else: mathematical models of domestic political conflict. Defence and Peace Economics, 3(4):341–357, 1992. [25] Dahlard L Lukes. Differential equations: classical to controlled. 1982. [26] Charles J Lumsden and Edward O Wilson. Genes, mind, and culture-the coevolutionary process. World Scientific, 2005. [27] Charles Z Marshak, M Puck Rombach, Andrea L Bertozzi, and Maria R D’Orsogna. Growth and containment of a hierarchical criminal network. Physical Review E, 93(2):022308, 2016. [28] Jie Pan and Anna Nagurney. Using markov chains to model human migration in a network equilibrium framework. Mathematical and computer modelling, 19(11):31–39, 1994. [29] Lev Semenovich Pontryagin. Mathematical theory of optimal processes. CRC press, 1987. [30] Laura M Smith, Andrea L Bertozzi, P Jeffrey Brantingham, George E Tita, and Matthew Valasik. Adaptation of an ecological territorial model to street gang spatial patterns in los angeles. Discrete and Continuous Dynamical Systems, 32(9):3223–3244, 2012. [31] Baojun Song. Dynamical epidemic models and their applications. Cornell University, 2002. [32] Richard S Varga. Matrix iterative analysis. prentice hall, englewood cliffs. New Jersey, 1962. [33] RS Varga. Factorization and normalized iterative methods, boundary problems in differential equations (r. langer, ed.), 1960. [34] Hamzeh Zureigat, Mohammad A Tashtoush, Ali F Al Jassar, Emad A Az-Zo’bi, Moham- mad W Alomari, et al. A solution of the complex fuzzy heat equation in terms of com- plex dirichlet conditions using a modified crank–nicolson method. Advances in Mathematical Physics, 2023, 2023.