An accelerated double-step-length iterative method for nonlinear systems of equations with applications to motion control

1. Introduction

Nonlinear systems of equations (SNE) are fundamental in a wide range of fields, from scientific research to various engineering applications. These systems are distinguished by their complexity and by the intricate relationships among variables, which often necessitate the use of advanced solution techniques. Consequently, researchers are motivated to develop effective iterative methods to address the challenges presented by these nonlinear equations [2– 6]. The primary aim of these methods is to achieve accurate solutions while minimizing the computational resources required for their determination.

A nonlinear system of equations is typically represented in the following mathematical form: \[\label{a1} \mathfrak{F}(s)=0. \tag{1}\] In this expression, \(\mathfrak{F}:\mathbb{R}^n\to\mathbb{R}^n\) denotes a nonlinear mapping that connects inputs from the \(n\)-dimensional real space to their corresponding output values. In this paper, we refer to \(\mathbb{R}^n\) as the \(n\)-dimensional real space, and we use the notation \(\|\cdot\|\) to denote the Euclidean norm, which serves as a standard measure of distance in this context. By investigating these systems and the methods for solving them, we aim to contribute meaningful insights to the field.

To solve problem (1), various iterative techniques can be employed. These techniques include double-direction methods [7– 10], double-step-length methods [11– 16], Newton’s method (NM) [2], and quasi-Newton methods (QNM) [6, 17]. Among these techniques, NM stands out due to its favorable characteristics, such as ease of implementation and rapid convergence. However, the Jacobian matrix (JM) must be computed and stored during each iteration. These techniques use the recursive formula to generate a sequence of points: \[\label{h3} s_{k+1}=s_k+\alpha_k q_k,\qquad k=0,1,\ldots, \tag{2}\] where \(\alpha_k\) is the step size, \(q_k\) is the search direction, and the current and next iterates are \(s_k\) and \(s_{k+1}\), respectively. The NM search direction \(q_k\) is determined by \[\label{ND} q_k=-(\mathfrak{F}^{\prime}(s_k))^{-1}\mathfrak{F}(s_k).\] Here, \(\mathfrak{F}^{\prime}(s_k)\) is the JM of \(\mathfrak{F}\) at \(s_k\). However, one limitation of Newton’s method is that it requires the derivative, namely the JM, to be computed at each iteration, which may not always be available or may not be precisely determined. In such cases, NM cannot be applied directly.

The global optimization problem can be derived from (1) [1]. Let \(\mathfrak{h}\) be a merit function given by \[\label{a7} \mathfrak{h}(s)=\frac{1}{2}\|\mathfrak{F}(s)\|^{2}. \tag{3}\] Consequently, problem (1) can be equivalently expressed as the following global optimization problem: \[\label{global} \min \mathfrak{h}(s),\qquad s\in\mathbb{R}^n,\qquad \mathfrak{h}:\mathbb{R}^n\to\mathbb{R}. \tag{4}\]

In order to achieve global convergence, it is essential to employ a step length in iterative methods. This is because the use of a step length reduces the function value, thereby helping to ensure global convergence. For instance, NM is locally convergent. However, NM becomes globally convergent when equipped with a step length. In addition, the step length can be computed using any suitable line-search technique. The exact line search is appropriate for use in practice [18] and satisfies \[\label{ab1} \mathfrak{h}(s_k+\alpha_k q_k) = \min_{\alpha>0}\mathfrak{h}(s_k+\alpha q_k). \tag{5}\]

In practical applications, determining the precise value of the step length is often challenging. Consequently, inexact line-search methods are the predominant approach used in practice [16, 19, 20]. In this regard, the authors in [21] introduced a specific line-search rule that effectively identifies an appropriate step length: \[\label{bb1} \mathfrak{h}(s_k+\alpha_k q_k)-\mathfrak{h}(s_k) \leq \sigma\alpha_k\nabla \mathfrak{h}(s_k)^Tq_k, \qquad \sigma\in(0,1). \tag{6}\]

