Optimized Computational Schemes for Finding Multiple Roots in Nonlinear Biomedical Engineering Problems

The capacity to solve nonlinear equations is an essential notion in computer science, engineering, and applied mathematics. In several real-world applications, the nonlinear equations are utilized to model various systems that are influenced by biology, chemistry, and physics [1], [2], [3]. Some of the most computationally intensive subcategories of these problemes are several-rooted nonlinear equations, which present major challenges for numerical analysis [4], convergence, and stability. Classical root-finding schemes will be ineffectual or, in some cases, unsuccessful, highlighting the critical need for increasingly complex numerical methodologies to solve

Nonlinear equations having multiple roots having multiplicity “$\sigma$” appear in a variety of scientific disciplines, including biomedical engineering [5], control systems [6], chemical kinetics [7], and epidemiology [8]. The system being described usually has numerous steady states or equilibrium points that these roots correspond to. For instance, a biological system may have multiple physiological states, while an engineering system may have multiple operational modes. Accurately identifying each of these roots is crucial to obtaining the complete behavior of such systems as well as to creating trustworthy simulations and forecasts. In biomedical models like glucose-insulin dynamics, enzyme kinetics, and disease propagation models, this is a particularly important example [9]. Multiple roots may represent disease-free and endemic equilibria (in epidemiological models) or, in glucose regulation models, normoglycemia, hyperglycemia, and hypoglycemia. It is therefore very useful to develop efficient numerical techniques to address such nonlinear situations.

Several iterative approaches for solving nonlinear equations have been developed over the years. It has been common practice to employ classical techniques like Newton-Raphson [10], Secant [11], and Bisection [12] for single roots. Their convergence characteristics, however, worsen significantly when expanded to multiple roots. To address this drawback, various scholars have offered new approaches and modifications. Chun and Ham [13], for example, developed modified Newton methods based on multiplicity data in 2006. Sharma and Gupta [14] in 2014 proposed higher-order iterative algorithms, specifically those with multiple roots. Fried [15] in 2016 proposed a chord iterative method for roots of uncertain multiplicity. Based on the weight function procedure, Chicharro et al. [16] in 2020 proposed a family of multiple root-finding schemes. Behl et al. [17] in 2023, proposed two-step optimal multiple root-finding methods and solved some problems in physics and chemistry. There are other scholars who used different techniques to find the multiple roots using various strategies see, e.g., and the references cited therein [18], [19], [20], [21], .

The literature review highlights the following significant gaps:

• The majority of the literature's iterative strategies for multiple roots are based on an explicit awareness of root multiplicity.

• In the absence of a priori multiplicity information, there aren't many high-efficiency optimal-order methods for numerous roots.

• Only a few articles have carefully explored the dynamical properties of such approaches on the complex plane, particularly the basins of attraction.

• Such advanced techniques are rarely applied to real-world nonlinear biomedical engineering models.

In this paper, we present a novel two-step optimal iterative approach for approximating nonlinear equations with multiple roots. The Kung-Traub conjecture states that the approach is optimal, and it does not require prior knowledge of the root's multiplicity. The important contributions of this study are:

• A two-step iterative approach with optimal order of convergence for multiple roots is constructed with the fewest number of function evaluations.

• The stability and convergence behavior of the suggested strategy are visualized through a thorough dynamical analysis that makes use of basins of attraction.

• Utilizing the approach to real biomedical models that produce multi-root nonlinear equations, such as models for disease, enzyme kinetics, and glucose regulation.

• Comparative analysis with existing approaches demonstrates increased accuracy, convergence rate, and robustness.

The primary goal of this study is to advance iterative methods by developing robust, high-efficiency algorithms for real-world, high-stakes biomedical modeling problems involving nonlinear equations with multiple roots where higher accuracy is essential. The remainder of the paper is structured as follows:

After the introductry section, Section 2 presents the mathematical formulation of the proposed two-step optimal strategy, convergence theory, and error analysis after the introductory sections. Section 3 discusses several biomedical engineering applications, developing nonlinear models from enzyme kinetics, glucose regulation, and epidemiological models and demonstrating how they have multiple roots. The accuracy, stability, and efficiency of the approach are evaluated through extensive numerical experiments on two real-world biomedical problems in Section 3. Section 4 concludes the study by presenting findings and proposed future research topics.

