Theoretical treatment and implementation of the SCM included Appell-Changhee polynomials for the fractional delayed carbon absorption-emission model

Results in Physics 58 (2024) 107459 Available online 15 February 2024 2211-3797/© 2024 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/bync-nd/4.0/). Contents lists available at ScienceDirect Results in Physics journal homepage: www.elsevier.com/locate/rinp Theoretical treatment and implementation of the SCM included Appell-Changhee polynomials for the fractional delayed carbon absorption-emission model M.M. Khader a,b, M. Adel c,∗, Muhammad Bilal Riazd,e, Hijaz Ahmad c,f,e,g aDepartment of Mathematics and Statistics, College of Science, Imam Mohammad Ibn Saud Islamic University (IMSIU), Riyadh, Saudi Arabia bDepartment of Mathematics, Faculty of Science, Benha University, Benha, Egypt cDepartment of Mathematics, Faculty of Science, Islamic University of Madinah, Medina, Saudi Arabia dIT4Innovations, VSB-Technical University of Ostrava, Ostrava, Czech Republic eDepartment of Computer Science and Mathematics, Lebanese American University, Byblos, Lebanon fNear East University, Operational Research Center in Healthcare, TRNC Mersin 10, Nicosia, 99138, Turkey gCenter for Applied Mathematics and Bioinformatics, Gulf University for Science and Technology, Mishref, Kuwait ARTICLE INFO MSC: 34A12 41A30 47H10 65N20 Keywords: Fractional delayed carbon absorption-emission model Caputo fractional derivative ACPs SCM Hopf bifurcation RKM4 ABSTRACT Some researchers have started to perform in-depth studies on carbon dioxide emissions phenomena because global warming has caused immense damage. In this work, as a generalization of these studies, we will investigate and describe the fractional delayed Carbon absorption-emission model which consists of two fractional (Caputo sense) differential equations. From this point, we study the behavior of the solution for this model by using an approximate technique based on Appell-type Changhee polynomials (ACPs). We present an approximation of the fractional-order derivative by using ACPs. We implement the spectral collocation method (SCM) and use the properties of the ACPs to convert the provided model into a set of algebraic equations. A special focus is placed on the examination of the existence and stability of the equilibrium point for the model according to the characteristic roots with the help of the Hopf bifurcation technique. We verify the effectiveness and efficiency of the presented numerical scheme by comparing its results to those of the fourth-order Runge–Kutta method (RKM4). The results confirm that the current method is a straightforward and efficient tool for simulating the solutions to such models. Introduction The prevailing trend is to safeguard the environment and conserve energy because unwise development has resulted in environmental destruction. Global cooperation encourages the rapid expansion of the world, but excessive use of natural resources has led to biodiversity loss, deforestation, and other environmental problems like land degradation. The world’s secondary sector is expanding quickly, which promotes a massive use of fossil fuels, three industrial wastes pollute the air, water, and land, and the greenhouse gases released from these processes intensify the greenhouse effect and contribute to global warming. We know that in 2020, the average worldwide temperature was around 1.2 degrees Celsius greater than it was before the Industrial Revolution. In September 2019, 66 nations came to the consensus that clear waterways and lush mountains are priceless resources and established the Climate Ambition Alliance. Global carbon neutrality has advanced as a ∗Corresponding author. E-mail addresses: [email protected] (M.M. Khader), [email protected] (M. Adel), [email protected] (M.B. Riaz), [email protected], [email protected] (H. Ahmad). result of all these actions. The legal requirements of countries include the objective of carbon neutrality. The EU said that it would become the first ‘‘carbon neutral’’ country ever in 2050. It is crucial to encourage peaking CO2emissions and achieve carbon neutrality since China is a significant carbon emitter. The secondary industry, which is mostly made up of manufacturing, has developed quickly due to the availability of abundant labor resources. The traditional manufacturing sector in China, with its high energy consumption and carbon emissions, continues to be the country’s main industry [1]. China is one of the biggest countries and plays a major role in the world economy, so created the Energy Conservation Law to stop the trend of global warming and take the lead among emerging nations. Therefore, there are many studies regarding China’s urban carbon dioxide emissions [1,2]. Researchers discovered that the realization of carbon neutralization in China will be considerably aided by the research [3]. https://doi.org/10.1016/j.rinp.2024.107459 Received 4 October 2023; Received in revised form 4 February 2024; Accepted 6 February 2024 Results in Physics 58 (2024) 107459 2 M.M. Khader et al. The third energy transformation’s leading character is now new energy, and it will steer the world towards carbon neutrality [4]. For China to reach peak CO2emissions and achieve carbon neutrality earlier than expected, it is crucial to build low-carbon cities, optimize industrial structure, cut carbon emissions from the steel industry, advance technology for reducing emissions, and lower the cost of carbon sequestration [5,6]. There are a few mathematical models on China’s carbon neutrality, for example, you can see [7–9]. We can use continuous differential equations to represent the timevarying processes of carbon release and absorption. Since carbon emission and absorption (CEA) are related to both the present and the past, we can also use the delayed differential equation (DDE) model to more completely and precisely represent the phenomenon of the dynamic system of the CEA and there are many research articles studied complex neural networks and DDEs (with fractional order) models [10– 12]. Because there are currently few research successes in this field, the main aim of this research is to use DDEs to characterize the CEA models. In this study, we will approximate the fractional derivative using the ACPs and use it to solve the proposed model using the SCM [13–20]. Additionally, we validate the effectiveness of the suggested algorithm by contrasting the outcomes with the RK4 results. Problem formulation Han et al. [3] obtained the following system of differential equations which depends on the conflicting nature of carbon output and absorption. They obtained that: 𝑑𝜓1(𝑡) 𝑑𝑡 =𝜓1(𝑡)(𝛼1−𝛼1𝛾1 𝑁2 𝜓2(𝑡)), 𝑑𝜓2(𝑡) 𝑑𝑡 =𝜓2(𝑡)(𝛼2−𝛼2𝛾2 𝑁1 𝜓1(𝑡)). (1) The entire parameters will be explained later. In [3], the authors postulated that the model’s realism would rise by including a carbon adsorption saturation component. The technology for reducing carbon emissions in 2024 would be more advanced than it would be when it reaches its full potential. As a result, the distance from the peak to the emissions of carbon before 2023 will have the biggest influence on the carbon emissions in 2024. Government investment in carbon emission reduction technologies will increase as these technologies are developed quickly, which will speed up the transformation of the industrial and energy structures and improve efficiency while they are being implemented. As a result, the authors in [3] established the following mathematical model: 𝜓1(𝑡) = 𝜓1(𝛼1−𝜖1𝜓2)+𝜆(M−𝜓1(𝑡−𝛿)), 𝜓2(𝑡) = 𝜓2(𝛼2−𝛽2𝜓1−𝜖2𝜓2),(2) where the dot refers to derivative w.r.t the time 𝑡; and 𝜖1=𝛼1𝛾1 𝑁2 , 𝛽2=𝛼2𝛾2 𝑁1 , 𝜖2=𝛼2 𝑁2 . The following is the description of the entire parameters: 𝜓1, 𝜓2: Quantity of carbon emissions and absorption, respectively. 𝛼1: Carbon emissions’ rate of annual growth. 𝛼2: Annual rate of carbon absorption growth. 𝜆: The energy structure reform’s impact factor on. M: Carbon emissions’ highest point. 𝑁1: The highest level of carbon emissions. 𝑁2: The maximum amount of carbon that can be absorbed. 𝛾1: The ratio of CEs to CAs that is competitive. 𝛾2: The competitive ratio of CEs to CAs. 𝛿: The delay in implementing technology to reduce the CEs in real production. Analysis of stability and the Hopf bifurcation Existence of a point of equilibrium In the case that the parameters of model (2) satisfy the hypotheses listed below [3]: 𝛶=√(𝜖1𝛼2+(𝜆−𝛼1)𝜖2)2− 4𝜖1𝜖2𝛽2𝜆M>0, 𝜖1𝛼2+(𝜆−𝛼1)𝜖2+𝛶 > 0, 𝜖1𝛼2−(𝜆−𝛼1)𝜖2−𝛶 > 0, (3) then the system (2) has a positive equilibrium point E=(𝜓1, 𝜓2), which is 𝜓1=𝜖1𝛼2+(𝜆−𝛼1)𝜖2+𝛶 2𝜖1𝛽2 , 𝜓2=𝜖1𝛼2−(𝜆−𝛼1)𝜖2−𝛶 2𝜖1𝜖2 . (4) Hopf bifurcation for Eand its stability There are many studies regarding bifurcation analysis and chaos control in one and two-dimensional spaces for fractional models [21– 23]. The stability of the model (2) at Ecan be analyzed by constructing the characteristic equation of (2), evaluated at E=(𝜓1, 𝜓2), as follows: 𝛬2+(−𝛼1+𝜖1𝜓2+𝜖2𝜓2)𝛬+𝜖2𝜓2(−𝛼1+𝜖1𝜓2)−𝜖1𝛽2𝜓1𝜓2 +(𝜆 𝛬 +𝜆 𝜖2𝜓2)𝑒−𝛬 𝛿 = 0.(5) When 𝛿= 0,(5) becomes 𝛬2+𝐴2𝛬+𝐵2= 0, where 𝐴2=(𝜆−𝛼1+𝜖1𝜓2+𝜖2𝜓2), 𝐵2=𝜖2𝜓2(𝜆−𝛼1+𝜖1𝜓2)−𝜖1𝛽2𝜓1𝜓2. When E=(𝜓1, 𝜓2)and the parameters of (2) satisfy (3), we can prove [3] 𝐴2>0, 𝐵2=𝜓2𝛶 > 0. Thus, Eis locally asymptotically stable (LAS) when 𝛿= 0. When 𝛿 > 0, we attempt to show whether Hopf bifurcation exists or not. We suppose that the Eq. (5) has a pure imaginary root, 𝛬=𝑖 𝜉 (𝜉 > 0, 𝑖 =√−1) is a real number. We get the following by substituting it into (5) and dividing the real and imaginary components: 𝜉2+𝜖2𝜓2(𝛼1−𝜖1𝜓2)+𝜖1𝛽2𝜓1𝜓2=𝜆 𝜉 sin(𝜉 𝛿) + 𝜆 𝜖2𝜓2cos(𝜉 𝛿), 𝜉(𝛼1−𝜖1𝜓2−𝜖2𝜓2)=𝜆 𝜉 cos(𝜉 𝛿) − 𝜆 𝜖2𝜓2sin(𝜉 𝛿).(6) Eq. (6) derives that equation given in Box I. Adding the square of the two equations, we obtain 𝜉4+𝐴3𝜉2+𝐵3= 0,(8) where 𝐴3= 2𝜖2𝜓2(𝛼1−𝜖1𝜓2)+ 2𝜖1𝛽2𝜓1𝜓2+(𝛼1−𝜖1𝜓2−𝜖2𝜓2)2−𝜆2, 𝐵3= 2𝜖1𝜖2𝛽2𝜓1𝜓2 2(𝛼1−𝜖1𝜓2)+(𝜖2𝜓2(𝛼1−𝜖1𝜓2))2 +(𝜖1𝛽2𝜓1𝜓2)2−(𝜆 𝜖2𝜓2)2. For convenience, we let 𝜉2=𝜙, then (8) becomes 𝑓(𝜙) = 𝜙2+𝐴3𝜙+𝐵3= 0.(9) When the parameters of (2) meet the assumption −𝐵3<0, then (9) has one root 𝜙2. If 𝐴3>0, 𝐵3>0hold, then (9) has no positive root. If 𝐴3<0, 𝐵3>0hold, then (9) has two positive roots 𝜙3, 𝜙4. Results in Physics 58 (2024) 107459 3 M.M. Khader et al. 𝑃2≜sin(𝜉 𝛿) = =𝜆 𝜉 (𝜉2+𝜖2𝜓2(𝛼1−𝜖1𝜓2)+𝜖1𝛽2𝜓1𝜓2)−𝜆 𝜖2𝜓2𝜉(𝛼1−𝜖1𝜓2−𝜖2𝜓2) 𝜆2(𝜉2+𝜖2 2𝜓2 2), 𝑆2≜cos(𝜉 𝛿) = =𝜆 𝜉2(𝛼1−𝜖1𝜓2+𝜖1𝛽2𝜓1𝜓2)+𝜆 𝜖2𝜓2(𝜉2+𝜖2𝜓2(𝛼1−𝜖1𝜓2)+𝜖1𝛽2𝜓1𝜓2) 𝜆2(𝜉2+𝜖2 2𝜓2 2). (7) Box I. We hypothesize that (9) has positive roots 𝜙𝑘(𝑘= 2,3,4), then 𝜉𝑘=√𝜙𝑘. From (7), we can address the important time delay value 𝛿(𝑗) 𝑘=⎧ ⎪ ⎨ ⎪ ⎩ 𝜉−1 𝑘(sin−1 (𝑃2)+ 2𝑗𝜋), 𝑆2>0, 𝜉−1 𝑘(− sin−1 (𝑃2)+ 2(𝑗+ 1)𝜋), 𝑆2<0, 𝑘 = 2,3,4, 𝑗= 0,1,2,…. (10) Let 𝛬=𝛬(𝛿)be the root of (5), satisfying 𝛬(𝛿(𝑗) 𝑘)=𝑖 𝜉𝑘(𝑘= 2,3,4). Differentiating both sides of (5) w.r.t 𝛿gives us the following Re (𝑑𝛬 𝑑𝛿 )−1|||||𝛿=𝛿(𝑗) 𝑘 = 2𝜉2 𝑘+ 2𝜖2𝜓2(𝛼1−𝜖1𝜓2)+ 2𝜖1𝛽2𝜓1𝜓2+(𝛼1−𝜖1𝜓2−𝜖2𝜓2)2−𝜆2 𝜆2(𝜉2 𝑘+𝜖2 2𝜓2 2) =𝑓′(𝜙𝑘)≠0. (11) Theorem 1 ([3]).Considering the stability of the point Efor system (2), we get: When (3) holds, we have: 1. If 𝐴3>0, 𝐵3>0, then Eis LAS for any 𝛿⩾0, and (9) has no positive root; 2. If 𝐵3<0holds, (9) has one positive roots 𝜙2, thus Eis LAS when 𝛿∈[0, 𝛿(0) 2)and unstable when 𝛿 > 𝛿(0) 2; 3. If 𝐴3<0, 𝐵3>0hold, (2) undergoes a Hopf bifurcation at E when 𝛿=𝛿(𝑗) 𝑘(𝑘= 3,4; 𝑗= 0,1,2,…). Then, ∃𝑛∈Nsuch that 0< 𝛿(0) 4< 𝛿(0) 3< 𝛿(1) 4<⋯< 𝛿(𝑛−1) 3< 𝛿(𝑛) 4< 𝛿(𝑛+1) 4. When 𝛿∈[0, 𝛿(0) 4]∪⋃𝑛 𝓁=1 (𝛿(𝓁−1) 3, 𝛿(𝓁) 4),Eof the model (2) is locally asymptotically stable, and when 𝛿∈⋃𝑛−1 𝓁=0 (𝛿(𝓁) 4, 𝛿(𝓁) 3)∪ (𝛿(𝑛) 4,+∞), the equilibrium Eis unstable. In this work, we will study the model (2) in its fractional (with order 𝜈) form as follows: 𝐷𝜈𝜓1(𝑡) = 𝜓1(𝛼1−𝜖1𝜓2)+𝜆(M−𝜓1(𝑡−𝛿)), 𝐷𝜈𝜓2(𝑡) = 𝜓2(𝛼2−𝛽2𝜓1−𝜖2𝜓2),(12) with the following initial conditions: 𝜓1(0) = 𝜓0 1, 𝜓2(0) = 𝜓0 2.(13) General concepts Fractional integration and derivative Here, the Caputo definition, which is one of the basic definitions in the fractional analysis field is given as follows. Definition 1. Following is the fractional derivative of the function 𝜓(𝑡), or 𝐷𝜈in Caputo’s interpretation [24] 𝐷𝜈𝜓(𝑡) = 1 𝛤(𝑛−𝜈)∫𝑡 0 𝜓(𝑛)(𝑠) (𝑡−𝑠)𝜈−𝑛+1 𝑑𝑠, 𝑛 − 1 < 𝜈 < 𝑛, 𝑛 ∈N.(14) The following characteristics apply to the Caputo fractional derivative 𝐷𝜈 𝐷𝜈𝐶= 0,where C is a constant, 𝐷𝜈𝑡𝓁=𝛤(𝓁+ 1) 𝛤(𝓁+1−𝜈)𝑡𝓁−𝜈,if 𝓁∈N∪ {0},and 𝓁≥⌈𝜈⌉.(15) It is a linear operator [25], i.e. 𝐷𝜈(𝑐1𝜙1(𝑡) + 𝑐2𝜙2(𝑡)) = 𝑐1𝐷𝜈𝜙1(𝑡) + 𝑐2𝐷𝜈𝜙2(𝑡), for some constants 𝑐1and 𝑐2. For more properties concerning this operator, the reader can be found in [24,25]. Changhee polynomials Changhee polynomials 𝐶ℎ𝑚(𝑡)are given by [26,27] 2 𝑧+ 2 (1 − 𝑧)𝑡= ∞ ∑ 𝑚=0 𝐶ℎ𝑚(𝑡)𝑧𝑚 𝑚!,(16) where 𝐶ℎ𝑚=𝐶ℎ𝑚(0) are the Changhee numbers [26]. The ACPs of degree 𝑚are given by 𝐶ℎ∗ 𝑚(𝑡) = 𝑚 ∑ 𝑗=0 (𝑚 𝑗)𝐶ℎ∗ 𝑚−𝑗𝑡𝑗.(17) It is easy to prove the following relationship in light of this formula (17) 𝑑 𝑑𝑡 𝐶ℎ∗ 𝑚(𝑡) = 𝑚 𝐶ℎ∗ 𝑚−1(𝑡).(18) From (18), we have 𝐶ℎ∗ 𝑚(𝑡) = ∫𝑡 0 𝑚 𝐶ℎ∗ 𝑚−1(𝑦)𝑑𝑦 +𝐶ℎ∗ 𝑚.(19) 𝐶ℎ∗ 0= 1 and for all 𝑚≥1we find 2𝐶ℎ∗ 𝑚+𝑚 𝐶ℎ∗ 𝑚−1 = 0. Since 𝑢(𝑡)∈L2[0,1], then 𝑢(𝑡) ≈ 𝑢∗(𝑡) = 𝑚 ∑ 𝑖=0 𝑐𝑖𝐶ℎ∗ 𝑖(𝑡).(20) ACPs-based approximation of the 𝑫𝝂 The following theorem to approximate 𝐷𝜈𝑢∗(𝑡)for the function 𝑢∗(𝑡) defined in (20). Theorem 2. Let 𝑢∗(𝑡)be given as shown in (20), then 𝐷𝜈𝑢∗(𝑡) = 𝑚 ∑ 𝑖=⌈𝜈⌉ 𝑖 ∑ 𝑗=⌈𝜈⌉ 𝑐𝑖𝜅𝑖,𝑗,𝜈 𝑡𝑗−𝜈,(21) Results in Physics 58 (2024) 107459 4 M.M. Khader et al. Fig. 1. The solution with various values of 𝜈. where 𝜅𝑖,𝑗,𝜈 is defined in terms of the Changhee number 𝐶ℎ∗ 𝑖−𝑗as follows 𝜅𝑖,𝑗,𝜈 = (𝑖)! 𝐶ℎ∗ 𝑖−𝑗 (𝑖−𝑗)! 𝛤(𝑗+1−𝜈). Proof. Consider the ACP, 𝐶ℎ∗ 𝑖(𝑡)of degree 𝑖, using (17), and (20) we get 𝐷𝜈𝑢∗(𝑡) = 𝑚 ∑ 𝑖=0 𝑐𝑖𝐷𝜈𝐶ℎ∗ 𝑖(𝑡) = 𝑚 ∑ 𝑖=⌈𝜈⌉ 𝑖 ∑ 𝑗=⌈𝜈⌉ 𝑐𝑖 (𝑖!)𝐶ℎ∗ 𝑖−𝑗 (𝑗!)(𝑖−𝑗)! 𝐷𝜈𝑡𝑗 = 𝑚 ∑ 𝑖=⌈𝜈⌉ 𝑖 ∑ 𝑗=⌈𝜈⌉ 𝑐𝑖 (𝑖!) 𝐶ℎ∗ 𝑖−𝑗 (𝑖−𝑗)! 𝛤(𝑗+1−𝜈)𝑡𝑗−𝜈 = 𝑚 ∑ 𝑖=⌈𝜈⌉ 𝑖 ∑ 𝑗=⌈𝜈⌉ 𝑐𝑖𝜅𝑖,𝑗,𝜈 𝑡𝑗−𝜈, (22) where 𝜅𝑖,𝑗,𝜈 defined in (21), and this finishes the proof. □ Numerical implementation Now, we give an outline of the implementation of the technique that will solve the model under study. Thus, the fractional delayed carbon absorption-emission model with Caputo fractional-order derivative 0< 𝜈≤1, is given in (12)–(13).𝜓1(𝑡), 𝜓2(𝑡)can be approximated using ACPs by 𝜓1,𝑚(𝑡),𝜓2,𝑚(𝑡), as follows: 𝜓1,𝑚(𝑡) = 𝑚 ∑ 𝑖=0 𝑎𝑖𝐶ℎ∗ 𝑖(𝑡), 𝜓2,𝑚(𝑡) = 𝑚 ∑ 𝑖=0 𝑏𝑖𝐶ℎ∗ 𝑖(𝑡).(23) From (21) and (23) in (12), we get 𝐷𝜈𝜓1(𝑡) = (𝜓1,𝑚(𝑡))(𝛼1−𝜖1𝜓2,𝑚(𝑡))+𝜆(M−𝜓1,𝑚(𝑡−𝛿)),(24) 𝐷𝜈𝜓2(𝑡) = (𝜓2,𝑚(𝑡))(𝛼2−𝛽2𝜓1,𝑚(𝑡) − 𝜖2𝜓2,𝑚(𝑡)).(25) By collocation (24)–(25) at 𝑚−1 points 𝑡𝑟=𝑟 𝑚−1 +1, 𝑟 = 1,2,…, 𝑚−1, we have: 𝐷𝜈𝜓1(𝑡𝑟) = (𝜓1,𝑚(𝑡𝑟))(𝛼1−𝜖1𝜓2,𝑚(𝑡𝑟))+𝜆(M−𝜓1,𝑚(𝑡𝑟−𝛿)),(26) 𝐷𝜈𝜓2(𝑡𝑟) = (𝜓2,𝑚(𝑡𝑟))(𝛼2−𝛽2𝜓1,𝑚(𝑡𝑟) − 𝜖2𝜓2,𝑚(𝑡𝑟)).(27) Substitute Eq. (23) into (13) and since 𝐶ℎ∗ 𝑖(0) = 𝐶ℎ∗ 𝑖, then (13) becomes the following algebraic equations: 𝑚 ∑ 𝑖=0 𝑎𝑖𝐶ℎ∗ 𝑖=𝜓0 1, 𝑚 ∑ 𝑖=0 𝑏𝑖𝐶ℎ∗ 𝑖=𝜓0 2.(28) To get 𝑎𝑖, 𝑏𝑖, 𝑖 = 0,1,…, 𝑚, we solve the nonlinear system of equations (26)–(28) using the Newton iteration method. In turn, this prompts us to construct the rough solution by replacement in (23). Numerical simulation By presenting simulations on a test case in a closed interval, we will be able to confirm the validity and quality of the provided scheme. We study the model (12)–(13) at various values of 𝜈, 𝑚, 𝛿, 𝜆, 𝛾1, 𝛾2. In all Figures, we take the same values for the constants of [3] as follows 𝛼1= 0.14, 𝛼2= 0.2,M= 138.135, 𝑁1= 276.27, 𝑁2= 151.9485. We consider the following initial conditions [3]: 𝜓0 1= 132, 𝜓0 2= 50. To assess the precision and caliber of the suggested scheme, we also compare the results given using the presented approach with those given by using the RKM4. We also evaluate the residual error function (REF) [28] to check the accuracy and quality of the proposed scheme. Figs. 1–7provide the numerical findings for the investigated model that were produced by using the suggested technique. 1. In Fig. 1, we give the solution with 𝜈= 1.0,0.95,0.85,0.75, with 𝑚= 7, 𝛿 = 0.2, 𝜆 = 0.5, 𝛾1= 2, 𝛾2= 0.5. 2. In Fig. 2, we give the solution with 𝜆= 0.4,0.6,0.8,1, with 𝑚= 8, 𝜈 = 0.95, 𝛿 = 0.2, 𝛾1= 2, 𝛾2= 0.5. 3. In Fig. 3, we present a comparison between the solutions with those that were obtained by the RKM4 at (𝜈= 1) with 𝑚= 8, 𝛿 = 0.2, 𝛾1= 2, 𝛾2= 0.5. 4. In Fig. 4, we study the solution with 𝑚= 5 (a,c) and 𝑚= 10 (b,d) by computing the REF; with 𝜈= 0.95, 𝜆 = 0.5,𝛿= 0.2, 𝛾1= 2, 𝛾2= 0.5. 5. In Fig. 5, the solution with 𝛿= 0.2,0.4,0.6,0.8, and 𝜈= 0.95, 𝑚 = 8, 𝜆 = 0.5,𝛿= 0.2, 𝛾1= 2, 𝛾2= 0.5is presented. 6. In Fig. 6, the solution at (𝜈= 0.96) with 𝛾1= 2,2.5,3with 𝑚= 8, 𝜆= 0.5,𝛿= 0.2, 𝛾2= 0.5. 7. In Fig. 7, the solution at (𝜈= 0.94) with 𝛾2= 0.5,0.75,1with 𝑚= 8,𝜆= 0.5,𝛿= 0.2, 𝛾1= 0.5. These findings show that the recommended technique can solve the suggested model in its fractional form with the Caputo sense and that the behavior of the numerical solution generated by the suggested method depends on the values of 𝜈, 𝑚, 𝛿, 𝜆, 𝛾1, 𝛾2. In the end, the presented method improves the efficiency of the method and its outputs significantly. Conclusions Our target is to investigate the dynamical behavior of the delayed carbon absorption-emission model in its fractional form. We studied the behavior of the solution for this model by using the SCM based on ACPs. We presented an approximation for the fractional derivative by utilizing ACPs. A special focus is placed on the examination of the existence and stability of the equilibrium point for the system according to the Results in Physics 58 (2024) 107459 5 M.M. Khader et al. Fig. 2. The solution with various values of 𝜆. Fig. 3. The solution by ACPs and RK4 methods 𝜈= 1. Fig. 4. The REF of 𝜓1(𝑡), and 𝜓2(𝑡)for various values of 𝑚= 5(𝑎, 𝑐)and 𝑚= 10(𝑏, 𝑑). characteristic roots with the help of the Hopf bifurcation technique. We determined the numerical solutions of this model for various values of 𝜈, and 𝑚with the help of the proposed technique. Then it is possible to confirm the suitability of the proposed method to successfully study this model. In addition, by including more terms from the series of approximation solutions, we can control and reduce the error. So, we conclude that the numerical simulation using the Caputo operator is the most convenient for the model under consideration in this paper. The results obtained using the RK4 approach are consistent with the results obtained in graphical form from the presented method. These results also show the accuracy and effectiveness of the proposed method in computing. Through this study and numerical simulation, we can confirm that when the natural growth rate of carbon absorption increases, the time it takes to achieve the carbon peak in CO2emissions will be shortened and then the peak value will also decrease. Since the natural growth rate of carbon absorption is very high, to reduce this rate we need to deepen industrial reform and improve the energy structure. Therefore, Results in Physics 58 (2024) 107459 6 M.M. Khader et al. Fig. 5. The solution with various values of 𝛿. Fig. 6. The solution with various values of 𝛾1. Fig. 7. The solution with various values of 𝛾2. based on these results, we must recommend increasing tree planting and improve the level of carbon storage technology, as well as improving the development and application of new energy technology to achieve and adjust the depth of the industrial structure and improve the energy structure. These recommendations and procedures are of utmost importance to avoid the risks resulting from global warming resulting from increased carbon dioxide. We will try to deal with the problem under study in the future utilizing another type of fractional derivative or another kind of polynomial that may be used as a generalization of our findings. In the end, it is possible to address some other parameters that affect the model and add to it to make it deeper and more realistic; and it is considered a generalization of this model and this study. CRediT authorship contribution statement M.M. Khader: Supervision, Validation, Visualization, Writing – original draft. M. Adel: Conceptualization, Data curation, Formal analysis, Funding acquisition. Muhammad Bilal Riaz: Funding acquisition, Investigation, Methodology, Project administration. Hijaz Ahmad: Writing – review & editing, Data curation, Conceptualization. Declaration of competing interest The authors declare that there is no conflict of interests regarding the publication of this paper. Data availability Data will be made available on request. Acknowledgment The Author Muhammad Bilal Riaz is highly obliged to Ministry of Education, Youth and Sports of the Czech Republic for their support through the e-INFRA CZ (ID: 90254). References [1] Zhang GX, Zhang PD, Xiu J, Chai J. Chin J Popul Resour 2018;1:12–27. [2] Cai BF, Wang JN, Yang SY, Mao YQ, Cao LB. Chin J Popul Resour Environ 2017;15:58–70. [3] Han L, Sui H, Ding Y. Mathematics 2023;10:1–21. [4] Zou CN, Xiong B, Xue HQ, Zheng DW, Ge ZX, Wang Y, et al. Pet. Explor Dev 2021;48:480–91. Results in Physics 58 (2024) 107459 7 M.M. Khader et al. [5] Zhang YX, Yin JP, Fu Y, Wang JX, Wang GS. Meteorol Env Res 2021;12:24–31. [6] Feng R, Dai BZ, Su H. Ecol Econom 2017;2:117–24. [7] Wang SJ, Gao S, Huang YY, Shi CY. Acta Geogr Sin 2020;30:757–74. [8] Muhammad KK, Muhammad IK, Muhammad R. Financ Innov 2020;6:56–68. [9] Liu K, Tao YM, Wu Y, Wang CX. Chin J Popul Resour 2020;18:97–102. [10] Panda SK, Vijayakumar V, Nagy AM. Chaos Solitons Fractals 2023;177:114263. [11] Panda SK, Nagy AM, Vijayakumar V, Hazarika B. Chaos Solitons Fractals 2023;175:114045. [12] Panda SK, Abdeljawad T, Jarad F. Results Phys 2023;177:106313. [13] Khader MM, Abualnaja KM. J Appl Anal Comput 2019;9:261–70. [14] Khader MM, Sharma RP. Math Comput Simulation 2021;181:333–50. [15] Adel M, Khader MM, Assiri TA, Kaleel W. Symmetry 2023;15(4):931. [16] Khader MM, Adel M. Fractal Fract 2022;6(363):1–19. [17] Saad KM, Khader MM, Gomez-Aguilar JF, Baleanu D. Chaos 2019;29:1–5. [18] Ibrahim YF, Abd El-Bar SE, Khader MM, Adel M. Fractal Fract 2023;7(4):307. [19] Adel M, Srivastava HM, Khader MM. Math Methods Appl Sci 2023;46:8362–71. [20] Khader MM, Saad KM. Int J Biomath 2018;11(8):1–15. [21] Almatrafi MB, Berkal M. Int J Anal Appl 2023;21:131. [22] Berkal M, Almatrafi MB. Fractal Fract 2023;7:344. [23] Almatrafi MB. Fractals 2023;31(10):2340160. [24] Kilbas AA, Srivastava HM, Trujillo JJ. Theory and applications of fractional differential equations. Elsevier; 2006, p. 204. [25] Saadatmandi A, Dehghan M. Comput Math Appl 2010;59(3):1326–36. [26] Kim DS, Kim T, Seo JJ. Adv Stud Theor Phys 2013;7(20):993–1003. [27] Lee JG, Jang LC, Seo JJ, Choi SK, Kwon HI. Adv Difference Equ 2016;1:1–10. [28] El-Hawary HM, Salim MS, Hussien HS. J Global Optim 2003;25(3):283–303.