The technique outlined in (6) indicates that the JM must be calculated at every iteration. This requirement has the potential to substantially increase computational demands, especially in the context of large-scale problems or when dealing with complex matrix calculations. To address these concerns, Yuan\(^*\) and Lu [17] proposed a new backtracking inexact technique for determining the step length \(\alpha_k\) as follows: \[\label{bb2} \|\mathfrak{F}(s_k+\alpha_k q_k)\|^2 \leq \|\mathfrak{F}(s_k)\|^2+\phi\alpha_k^2\mathfrak{F}(s_k)^Tq_k,\] where \(\phi\in(0,1)\). Another inexact line search was introduced in [1], and the step length can be calculated as follows.

Let \(r\in(0,1)\), \(\omega_1,\omega_2>0\), and let \(\{\chi_k\}\) be a given positive sequence such that \[\label{ab29} \sum_{k=0}^{\infty}\chi_k<\chi<\infty.\] Then the step length is chosen to satisfy \[\label{a30} \mathfrak{h}(s_k+\alpha_k q_k)-\mathfrak{h}(s_k) \leq -\omega_1\|\alpha_k\mathfrak{F}(s_k)\|^2 -\omega_2\|\alpha_k q_k\|^2 +\chi_k\mathfrak{h}(s_k).\]

Despite the favorable properties of NM, it requires the computation of the Jacobian or its approximation at each iteration. This requirement makes it less suitable for large-scale problems, as it necessitates extensive matrix storage, which can be costly in numerical experiments. To address these concerns, matrix-free approaches have been proposed. The double-direction method (DDM) is one such successful matrix-free approach [10], which generates the iterates using the following formula: \[\label{dd} s_{k+1}=s_k+\alpha_k q_k+\alpha_k^2 b_k, \tag{7}\] where \(q_k\) and \(b_k\) are the respective search directions. The essential idea underlying the DDM is that it incorporates two adjustments in the system represented by (7). If one correction fails, the other serves as a backup to correct the system.

In their work, Petrovi\(\acute{c}\) and Stanimirovi\(\acute{c}\) [10] made significant contributions to the field of unconstrained optimization by proposing a DDM. This innovative approach utilizes an acceleration parameter \(\delta_k>0\) to construct an approximation of the Hessian matrix, which is expressed as \[\nabla^2\mathfrak{h}(s_k)\approx \delta_k I,\] where \(I\) denotes the identity matrix and the acceleration parameter is given by \[\delta_{k+1}^{ADD} = 2 \frac{ \mathfrak{h}(s_{k+1})-\mathfrak{h}(s_k) -\alpha_k\nabla\mathfrak{h}(s_k)^T \bigl(\alpha_kq_k-\delta_k^{-1}\nabla\mathfrak{h}(s_k)\bigr) }{ \bigl(\alpha_kq_k-\delta_k^{-1}\nabla\mathfrak{h}(s_k)\bigr)^T \bigl(\alpha_kq_k-\delta_k^{-1}\nabla\mathfrak{h}(s_k)\bigr) }.\]

To enhance the efficiency of the DDM presented in [10], Petrovi\(\acute{c}\) further refined it in [16] by introducing a double-step-length scheme (DSM). This was done by replacing \(\alpha_k^2\) in (7) with \(\beta_k\), which is defined mathematically as \[\label{dsl} s_{k+1}=s_k+\alpha_k q_k+\beta_k b_k. \tag{8}\] In the above equation, \(\alpha_k\) and \(\beta_k\) represent two independently chosen step sizes, allowing for a more flexible and adaptive approach to convergence. In addition, the proposed acceleration parameter is given by \[\delta_{k+1}^{ADSS} = 2 \frac{ \mathfrak{h}(s_{k+1})-\mathfrak{h}(s_k) + (\alpha_k\delta_k^{-1}+\beta_k)\|\nabla\mathfrak{h}(s_k)\|^2 }{ (\alpha_k\delta_k^{-1}+\beta_k)^2\|\nabla\mathfrak{h}(s_k)\|^2 }.\]

The numerical results in [16] illustrate that this refined method demonstrates superior performance compared with the original DDM, achieving a reduction in the number of iterations required as well as a decrease in processor time.