Iterative techniques are important for precisely locating the roots of many nonlinear equations, which are prevalent in challenging scientific and engineering problems. They increase computational efficiency by gradually improving estimation with repeated calculations. They offer reliable approaches for determining the equilibrium points required to describe physiological processes in biomedical models. They are adaptable and flexible techniques that are useful in accurately and consistently resolving real-world nonlinear systems. Here, we discuss some well-known existing methods in the literature. Li et al. [22] proposed the following optimal scehems for finding multiple root of Eq. (1) as

where,

$\phi {11}=-\frac{1}{2}\sigma \left( \sigma -2\right) ,\phi {12}=-\frac{1}{\sigma },\phi _{13}=\frac{1}{\sigma }\left( \frac{2+\sigma }{\sigma }\right) ^{\sigma }.$

Zhou et al. [23] proposed the following schemes for finding multiple roots as

where,

$\phi _{21}=\left( \frac{\sigma +2}{\sigma }\right) ^{\sigma }\frac{f^{\prime }(y^{\left[ t\right] })}{f^{\prime }(x^{\left[ t\right] })};\phi _{22}=- \frac{\sigma \left( \sigma +2\right) ^{3}}{8}\left( \frac{\sigma }{\sigma +2} \right) ^{\sigma };\phi _{23}=\sigma ^{3}+3\sigma ^{2}+2\sigma -4.$

Liu and Zhou [24] proposed the following optimal two fmailes for solving multiple roots of nonlinear equations as

and

Soleymani et al. [25] proposed the following family of fourth order methods as

where,

$\phi _{31}=\frac{1}{16}\sigma ^{3-\sigma }\left( 2+\sigma \right) ^{\sigma },$

$\phi _{32}=\frac{8-\sigma \left( 2+\sigma \right) \left( -2+\sigma ^{2}\right) }{8\sigma },$

$\phi _{33}=\frac{1}{16}\left( -2+\sigma \right) \sigma ^{-1+\sigma }\left( 2+\sigma \right) ^{3-\sigma }.$

Here, we proposed the following family of fourth-order method (MSN$^{[\ast ]}$) for finding multiple roots of nonlinear equations

where, $\beta \in \mathbb{R}.$

Convergence analysis is required for numerical systems that solve nonlinear equations to ensure reliability and efficiency. This is done by determining the conditions under which the approach converges to the solution. Besides, convergence analysis can be considered a basis for evaluating the methods and adjusting them to address specific problem types or restrictions. The local convergence of the developed method is discussed in Theorem 1.

Theorem 1. Assume that $f: I\rightarrow R$ is a sufficiently differentiable function defined on an open interval, and that $\zeta $ is a root of multiplicity $\sigma $. If the initial guess is close enough to the exact root, the two-step iterative technique

will converge with order four, satisfying the error relation

where, $d_{j}=\frac{\sigma !}{j!}\frac{f^{\left( i\right) }\left( \zeta\right) }{f^{\left( \sigma \right) }\left( \zeta \right) };j\geq 1.$

Proof: Let the simple solution $\zeta $ of f and $e^{\left[ t\right] }=x^{\left[ t\right] }-\zeta ,e^{\left[ t+1\right] }=z-\zeta $. Expanding $f(x^{\left[ t\right] })$ and $f^{\prime }(x^{\left[ t\right] })$ at $x=$ $\zeta$ utilizing Taylor series expansion, we yeild the following

Dividing $f(x^{\left[ t\right] })$ by $f^{\prime }(x^{\left[ t\right] }),$ we have

where,

$\tau {11}=\frac{d{1}}{\sigma }\left( \frac{d_{1}^{2}}{\sigma }+2d_{1}^{2}+\sigma d_{1}^{2}-2d_{2}-\sigma d_{2}\right),$

$\tau {12}=-d{1}^{3}-\frac{d_{1}^{3}}{\sigma ^{3}}-\frac{3d_{1}^{3}}{\sigma^{2}}-\frac{3d_{1}^{3}}{\sigma }+2d_{1}d_{2}+\frac{4d_{1}d_{2}}{\sigma ^{2}}+\frac{6d_{1}d_{2}}{\sigma }-d_{3}-\frac{3d_{3}}{\sigma }.$

