## 1 Introduction  
Support vector regression [1; 2] is an extension of the Support Vector Machine (SVM), which was initially introduced to solve classification problems [3; 4]. In SVR, the goal is to find a function that best fits the data, providing a good generalization as well as being robust against noise by considering a regularized \(\varepsilon\)-insensitive loss function. The SVR typically solves the Lagrangian dual problem with \(2L\) Lagrange multipliers, twice the number of training samples, by taking advantage of the kernel trick that makes it possible to implicitly transform the input data into a higher dimensional space. The landscape of SVM/SVR optimization methods is rich and diverse. Cutting plane algorithms with improved line search techniques have been used to solve nonsmooth convex linear SVM problems in the primal [5; 6]. Proximal gradient-based methods have found success in solving linearly constrained Quadratic Programming (QP) problems with box constraints, applicable to Huberized SVMs and various SVM variants [7; 8]. Both subgradient- and gradient-based methods have been explored for solving strongly convex optimization problems, with numerical experiments conducted on soft-margin linear SVMs [9]. Recent studies have reformulated SVR with the \(l_{2}\) loss function, leading to nonsmooth unconstrained dual problems with \(L\) Lagrange multipliers. These studies investigated various approaches, including solving smooth approximations, applying generalized derivatives, and employing functional iterative and Newton methods [10; 11; 12]. [PERSON] et al. [13] addressed nonsmooth dual problems in nonparallel support vector ordinal regression using an alternating direction method of multipliers, requiring kernel computations at each iteration. [PERSON] and [PERSON] [14] introduced a semismooth Newton method for solving both support vector classification and regression problems with the \(l_{2}\) loss function in the primal. As the number of training samples in SVM/SVR grows, computational efficiency becomes a major challenge due to the appearance of the massive kernel matrix. To address this issue, the SMO algorithm is used to decompose the problem into smaller subproblems where only two Lagrange multipliers are updated in each iteration. Originally developed by [PERSON] [15] for smooth dual QP problems with \(2L\) Lagrange multipliers, the SMO algorithm has undergone numerous developments. [PERSON] et al. [16] and [PERSON] et al. [17] introduced first-order and second-order information, respectively, into the WSS procedure, a key component of SMO. [PERSON] and [PERSON] [18] introduced an SMO algorithm for solving nonsmooth, indeed piecewise quadratic, optimization problems by dealing with \(L\) optimization parameters. Other studies, such as [PERSON] et al. [19] and [PERSON] et al. [20] extended this approach by using first-order (FO) information for WSS. Additionally, [PERSON] [21; 22] further extended the WSS procedure by integrating SO-like information. This extension involved the innovative concept of minimizing an upper bound on the difference between consecutive loss function values, effectively addressing the challenges of solving the piecewise quadratic dual optimization problem. In [23], a WSS procedure was also developed by combining the advantages of the methods in [19; 20] based on the FO-like information and the method in [24] for solving piecewise quadratic problem arising in nonparallel SVR. In [25], an SMO algorithm for solving QP problems that arise in LSSVM was developed, taking advantage of handling \(L\) variables. This algorithm utilizes the WSS procedure with FO information. Later, in [26], the SMO algorithm for LSSVM was extended by comparing the performance of WSS procedures that employ both first-order and second-order information. In particular, studies [17; 21; 26] have consistently demonstrated the advantages of SO-based WSS over FO-based approaches in terms of efficiency.