The existing literature on derivative-free DDM for solving SNE is notably limited. This lack of comprehensive research inspired Halilu and Waziri [8] to develop a derivative-free method that employs a DDM for addressing (1), as described in Eq. (7). The approach proposed an acceleration parameter \(\delta_k>0\) to create an approximation of the Jacobian matrix, expressed as \[\label{AC} \mathfrak{F}^{\prime}(s_k)\approx \delta_k I. \tag{9}\] The acceleration parameter is given by \[\delta_{k+1}^{IDFDD} = \frac{v_k^Tv_k} {(\alpha_k+\alpha_k^2\delta_k)v_k^Tq_k},\] where \[v_k=\mathfrak{F}(s_{k+1})-\mathfrak{F}(s_k).\] Subsequently, in [22], a modified version of the algorithm in [8] was developed. They established that the proposed method exhibits global convergence, contingent on the condition that the JM is both positive definite and bounded. This finding is significant, as it lays the groundwork for reliable solutions to nonlinear equations without requiring derivatives.

Building upon this framework, Abdullahi et al. [7] made substantial improvements to the performance of the DDM. They modified the initial methodology from [8] by integrating a conjugate-gradient approach for symmetric nonlinear equations. Their enhancement not only retained the derivative-free nature of the method but also facilitated global convergence, achieved through the innovative line-search strategy proposed by Li and Fukushima [1]. This advancement marks notable progress in the efficiency and applicability of DDM for solving SNE.

In another contribution, Halilu and Waziri [9] further refined the DDM discussed in (7) for solving SNE. Their latest method demonstrates global convergence under some mild assumptions, which broadens its applicability. The research includes extensive numerical experiments that validate the efficiency and effectiveness of their proposed method, showcasing promising results that contribute to the ongoing development of robust techniques for solving SNE in a derivative-free manner.

In their insightful study presented in [13], the authors effectively built upon the foundational concepts in Eq. (8) by introducing the DSM. Here, the acceleration parameter is given by \[\label{a16} \delta_{k+1} = \frac{v_k^Tv_k} {(\alpha_k+\beta_k)v_k^Tq_k}.\] This innovative method aims to enhance the efficacy of the DDM described in [8]. The results from their research are particularly encouraging, as they indicate that the DSM leads to significantly faster convergence than the previous method. These findings not only highlight the potential of their approach to optimize numerical algorithms but also pave the way for further advancements in this area of study.

It has been observed that the DSM for unconstrained optimization problems [16] and the DSM for SNE [13] use a line-search approach to determine the values of the two step lengths. Motivated by these concepts, our objective is to develop an improved DSM for solving SNE, where the two step lengths can be computed using a different line-search approach in order to improve both the efficiency and accuracy of the DDM in finding solutions to Problem (1).

The structure of the paper is as follows. The derivation of the suggested approach is explained in the next section. A convergence analysis of the suggested algorithm is given in §3. A number of numerical tests and motion-control applications of the suggested method are presented in §4. The article is concluded in §5.

2. Motivation and the Proposed Method

In this section, we present a derivative-free approach for solving SNE using an acceleration parameter. We begin by recalling the derivative-free spectral method (DSM) introduced in [13], where the iterative sequence \(\{s_k\}\) is generated by \[\label{schm} s_{k+1}=s_k+\alpha_k q_k+\beta_k q_k, \tag{10}\] and the search direction is defined as \[\label{ash18-old} q_k=-\delta_k^{-1}\mathfrak{F}(s_k). \tag{11}\]

Subsequently, in [15], a double-step-length strategy was proposed by defining two search directions in the double-step-length iteration (7) as \[\label{ash18a} q_k=-\eta\delta_k^{-1}\mathfrak{F}(s_k), \tag{12}\] and \[\label{ash18b} b_k=-\mathfrak{F}(s_k), \tag{13}\] where \(\eta\in(1,2)\) is a correction parameter obtained by hybridizing the double-step-length iteration with the Picard–Mann process [23].

Although the DSM proposed in [13] is fully derivative-free and computationally efficient, its reliance on a single acceleration parameter \(\delta_k\) may limit its effectiveness. In particular, when \(\delta_k=1\), the method reduces to a steepest-descent-type iteration, which may result in slow convergence. The double-step-length strategy introduced in [15] incorporates a correction parameter \(\eta\) through hybridization with the Picard–Mann process, leading to improved numerical performance. However, the method employs two distinct search directions, one of which coincides with the classical steepest-descent direction, which is known to converge slowly in certain situations. This may affect the consistency of the update scheme and prevent the method from fully exploiting the acceleration effect of the correction parameter across all directions.