Using the expression provided in the first step of Eq. (7), we obtain the following result in Maple 18:

Expanding $f(y^{\left[ t\right] })$ at $y=$ $\zeta $ utilizing Taylor series expansion, we have

Dividing $f\left( y^{\left[ t\right] }\right) $ by $f(x^{\left[ t\right] })$ and computing the term

where,

$\tau _{21}=\frac{\left( -9+2\beta -\sigma +2\sigma ^{2}\right)d_{1}^{3}+\left( 2\sigma -6\sigma ^{2}\right) d_{1}d_{2}+6\sigma ^{2}d_{3}}{2\sigma ^{2}}.$

Using $\left( \frac{f\left( y^{\left[ t\right] }\right) }{f\left( x^{\left[ t\right] }\right) }\right) ^{\frac{1}{\sigma }}$ to compute error of

Now using Eq. (16) and Eq. (17), in the last-setp of the iterative scheme to obtaine its order of convergence as:

Hence prove the theorem.

Numerical solutions for fractional nonlinear equations help validate theoretical models in practice. They provide an approximation when analytical solutions are hard or impossible to obtain. This would enable the study of complex systems with fractional dynamics, providing insight into real-world applications in physics, engineering, operations research, and finance. Furthermore, numerical findings can be used to compare and evaluate different iteration techniques to improve convergence behavior efficiency. Providing more control over modeling and support in evaluating the potential effects of fractional parameters on system behavior are among the other advantages. Given these, a numerical approach is required for solving fractional nonlinear equations because obtaining exact solutions is often inapplicable. This section provides four practical examples from thermodynamics, computational chemistry, engineering control systems, and operations research.

where, $e^{[t]}$ absolute error. To check the efficiency and stability, we compare our scheme with the following well known numerical schemes exists in the literature as:

1-Zafar et al. [26] method (FZM$^{[\ast ]}$) of convergence order four

where, $\alpha \in \mathbb{R}.$

2-Shengguo et al. [27] method (SHM$^{[\ast ]}$) of fourth order convergence

3-Sharma and Sharma [28] method (SSM$^{[\ast ]}$) of convergence order four for finding multiple roots of Eq. (1) as given as

where,

$\vartheta _{11}=\sigma ^{4}-4\sigma ^{2}+8\sigma,$

$\vartheta _{12}=\left( \sigma +2\right) ^{2}\left( \frac{\sigma }{\sigma +2}\right) ^{\sigma },$

$\vartheta _{13}=2\left( \sigma -1\right) -\left( \left( \sigma +2\right)\left( \frac{\sigma }{\sigma +2}\right) ^{\sigma }\frac{f^{\prime }(y^{\left[t\right] })}{f^{\prime }(x^{\left[ t\right] })}\right).$

4-Zhou et al. [29] method (ZHM$^{[\ast ]}$) for finding multiple roots with known multiplicity as given

where,

$\vartheta _{21}=\sigma ^{3}\left( \frac{\sigma +2}{\sigma }\right) ^{2\sigma},$

$\vartheta _{22}=-2\sigma ^{2}\left( \sigma +3\right) \left( \frac{\sigma^{4}+6\sigma ^{3}+8\sigma ^{2}+8\sigma }{\sigma }\right).$

The order of convergence of the method Eq. (23) is four. Using the flowchart Figure 1 and Algorithm 1, we compute the multiple roots of nonlinear equations in Maple 2020 using a 64-bit operating system and a 128-digit floating-point algorithm.

Example 1. Glucose-Insulin Regulatory Model [30]

Blood glucose regulation is a straightforward physiological mechanism that involves interactions between insulin and glucose. Diabetes patients have reduced regulation of this process, therefore understanding the underlying dynamics is crucial for management and treatment. Insulin and glucose interactions are explained bymathematical modeling, especially in various physiological and pathological conditions.

Principles of Basic Physiology

• Food-derived glucose enters the bloodstream and is absorbed by cells, generating energy.

