Mathematical modeling and computational analysis of hepatitis B virus transmission using the higher-order Galerkin scheme

1 Introduction

Hepatitis B, a severe and potentially fatal viral infection of the liver, represents a considerable challenge to global health. The hepatitis B virus (HBV) is the primary cause for this disease, which can result in both acute and chronic liver damage, leading to significant morbidity and mortality. HBV transmission predominantly occurs through the exchange of bodily fluids, such as blood, semen, and vaginal secretions, with adults generally encountering a self-limiting infection. However, a considerable number of individuals may progress to chronic hepatitis B, potentially resulting in liver cirrhosis, hepatocellular carcinoma, and various other serious complications. The precision of HBV diagnosis and screening is critical for identifying individuals with acute or chronic infections, as well as facilitating prompt interventions and preventive measures. Mathematical modeling is utilized as a fundamental tool for elucidating the dynamics and attributes of diverse diseases. Numerous groups of researchers have constructed models for various diseases, including human immunodeficiency virus (HIV) infection, hepatitis B, cancer, and infectious diseases, facilitating more comprehensive comprehension, prediction, and management of their dissemination and effects. Yavuz et al. [1] introduced a novel fractional model that explained the behaviors of the HBV. The study involved the modeling, equilibria, stabilities, and the analysis of existence and uniqueness pertaining to the model, which was afterward solved numerically employing the Adams–Bashforth scheme. Mustapha et al. [2] developed a systematic deterministic mathematical model that illustrates the transmission dynamics of HBV. The model parameters were estimated utilizing actual data from South Africa through the application of the nonlinear least-squares curve fitting method. The numerical analysis conducted confirmed that minimizing the transition of chronically infected children to adulthood through treatment is essential for the eradication of hepatitis B in South Africa. Naik et al. [3] carried out an investigation into the global dynamics of a fractional-order model for the hepatitis C virus (HCV), incorporating an effective treatment strategy. The model’s fixed points were identified, followed by a comprehensive stability analysis. The fractional Adams method was employed in obtaining the numerical solutions of the model. In the end, they validated the numerical simulations compared to the theoretical and numerical results. Din and Abidin [4] performed an investigation into a recently developed model for hepatitis B infection, utilizing the Atangana–Baleanu–Caputo fractional order derivative. The Ulam–Hyers-type stability was achieved, along with a qualitative analysis of the corresponding solution, through the application of a well-established principle from fixed point theory. The results obtained have been verified by numerical verification through MATLAB. Bolaji et al. [5] implemented a deterministic compartmental co-infection model to elucidate the effects of tuberculosis (TB) infection on the co-infection dynamics of the two diseases, particularly in environments where treatment for TB is readily accessible. The four-dimensional systems of ordinary differential equations (ODEs) constructed by Joshi et al. [6] quantify the impact of burnt and recycled plastic on air pollution. The well-posedness and qualitative features were explored. The plastic waste model’s reproduction number and local/global stability were thoroughly examined. Naik et al. [7] conducted a thorough investigation and analysis of the fractional-order epidemic model concerning the mutual influence in HIV/HCV co-infection. Farman et al. [8] explored an unconditionally convergent semi-analytical method utilizing recent evolutionary computational techniques, alongside evolutionary Padé approximation (EPA), for addressing the complexities of a nonlinear hepatitis B model. Naik et al. [9] performed an analysis of a fractal-fractional operator within the context of an epidemic model, specifically pertaining to COVID-19 modeling, and examined the variations in the infection rate across society. Attaullah et al. [10] formulated a higher-order Galerkin time discretization and conducted numerical comparisons for two models of HIV infection. In order to substantiate the findings, the model was assessed through the application of the Runge–Kutta fourth-order (RK4) method. Furthermore, they examined the impact of diverse physical attributes by altering their values and analyzing them through graphical representations. Jan et al. [11] developed a model for the dynamics of asthma, incorporating smoking and environmental factors within the fractional Caputo–Fabrizio (CF) framework to elucidate its dynamic behavior. Attaullah et al. [12] examined the dynamical behavior of HIV infection, incorporating the influence of a variable source term via the Galerkin method, and elucidated the characteristics of HIV infection. HBV remains a pressing global health challenge, prompting researchers to investigate treatment responses and identify susceptibility factors to inform prevention and intervention strategies. Several researchers looked at how well hepatitis B patients responded to treatment and tried to figure out who was most likely to get infected. They put a lot of effort into making a mathematical model to help them understand how HBV spreads. They attempted to enhance the existing knowledge of hepatitis B by providing detailed mathematical frameworks that illustrate the disease’s transmission dynamics, guide public health policy, and refine treatment approaches. They discussed comprehension of the intricate relationships among HBV transmission, host immunity, and treatment outcomes, thereby addressing a critical public health challenge and enhancing patient outcomes. Mohideen and Lakshmi [13] presented a stochastic gradient algorithm to find the co-connection and to find the resulting accuracy of acute as well as chronic prediction of patients. Din et al. [14] developed a model of HBV under computational fluid dynamics and used domain polynomials to find out the numerical results of the proposed model. Reza et al. [15] suggested the dynamics of patients with HBV vaccination and treatments. The spread of the disease may be reduced by vaccination and spreading awareness. Shang and Kao [16] discussed and summarized the relevant issues of HBV. Nana et al. [17] focused on the treatment model of liver disease and critically investigated the dynamics of the proposed model. Okosun et al. [18] investigated the dynamics of local as well as global stability of HBV by using the Lyapunov function. Lashari et al. [19] investigated direct and indirect transmissions by using the backward bifurcation technique. Forde et al. [20] presented an optimal control strategy for controlling the immunity response against the HBV and showed that high vaccination induces immunity to control the infection. Ntaganda et al. [21] used the fuzzy logic strategy for hepatitis B by applying an optimal control technique and considering infants up to the age of 12 months. Zhang et al. [22] suggested a model for the transmission of HBV in high endemic areas with vertical transmission. Global stability for the HBV model was derived via the Lyapunov function. Anna and Susanna [23] summarized the HBV therapies with limitations, novel approaches and therapeutic points for the cure of HBV. Van den et al. [24] demonstrated the basic reproduction number R 0 for compartmental transmission of infection. Zou et al. [25] proposed a mathematical model to prevail hepatitis B transmission in China by taking R 0 = 2.406 . Reza et al. [26] discussed the vaccination and treatment for the hepatitis B dynamical infection. A R 0 was also derived and analyzed for the infection-free steady state. Hews et al. [27] considered the model as explicit time delay in the production of HBV by providing the lamivudine therapy which showed best agreement with clinical data. Ntaganda [28] presented the HBV dynamics control problem by using Pontryagin’s maximum principle as well as a direct approach. Armbruster and Brandeau [29] evaluated the effect and cost of the chronic infection for different screening levels and contact tracing. Hattaf et al. [30] aimed at the usage of optimal control for the system of the ODE model of hepatitis B disease. Optimal control was derived by using Pontryagin’s maximum principle. Ntaganda and Gahamanyi [31] employed the fuzzy logic technique to determine the optimal control solution. A numerical comparison was derived through the direct method. Bhattacharyya and Ghosh [32] suggested that hepatitis B infection vaccination and antiviral treatment are necessary to be in a proper ratio to take control of HBV. Li et al. [33] discussed the infection for both vertical and horizontal transmission of hepatitis B for the host population and the basic reproduction number R 0 for the vertical transmission. Zou et al. [34] proposed that immunization is very necessary to control the HBV infection, especially for those who are at high risk of this infection. Haq et al. [35] focused on vector-borne diseases by the direct method taking the population as constant. Haq et al. [36] studied the approximate solution with the help of Laplace transformation and used the RK-4 method to find the numerical solution of the proposed model. Pang et al. [37] presented the susceptible infectious removed (SIR) epidemic model by considering the effects of vaccination to control the infection of the HBV. Liu et al. [38] discussed the SIR model with a nonlinear incidence rate and used the Lyapunov function to derive global stability. Zaman et al. [39] used the RK-4 method to find the numerical solution for the proposed model for HBV. Zeb et al. [40] presented the smoking model by using the non-standard finite difference method. RK-4 and ODE45 were used to find out the numerical results. Rahman et al. [41] showed global stability at different equilibrium points, and the numerical solutions were shown graphically. Zaman et al. [42] proposed an optimal control technique to study infected individuals, and the dynamical behavior of HBV R 0 was derived. Khan and Zaman [43] presented Lyapunov function theory to derive global stability. The basic reproduction number R 0 was derived via the next-generation matrix method.