Motivated by these observations, we propose a simplified accelerated derivative-free scheme by reducing the two directions to a single search direction: \[\label{ash18} q_k=-\eta\delta_k^{-1}\mathfrak{F}(s_k). \tag{14}\] This choice integrates both the acceleration and correction parameters into a single direction. This reduction addresses the identified gaps by:

1. maintaining uniform acceleration and correction through the parameter \(\eta\) in all iterations;

2. avoiding the slow convergence associated with the steepest-descent direction;

3. reducing computational cost by eliminating the need to compute two separate directions.

To further reduce computational effort, we assume a relationship between the step lengths \(\alpha_k\) and \(\beta_k\) as \[\label{betak} \beta_k=\frac{\alpha_k}{\eta}. \tag{15}\]

Next, we update the acceleration parameter \(\delta_{k+1}\) using the Taylor series expansion below: \[\label{a23a} \mathfrak{F}(s_{k+1}) \approx \mathfrak{F}(s_k)+\mathfrak{F}^{\prime}(\xi)(s_{k+1}-s_k), \tag{16}\] where \(\xi\in[s_k,s_{k+1}]\) can be expressed as \[\xi=s_k+\varrho(s_{k+1}-s_k),\qquad 0\leq\varrho\leq1.\] By selecting \(\varrho=1\), we have \[\label{a23ab} \mathfrak{F}(s_{k+1}) \approx \mathfrak{F}(s_k)+\mathfrak{F}^{\prime}(s_{k+1})(s_{k+1}-s_k). \tag{17}\]

To make the computations more efficient, we use a diagonal matrix to approximate the JM. This approximation eliminates the need for complex computations and storage associated with the full JM, thereby simplifying the overall procedure. Consequently, we use the approximation \[\label{a25a} \mathfrak{F}^{\prime}(s_{k+1})\approx \delta_{k+1}I. \tag{18}\]

Thus, based on (17) and (18), it is straightforward to verify the following secant equation: \[\label{a23ac} v_k\approx \delta_{k+1}u_k, \tag{19}\] where \[u_k=s_{k+1}-s_k = \left(\alpha_k+\frac{\alpha_k}{\eta}\right)q_k\] and \[v_k=\mathfrak{F}(s_{k+1})-\mathfrak{F}(s_k).\] By multiplying both sides of (19) by the transpose of \(v_k\), we obtain the following update for the acceleration parameter: \[\label{AP} \delta_{k+1} \approx \frac{v_k^Tv_k} {\left(\alpha_k+\frac{\alpha_k}{\eta}\right)v_k^Tq_k}. \tag{20}\]

The procedure for computing the step length is summarized in Algorithm 1.

3. Convergence Analysis

In this section, we present the global convergence analysis of the proposed Algorithm 2, namely the Accelerated Double-Step-Length Method (ADSM). We begin by defining the level set \[\label{a46} \Omega=\{s\mid \|\mathfrak{F}(s)\|\leq \|\mathfrak{F}(s_0)\|\}. \tag{23}\]