• Insulin, secreted by the pancreas, promotes glucose absorption by cells and reduces blood glucose levels.

• The link is nonlinear and involves feedback processes in which glucose levels influence insulin secretion and vice versa.

1-Determine the State Variable: Let $x(t)$, glucose levels in blood at time t, and $\zeta =\sqrt{\pi }$ be the desired (normal) glucose level (normalized for modeling purposes).

2-Feedback Control of Insulin and Glucose: Consider a cubic feedback model for the insulin respon.

A large deviation from normal glucose levels (too high or too low) results in a nonlinear biological feedback as described above. It depicts both little and major deviations, with less input on the former.

3-Regulation by Enzymes or Hormonal Saturation: Some reactions (for example, insulin activity and glucose uptake by cells) become saturated at increasing levels. A logarithmic term is commonly used to model such behavior:

• The slower logarithmic growth occurs as x grows.

• The shift $+1$ prevents undefinable behavior when $x=0$ and guarantees that the argument is positive.

Modeling this as a resistive term (opposing the feedback) and scaling it by a tiny factor (e.g., 0.5) yields:

Consequently, the general behavior is:

The nonlinear equation Eq. (27) has exact root is 0.5210797742 with multiplicity 3. The numerical outcomes are presented in Table 1.

Table 2 clearly shows that in terms of percentage convergence, which shows that the stability of the schemes for solving biomedical problems, our method behaves much better than existing methods and consumes less computational time, demonstrating that it is more efficient and reliable on supercomputers to solve large-scale problems. The approximate order of convergence matched the theoretical order of convergence and outperformed other existing approaches for an initial starting value of 1.012.

Analysis of Basins of Attraction [31]

The basins of attraction depicted in this study provide information regarding the proposed method's stability and robustness. Areas of attraction that are broadly and evenly spread exhibit increased tolerance and robustness to higher initial estimates. Compared to alternative fractional-based approaches, the suggested system performs better in:

• The suggested algorithm features wider basins, indicating superior global convergence characteristics.

• Root stability is guaranteed by basin symmetry, particularly in multi-root designs.

• Reduced sensitivity to perturbations is indicated by smooth boundary curves in attraction zones, which suggests numerical stability.

To demonstrate this, we calculated and reported the percentage of initial points that converged to the exact roots in Table 2. Table 2 displays the total number of points (64000), the computed time (in seconds), the number of iterations, and the convergence percentage (Per-convergence). The created technique shows consistency with the scheme's CPU-time stability for solving nonlinear equations, however it requires less iterations and has a bigger convergence region when compared to the existing methods. Figure 2 illustrates our method's basins of attraction in comparison to the existing approaches. It is evident that, in comparison to FZM$^{[\ast ]},$ SHM$%^{[\ast ]},$ SSM$^{[\ast ]},$ and ZHM$^{[\ast ]}$, our newly developed method MSN$^{[\ast ]}$ convergence behavior and region are superior.

Our method performs better than existing methods in dynamical convergence analysis and residual error on both stopping criteria, as shown in Table 1 and Table 2.

Physical significance: We conclude the following from the whole discuss that:

• A model with cubic feedback for glucose variations from the norm level.

• A logarithmic model of resistance to glucose uptake or insulin action.

• The point of balance in the glucose-insulin system can be found by solving Eq. (27), and adding them together represents equilibrium.

• In terms of conistancy and stability our method perform better as compared to exsitng scheme for finidng multiple roots of Eq. (27).

Example 2. Pharmacokinetics-Biomedical Engineering Problem

Pharmacokinetics (PK) is a discipline of pharmacology that investigates how the body absorbs, distributes, metabolizes, and excretes drugs over time. To put it simply, it provides a response to:"What does the body do to the drug?"

PK is essential in:

• Drug development and approval.

• Dosing regimen design.

• Personalized healthcare.

A simple one-compartment model with first-order elimination looks like this:

where, $x\left( t\right)$ is drug concentration at time $t$, $v_{\sigma ax}$ is enzymatic saturation and $k$ is the elimination rate constant. This yields nonlinear differential equations, which, under steady-state assumptions, can be reduced to nonlinear algebraic equations with various roots indicating different system equilibrium states.