The effective technique for examining the dynamics of various diseases as they arise in real-world situations is mathematical modeling. Various extensively utilized models have been implemented by biologists and mathematicians to understand the dynamics of how infectious diseases propagate through populations. Mathematical modeling serves as a highly versatile technique, yielding insightful findings pertinent to the realm of health care. The following points outline the primary aims of the research:

1. Establish a novel mathematical model for HBV with vaccination and treatment rate. The model is solved using higher-order Galerkin Scheme.

2. To find the basic reproduction number R 0 by implementing the next-generation matrix approach.

3. To find the model stability by using the Lyapunov function.

4. To discuss the local and global stability to assess the most sensitive parameters of disease transmission.

5. In order to confirm the validity and accuracy of the proposed scheme, the model is solved using the classical RK and compared the results obtained from both schemes.

2 Proposed mathematical model for HBV infection

In this section, the HBV mathematical model is considered. We divide the system into seven compartments, i.e., susceptible individuals denoted by S ( t ) , exposed individuals termed L ( t ) , people with acute infection identified by A ( t ) , individuals with chronic HBV infection represented by B(t), carriers of HBV infection denoted by C(t), recovered population characterized as R ( t ) , and lastly the vaccinated individuals represented by V ( t ) . The parameters used in the presented model are as follows: Λ represents the birth rate, ψ represents the carrier transmission rate of HBV, γ 1 demonstrates rate of increase from acute to chronic, γ 2 is the rate at which individuals move to recovery from the carrier class, d 0 and d 1 represent the biological and provoked risks of mortality from hepatitis, respectively, while a describes the rate at which people shifts from latent to acute class, η represents the rate of vertically infected individuals, ξ represents the fraction of effectively immunization, and β represents the transmission rate of HBV infection. The graphical representation of the proposed model is given in Figure 1. The proposed model is as follows (Table 1):

