Approximate numerical differentiation methods have also been used to calculate eigenvalue sensitivities, including finite differentiation methods (FD) and complex variable methods (CVM). The FD computes derivatives by adding a real-valued perturbation to the objective function, in this case, the GEP equation. Computing an FD approximation of a derivative is straightforward. However, the difficulties are given by the step size dependency, computational cost, and subtraction errors inherent to the method [22; 23; 24]. In contrast, CVM uses complex variables to compute the sensitivities of eigenvalues. [PERSON] and [PERSON] [25] presented an accurate, robust, and easy-to-implement CVM to differentiate eigenvalues. This method is only available for first-order sensitivities, and it requires numerical solvers for the GEP capable of handling complex algebra. Furthermore, [PERSON] et al. [26] presented a CVM to compute arbitrary-order eigenvalue sensitivities in Phononic Materials using the QR decomposition method; however, the method becomes impractical for a large number of sensitivities in systems with a high number of DOF because it requires several iterations over large matrices. In addition, it calculates the derivatives for all the eigenvalues corresponding to all the DOF in the system without distinction, increasing the computational complexity. Recently, [PERSON] et al. [27] introduced a CVM to differentiate eigenvalues and eigenvectors using hyperdual numbers. Although arbitrary-order sensitivities can be computed with such approach, the method involves the solution of multiple linear systems of equations for computing the sensitivities of both the eigenvaluesand eigenvectors, which increases the computational cost. In addition, the method has problems when calculating the sensitivities of systems that exhibit repeated eigenvalues. The calculation of the sensitivities of the eigenvectors for systems with repeated eigenvalues has been reported as a challenging and open area of research [6]. A model with distinct eigenvalues has one unique eigenvector for each eigenvalue. This implies that the eigenvectors are orthonormal and continuously differentiable with respect to any model input parameters. In contrast, a model with repeated or close eigenvalues has many orthonormal eigenvectors associated with multiple eigenvalues. This adversely affects the differentiability and continuity of the eigenvectors with respect to the model input parameters [6]. In the literature, there are two commonly used methods for calculating sensitivities of eigenvectors: model superposition methods [17; 28] and Nelson-like methods [29]. These methods use semi-analytical expressions and work when repeated and distinct eigenvalues are available. The modal superposition methods calculate the sensitivities of the eigenvectors with a linear superposition of all the structure's eigenvalues and eigenvectors. As obtaining all the modes in large systems is computationally expensive, these methods become impractical. Nonetheless, improvements have been proposed where the eigenvector sensitivities can be approximated using only a smaller set of modes [30; 31; 32]. As the superposition is truncated with these modifications the computations become prone to numerical errors. In addition, the modal superposition methods are complicated to implement because of the complexity of the semi-analytical equations, especially for higher-order sensitivities [33]. Nelson-like methods were originally presented in [29], with the main advantage of requiring only the information from the modes of interest to compute the sensitivities of the eigenvectors in distinct eigenvalues problems. Several authors [34; 35; 36; 37] made extensions to the initial approach, including the capabilities to deal with repeated eigenvalues [38; 39], and to expand the method to obtain higher-order sensitivities [21; 40].