where,

$x$: Concentration of drug.

$x^{2}\left( x-2\right) ^{2}$: Nonlinear dosage response with thresholds.

$e^{-x}$: Natural decay or removal.

The solution to the problem of illustrating steady states or thresholds (such as toxicity and ineffectiveness).

The numerical outcomes are presented in Table 1.

Table 3 clearly shows that in terms of cumulative residual error and functional value error, which shows that the stability of the schemes for solving biomedical problems, our method behaves much better than existing methods and consumes less computational time, demonstrating that it is more efficient and reliable on supercomputers to solve large-scale problems. The approximate order of convergence matched the theoretical order of convergence and outperformed other existing approaches for an initial starting value of 1.012.

Analysis of Basins of Attraction

To demonstrate this, we calculated and reported the percentage of initial points that converged to the exact roots in Table 4. Table 4 displays the total number of points (64000), the computed time (in seconds), the number of iterations, and the convergence percentage (per convergence). The created technique shows consistency with the scheme's CPU-time stability for solving nonlinear equations; however, it requires fewer iterations and has a bigger convergence region when compared to the existing methods. Figure 3 illustrates our method's basins of attraction in comparison to the existing approaches. It is evident that, in comparison to FZM$^{[\ast ]},$ SHM$%^{[\ast ]},$ SSM$^{[\ast ]},$ and ZHM$^{[\ast ]}$, our newly developed method MSN$^{[\ast ]}$ convergence behavior and region are superior.

It is evident from Table 3 and Table 4 that our approach outperforms current approaches in terms of residual error on both stopping criteria and our methods shows better convergence behaviour as compared to existing methods.

Numerical simulations using iterative approaches are important for confirming theoretical findings and determining the practical effectiveness of suggested algorithms. When applied to real-world issues, they offer useful information on convergence properties, stability, and accuracy, particularly in intricate domains like biomedical engineering.

• In addition to its higher order convergence rate, the method can approximate nonlinear equation solutions more accurately, especially in complex, memory-based systems (see Table 1, Table 2, Table 3 and Table 4).

• Numerical results show that the proposed method significantly reduces residual error, computational order of convergence, processing time, and percentage convergence (Figure 1 and Figure 2) when compared to conventional methods, making it more suitable for real-world applications.

• Although standard strategies fail to address fractional derivative problems, our method outperforms them in fields such as biomedical engineering.

• The method has a wide range of applications, from technical control systems to operations research, making it a versatile solution for nonlinear optimization issues.

Despite their effectiveness, iterative techniques for locating multiple roots of nonlinear equations have limitations. Their convergence frequently hinges on the precision of the initial approximation and understanding of the multiplicity of the root, which are not always predictable in practical situations. In addition, if not carefully designed, these approaches may exhibit poor convergence or ill-conditioning at numerous roots. In some circumstances, the approaches will also be sensitive to round-off or diverge when used to very nonlinear or ill-conditioned problems.

In this study, we provide a new two-step iterative method designed specifically for rapidly determining multiple roots of nonlinear equations. Thorough theoretical convergence analysis confirms that the proposed method has better convergence behavior. To demonstrate its efficacy and use, we conducted extensive dynamical evaluations on a wide range of tough biomedical engineering problems see Figure 1, Figure 2 and Figure 3. The numerical results show the accuracy, speed of convergence, and stability of our method MSN$^{[*]}$ compared to the state-of-the-art approaches FZM$^{[\ast ]},$ SHM$^{[\ast ]},$ SSM$^{[\ast ]},$ and ZHM$^{[\ast ]}$, in the literature. In practical real-world applications, the results (Table 1, Table 2, Table 3 and Table 4) validate that the suggested method, MSN$^{[*]}$ provides a robust and effective substitute for multirooting nonlinear equations.

Future research can go in various promising directions. First, the suggested approach would become much more applicable to more intricate biological and engineering issues if it were applied to fractional-order models and systems of nonlinear equations. Second, it might be more computationally efficient if adaptive step-size or memory-based acceleration techniques are incorporated. Finally, using the approach for real-time biomedical data and uncertainty quantification would provide useful information on its robustness under realistic conditions.