and the initial conditions are given as

Figure 1

Pictorial representation of the proposed model.

S ≥ 0 , L ≥ 0 , A ≥ 0 , B ≥ 0 , C ≥ 0 , R ≥ 0 , V ≥ 0 .

3 Continuous Galerkin Petrov method for the HBV model

In this section, we implement the Galerkin scheme, particularly the cGP(2) scheme (see previous studies [44,45,46,47,48,49,50,51,52,53] for detailed information) to tackle numerically the system of ODEs that mathematically characterize the aforementioned model. By implementing the standard form of the ODEs, we implement the Galerkin scheme to approximate the solutions, leading to a thorough examination of the model’s dynamics and behavior. Let us suppose

initially at time t = 0.

For ϱ : K = [ P , Q ] → S & K = [ P , Q ] function for positive t > 0 ,

where d h ρ ( h ) is the time derivative. ρ ( t ) = ( ρ 1 ( t ) , ρ 2 ( t ) , ρ 3 ( t ) , ρ 4 ( t ) , ρ 5 ( t ) , ρ 6 ( t ) , ρ 7 ( t ) ) ∈ S at time t = 0, and ε = ( ε 1 , ε 2 , ε 3 , ε 4 , ε 5 , ε 6 , ε 7 ) is non-linear defined as ε : K ×   S → S .   → S

The formulation of problem of (3.12) is as follows: find ρ ∈ X ′ such that ρ ( 0 ) = ρ 0 and