Since \(\eta^{-1}\delta_k I\) approximates \(\mathfrak{F}'(s_k)\) along \(s_k\), the following assumption can be made.

4. Numerical Experiments

This section provides a comprehensive summary of the results obtained from numerical experiments on systems of nonlinear equations. The purpose of these tests is to assess the performance and effectiveness of the proposed strategy for solving (1). In addition, the discussion highlights the successful application of the proposed method to motion control problems, thereby demonstrating its relevance and usefulness in real-world scenarios.

4.1. Experiments on Systems of Nonlinear Equations

The efficiency and robustness of the proposed ADSM approach are compared with those of the following existing methods:

1. The double-step-length method (IDS) presented in [13].

2. The modified double-direction method (MDFDD) presented in [22].

3. The accelerated derivative-free method (ALG1) presented in [24].

4.2. Experimental Setup

All computer codes used in this study were implemented in Matlab version 9.4.0 (R2020a) and executed on a personal computer equipped with a 1.80 GHz CPU processor and 8 GB of RAM. To ensure a fair comparison, all four algorithms employed the same line-search strategy described in equation (21).

The parameter settings for the IDS, MDFDD, and ALG1 methods were selected according to the recommendations provided in [13], [22], and [24], respectively. For the proposed ADSM algorithm, the parameters were chosen as \[\omega_{1}=\omega_{2}=10^{-4}, \qquad r=0.3, \qquad \eta=1.9,\] with the sequence \[\chi_{k}=\frac{1}{(k+1)^{2}}.\]

The stopping criteria for all algorithms were defined as follows: the iterative process was terminated if either the total number of iterations exceeded 1000 or the residual norm satisfied \[\|\mathfrak{F}(s_{k})\|\leq 10^{-5}.\]

To evaluate the effectiveness and robustness of the proposed method, extensive numerical experiments were conducted on six benchmark test problems. Each problem was solved using various initial points and dimensions ranging from 10 to 100,000, thereby assessing the scalability and performance of the algorithms on both moderate- and large-scale instances. Table 1 presents the initial points employed in the experiments.

The following test problems are used in the numerical experiments.

Detailed results of the numerical experiments are reported in Tables 2–5. The performance of the methods was evaluated using four principal metrics: the number of iterations (Iter), the number of function evaluations (Fun), the CPU time in seconds (Time(s)), and the residual norm \(\|\mathfrak{F}(s_k)\|\). These metrics provide insight into the efficiency, computational cost, and accuracy of the competing algorithms.

The number of iterations (Iter) represents the total number of algorithmic cycles required to reach convergence. In general, a smaller number of iterations indicates faster convergence and improved algorithmic efficiency. The results presented in the tables show that ADSM consistently requires fewer iterations compared with IDS, MDFDD, and ALG1. This reduction in iteration count demonstrates the stronger convergence behavior of ADSM and its capability to approach the solution more rapidly.

The number of function evaluations (Fun) measures how often the objective function is computed during the line-search procedure. Since function evaluation typically constitutes the dominant computational expense in large-scale optimization, this quantity is crucial. The reported results indicate that ADSM performs significantly fewer function evaluations than the compared methods. This confirms that the proposed strategy reduces line-search overhead and improves overall computational efficiency.

CPU time (Time(s)) reflects the actual computational time required to solve each problem. This metric depends on both the iteration count and the computational work performed per iteration. From a practical standpoint, a shorter CPU time implies a more efficient and scalable algorithm. The numerical results show that ADSM consistently achieves lower CPU times than IDS, MDFDD, and ALG1, frequently solving the test problems in substantially less time.

The accuracy of the computed solutions is assessed through the residual norm \(\|\mathfrak{F}(s_k)\|\), which measures how well the obtained solution satisfies Eq. (1). A sufficiently small residual norm indicates successful convergence to a solution of the original problem. The tables reveal that all four methods achieve residual norms within the prescribed tolerance, confirming their correctness and reliability.

Overall, while all four methods are capable of solving Problem (1), ADSM clearly demonstrates superior performance. It converges more rapidly and requires significantly fewer iterations, fewer function evaluations, and less computational time than the competing methods. These results confirm the efficiency, robustness, and practical effectiveness of the proposed algorithm.

Building upon the performance profile established by Dolan and Mor\(\acute{e}\) [25], Figures 1–3 provide a comprehensive overview of the performance and efficiency of the ADSM, IDS, MDFDD, and ALG1 methods. For each solver \(s\) and problem \(p\), the performance ratio is defined as \[r_{p,s} = \frac{t_{p,s}}{\min_{s\in S}t_{p,s}}, \tag{50}\] where \(t_{p,s}\) denotes the computational effort required by solver \(s\) to solve problem \(p\). The cumulative distribution function \[p_s(\tau) = \frac{1}{n_p}\,\mathrm{card}\{p\in P:r_{p,s}\leq \tau\} \tag{51}\] represents the probability that solver \(s\) performs within a factor \(\tau\) of the best solver. In particular, \(p_s(1)\) gives the fraction of problems for which solver \(s\) is the best, while the behavior of \(p_s(\tau)\) as \(\tau\) increases reflects the robustness of the solver.

Figure 1 presents the performance profile with respect to the number of iterations. Within the interval \(0\leq \tau\leq 5\), ADSM clearly dominates the other methods. It attains the highest value of \(p(\tau)\) near \(\tau=1\), indicating that it solves the largest proportion of test problems using the minimum number of iterations. Its curve rises steeply and approaches full robustness at relatively small values of \(\tau\).

ALG1 ranks second and demonstrates competitive behavior, although its performance remains slightly below that of ADSM for small \(\tau\). MDFDD exhibits moderate efficiency, requiring larger performance ratios before reaching comparable probability levels. IDS shows the slowest growth of \(p(\tau)\) and therefore requires more iterations relative to the best-performing solver.

Figure 2 shows the performance profile based on the number of function evaluations, with \(\tau\) ranging from 0 to 7. ADSM again achieves the highest \(p(\tau)\) at \(\tau=1\), demonstrating superior efficiency in terms of function evaluations. Since function evaluation is often the dominant computational cost in large-scale optimization, this result highlights the practical advantage of ADSM. ALG1 remains highly competitive and closely follows ADSM across increasing values of \(\tau\). IDS ranks third, while MDFDD exhibits the weakest performance in this metric. Notably, MDFDD begins to show significant probability growth only at larger \(\tau\) values, indicating that it generally requires substantially more function evaluations.

Figure 3 presents the performance profile based on CPU time in seconds, with \(0\leq \tau\leq 6\). ADSM consistently attains the largest \(p(\tau)\) values for small performance ratios, indicating that it is most frequently the fastest solver. ALG1 again ranks second and shows competitive time efficiency.

IDS demonstrates moderate performance, whereas MDFDD requires larger \(\tau\) values before reaching comparable probability levels. Although all methods approach full robustness as \(\tau\) increases, ADSM achieves this level more rapidly.

Across all three performance measures, namely iterations (\(0\leq \tau\leq 5\)), function evaluations (\(0\leq \tau\leq 7\)), and CPU time (\(0\leq \tau\leq 6\)), ADSM consistently achieves the highest probability of being the best solver, \(p(1)\), exhibits the fastest growth in \(p(\tau)\), and reaches full robustness at smaller performance factors compared with the competing methods. ALG1 demonstrates stable and competitive behavior, ranking second overall. IDS provides moderate performance, while MDFDD shows comparatively weaker efficiency, particularly in function evaluations and CPU time. These results confirm that ADSM is the most efficient and robust solver among the four methods for the tested large-scale unconstrained optimization problems.

4.3. Application of the Proposed Method to Motion Control Models

This section discusses the motion control problem arising in two-joint planar robotic manipulators. Since (1) is equivalent to the global optimization Problem (4), Algorithm 2 can be applied to solve the resulting nonlinear system.

According to [26], the discrete-time kinematic equation of a two-joint planar manipulator at the position level is given by \[\label{mc1} \varphi(\theta_k)=u_k, \tag{52}\] where \(\varphi(\cdot)\) denotes the nonlinear forward kinematic mapping.

The kinematic mapping has the form \[\label{mc2} \varphi(\theta) = \begin{bmatrix} a_1+a_2\\ b_1+b_2 \end{bmatrix}, \tag{53}\] where \[a_1=l_1\cos(\theta_1), \qquad a_2=l_2\cos(\theta_1+\theta_2),\] and \[b_1=l_1\sin(\theta_1), \qquad b_2=l_2\sin(\theta_1+\theta_2).\] Here, \(l_1\) and \(l_2\) denote the lengths of the first and second links, respectively. The vector \[\theta_k=[\theta_1,\theta_2]^T\in\mathbb{R}^2\] represents the joint angles, while \[\varphi(\theta_k)\in\mathbb{R}^2\] denotes the position of the end effector.

At any time instant \(t_k\in[0,t_f]\), the motion control objective can be formulated as the minimization problem \[\label{mc3} \min_{u_k\in\mathbb{R}^2} \frac{1}{2}\|u_k-u_{dk}\|^2, \tag{54}\] where \(t_f\) denotes the total task duration and \(u_{dk}\) represents the desired trajectory at time \(t_k\).

Following [27], the link lengths are chosen as \(l_1=l_2=1\). A Lissajous curve is used as the desired trajectory of the end effector: \[\label{mc4} u_{dk} = \begin{bmatrix} \frac{3}{2}+\frac{1}{5}\sin\!\left(\frac{2\pi t_k}{5}\right)\\[2mm] \frac{\sqrt{3}}{2}+\frac{1}{5}\sin\!\left(\frac{3\pi t_k}{5}\right) \end{bmatrix}. \tag{55}\]

4.4. Derivation of the Nonlinear System for the Motion Control Problem

In the motion control task, the aim is to drive the end effector of the manipulator to track the desired trajectory \(u_{dk}\). This requirement leads to the optimization problem (54). However, the end-effector position \(u_k\) is determined by the forward kinematics (52). Substituting (52) into (54), we obtain

\[\label{opt_theta} \min_{\theta_k\in\mathbb{R}^2} \frac{1}{2}\|\varphi(\theta_k)-u_{dk}\|^2. \tag{56}\]

Define the merit function \[\mathfrak{h}(\theta_k) = \frac{1}{2}\|\varphi(\theta_k)-u_{dk}\|^2.\] Minimizing (56) is equivalent to solving the nonlinear equation \[\label{motion_F} \mathfrak{F}(\theta_k) = \varphi(\theta_k)-u_{dk} = 0, \tag{57}\] where \(\mathfrak{F}:\mathbb{R}^2\rightarrow\mathbb{R}^2\) represents the residual mapping. Therefore, Algorithm 2 is applied by setting \[s_k=\theta_k,\] and \[\mathfrak{F}(s_k) = \mathfrak{F}(\theta_k) = \varphi(\theta_k)-u_{dk}.\]

The algorithm iteratively updates \(\theta_k\) until the stopping criterion \[\label{stop_motion} \|\mathfrak{F}(\theta_k)\| \leq 10^{-5} \tag{58}\] is satisfied. The error plots reported in the numerical experiments correspond to the residual norm \(\|\mathfrak{F}(\theta_k)\|\), which represents the tracking error between the actual end-effector position \(\varphi(\theta_k)\) and the desired trajectory \(u_{dk}\).

We established the following parameters for this experiment: \[\omega_1=\omega_2=10^{-4}, \qquad q=10^{-4}, \qquad r=0.3, \qquad \chi_k=\frac{1}{(k+1)^2},\] with \(t_f=10\) seconds and the initial point \[\theta_0= \left[0,\frac{\pi}{3}\right]^T.\] Furthermore, the working interval \([0,10]\) is divided into 200 equal subintervals.

Figures 4–7 illustrate the robot trajectories generated by the ADSM algorithm. The results show that the proposed method successfully guides the end effector to follow the desired trajectory. In addition, Figures 6 and 7 present the tracking error of the end effector during the motion process. It can be observed that the error decreases rapidly and remains approximately at the level of \(10^{-5}\), which satisfies the prescribed stopping tolerance. These results confirm the effectiveness and applicability of the proposed ADSM method for solving motion control problems in robotic manipulators.

5. Conclusion

This study presents a hybridized step-size method specifically designed to tackle systems of nonlinear equations. A notable advantage of the proposed method is the use of two different line-search techniques and its derivative-free characteristic, which significantly enhances its efficiency for solving problem (1). This feature makes it particularly advantageous in scenarios where obtaining derivative information may be challenging or costly. To thoroughly evaluate the performance of the proposed method, we conducted extensive numerical comparisons using a set of large-scale test problems. The results are summarized in Tables 2–4 and visually represented in Figures 1–3. Our findings reveal that the proposed method consistently requires fewer iterations and less CPU time compared with several other established techniques in the field. Additionally, we applied the proposed method to a motion control problem involving two-joint planar robotic manipulators, demonstrating its practical applicability. The outcomes of this application are illustrated in Figures 4–7, underscoring the method’s effectiveness in real-world situations. In future work, we intend to explore the development of a two-step-size method for solving motion control problems with four degrees of freedom.

Author Contributions

Abubakar Sani Halilu: Methodology, Software, Writing–original draft. Arunava Majumder and Mohamad Afendee Mohamed: Validation, Writing–review and editing, Supervision. Mohammed Yusuf Waziri and Sulaiman M. Ibrahim: Conceptualization, Methodology. Kabiru Ahmed and Muhammad Abdullahi: Validation, Formal analysis, Investigation.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.