in which X′ suggests the optimal solutions, Y′ serves for experiment, when F = [ P , Q ] describes its time. Assuming the domain E ( F , S ) = E O ( J , S ) as the space of linear equations ρ : F → S interconnect with the norm to explain the operation t → ρ ( t ) .

Square the integral function Ł 2 ( F   , S ) discontinuous functions manifested as

So, attempting to discretize the Galerkin time, we prevalently partition the period interval J into N subphases. F n = [ t n − 1 , t ] , where n = 1 , 2 , 3 , … N , and 0 = t 0 < t 1 < … . t N − 1 < t n = S . The time linearization component is given the symbol, and so the extreme time will just be illustrated with this sign τ = max 1 ≤ n ≤ N τ n , where τ n = t n − t n − 1 time intervals of J n . Now, we will estimate the continuous answer off a function a τ : F → S , which is piecewise polynomial of a certain order w.r.t. time. Then,

where H l ( F n , S ) = { ρ : F n → S , ρ ( t ) = ∑ j = 0 l U j t j , for all t ∈ F n , U j ∈ S , ∀ f } , and the Y τ k , get

with ρ τ ∈   X ′ τ l , i.e., ρ τ ( 0 ) = 0

As the discrete test space, Y ′ τ l is discontinuous and complicated and can be solved in a time marching operation where successively local problems on the time intervals are resolved. So, we go for the test function

Similarly, for the cGP(k) method:

i.e.,

For estimating ρ τ | t n , replace

where other components of S and another real-valued mechanism comprising their weighting U n j . The Lagrange partial derivatives t n , j ∈ F n meet the conditions having ( l + 1 ) adequate nodal positions.

ϑ n , j ( t n , j ) = Δ i , j for i , = 0 , 1 , 2 , … l . and j = 0 , 1 , 2 , … l .

Δ i , j , Kronecker delta is given by

let t n , 0 = t n − 1 , stands for the equation’s fundamental assumptions, i.e.,

The points t n , 1 , t n , 2 , … , t t n , l are being taken from the GauB ℬ -Lobatto formula at F n interval for d t ρ τ , having

Thus, we have

In a straight forward manner, we arrive at the following

with the basic ϑ n , j ∈ H k ( F n ) function having

Designate basic functions having

Now, defining basic real J n the time interval by

ϑ n , j t = ϑ ˆ j ( t ) with t ˆ = α n − 1 ( t ) = 2 τ n t + t n − 1 − t n 2 ∈ k .

In a similar fashion, we create the most appropriate comparative wavelet coefficients for the sample basic functions   φ n , i , i.e., ψ ˆ i ∈ H l − 1 ( k ) ,

By using the Galerkin method, we get φ ∈ H l − 1 and ∀ ϑ ∈ S ,

This implies that

Here, α μ ^ are the weights and t ˆ μ ∈ [ 1 , − 1 ] are the integration prong with t ˆ 0 = − 1 and t ˆ l = 1

Let J ˆ n be cGP(k), and undetermined U n j ∈ S ∀ j = 1 , 2 , 3 , … , l , thus

For U n 0 = U n − 1 l and n greater than 1 and U 1 0 = ρ 0 for n equal to 1, with z i , j

Once we have broken this system, the initial value is set for the current time interval

cGP (1) Method:

According to GauB ℬ -Lobatto two-point formulae,

t n , 0 = t n − 1 , h n , 1 = t n & ρ 0 ^ = ρ 1 ^ = 1 ,

ψ ˆ 0 ( t ˆ ) = 01 , can have sign U n 1 = ρ τ ( t n ) ∈ S , i.e.,

cGP(2) Method:

Now, the three-point GauB-Lobatto method which is also known as the Simpson rule is used for defining the basic function with weights α 0 ^ = α 2 ^ = 1 / 3 , α 1 ^ = 4 / 3 and t 0 ^ = − 1 , t 1 ^ = 0 , t 2 ^ = 1 .

Then, we have

So, the system is to be solved for the U n 1 , U n 2 ∈ V for the known U n 0 = U n − 1 2 which becomes

4 Comparative analysis of Galerkin and RK schemes

In this section, we conduct a comprehensive comparative analysis of the numerical outcomes obtained using the RK4-method and cGP(2)-method for solving the HIV model. Both numerical schemes are implemented for the model using a computer code written in MATLAB, and simulations are performed over a time span of t ∈ [ 0 , 1 ] . The accuracy of the Galerkin method is assessed by comparing its solutions with those obtained using the RK4 method at various time points within the specified interval. The results are visualized in Figures 2–8, which depict the temporal evolution of susceptible, latent, acutely infected, chronically infected , carrier, recovered, and vaccinated populations, respectively. A thorough examination of the graphical output reveals that the Galerkin method exhibits excellent agreement with the RK method, demonstrating its reliability and accuracy in approximating solutions to the HIV model. Furthermore, the proposed approach is shown to be a trustworthy tool for finding approximate solutions to real-world problems, particularly in the context of mathematical modeling of infectious diseases.

Figure 2

Comparison of RK4 and cGP(2) results for susceptible population S ( t ) .

Figure 3

Comparison of RK4 and cGP(2) results for latent population L ( t ) .

Figure 4

Comparison of RK4 and cGP(2) results for acutely infected individuals A ( t ) .

Figure 5

Comparison of RK4 and cGP(2) results for chronically infected individuals B ( t ) .

Figure 6

Comparison of RK4 and cGP(2) results for carrier-infected individuals C ( t ) .

Figure 7

Comparison of RK4 and cGP(2) results for recovered population R ( t ) .

Figure 8

Comparison of RK4 and cGP(2) results for vaccinated population V ( t ) .

Furthermore, we investigated the influence of parameter ξ on the dynamical behavior of the model, which describes the temporal evolution of the susceptible, latent, acutely infected, chronically infected, carrier, recovered, and vaccinated populations. The results are presented graphically in Figures 9–15, respectively. These figures illustrate how the parameter ξ affects the dynamics of each population, providing valuable insights into the impact of this parameter on the spread of the disease and the effectiveness of vaccination strategies.

Figure 9

Dependence of susceptible population as a function of varying “ ξ ”.

Figure 10

Dependence of latent population as a function of varying “ ξ ”.

Figure 11

Dependence of acutely infected population as a function of varying “ ξ ”.

Figure 12

Dependence of chronically infected population as a function of varying “ ξ ”.

Figure 13

Dependence of carrier infected population as a function of varying “ ξ ”.

Figure 14

Dependence of recovered population as a function of varying “ ξ ”.

Figure 15

Dependence of vaccinated population as a function of varying “ ξ ”.

Furthermore, the analysis results in Figures 9–15 demonstrate a unifying pattern, indicating that the concentrations of susceptible, latent, acutely infected, chronically infected, carrier, recovered, and vaccinated individuals all exhibit a growing trend corresponding to an increase in the value of parameter ξ . In particular:

1. The population of susceptible (Figure 9) demonstrates a consistent increase in concentration as the parameter ξ escalates, signifying a greater number of individuals vulnerable to infection.

2. The latent population (Figure 10) demonstrates an increasing concentration as the parameter ξ increases, indicating a larger number of individuals in the latent phase of infection.

3. The population experiencing acute infection (Figure 11) exhibits a significant rise in concentration as the parameter ξ increases, suggesting a greater number of individuals affected by acute infection.

4. The population with chronic infections (Figure 12) exhibits a comparable trend, where an increase in parameter ξ corresponds to a greater concentration of individuals affected by chronic infections.

5. The carrier population (Figure 13) demonstrates an increase in concentration as the parameter ξ increases, indicating a larger number of individuals harboring the infection.

6. The data presented in Figure 14 illustrate a notable increase in concentration corresponding to the increase of parameter ξ , suggesting a greater number of individuals successfully recovering from infection.

7. The data regarding the vaccinated population (Figure 15) reveal a significant increase in concentration as the parameter ξ increases, indicating a greater number of individuals receiving vaccination.

The empirical findings indicate that parameter ξ significantly influences the dynamics of the model, as higher values correspond to an increased number of individuals within each category.