#### Abstract

The scope of this study is to present a contribution to the geometrically nonlinear free and forced vibration of multiple-stepped beams, based on the theories of Euler–Bernoulli and von Karman, in order to calculate their corresponding amplitude-dependent modes and frequencies. Discrete expressions of the strain energy and kinetic energies are derived, and Hamilton’s principle is applied to reduce the problem to a solution of a nonlinear algebraic system and then solved by an approximate method. The forced vibration is then studied based on a multimode approach. The effect of nonlinearity on the dynamic behaviour of multistepped beams in the free and forced vibration is demonstrated and discussed. The effect of varying some geometrical parameters of the stepped beams in the free and forced cases is investigated and illustrated, among which is the variation in the level of excitation.

#### 1. Introduction

In different engineering fields, such as the automotive, aeronautical, civil, or mechanical engineering, many structural components that permit the reinforcement of the whole system, a good resistance to the working loads, or the transmission of motions, can be modeled as stepped beams. In their operation, this type of beams is often subjected to a high dynamic load that needs to be investigated in order to be able and reliably predict the dynamic response of the beam. Many research works have examined the linear behaviour of this type of beams, making the analysis very easy, but leading generally to unreliable results when large vibration amplitudes are encountered. Since, strictly speaking, no system behaves always linearly, the linear results are often misleading. As the nonlinear effects arise in the case of the large displacements induced in many dynamic operating conditions, this led us to include the influence of geometrical nonlinearity in the present analysis.

A literature survey on the transverse linear vibration of stepped beams goes backwards to Taleb and Suppiger [1], who presented some of the results of the integral equation theory, applied the Cauchy function method to achieve an approximate estimation of the fundamental frequency and modal configuration of the transverse vibration of beams with one step, and compared it with the exact solution. The study of the performance of a four-degree-of-freedom-per-node element model for the vibration analysis of uniform and multistepped beams was carried out by Balasubramanian and Subramanian [2]. Laura et al. [3] used the finite element method and compared their results with the experimental ones. Popplewell and Daqing Chang [4] used the Rayleigh–Ritz method to determine the unknown constants of the transverse deflection of a stepped beam. Krishnan et al. [5] discussed the use of the finite difference method in the study of stepped beams. Naguleswaran [6–9] studied the vibration of stepped beams in different configurations involving different boundary conditions, various numbers of steps, and types of cross sections. The natural frequencies were determined by iterative procedures based on linear interpolations. The vibration of a stepped beam with elastic supports at the ends was studied by De Rosa [10]. Sato [11] presented the analytical results relative to the effect of an abrupt change of cross section on the natural frequency of a beam with a groove, which were very well correlated with the experimental results. Sairigul and Aksu [12] used a method based on variational principles in conjunction with the finite difference technique to investigate the free vibration analysis of stepped Timoshenko beams. The Adomian decomposition method (ADM) was applied by Mao in [13, 14] in order to examine the free vibration of stepped and multistepped beams. A new method for calculating the frequencies and mode shapes was presented by Lin and Ng [15], who considered this method to be well adapted to forced-vibration applications.

Although most work on stepped beams has focused on linear vibration analysis, few papers have included the effect of geometric nonlinearity in their studies. Using the transfer matrix method, Sato [16] exposed a study on the nonlinear free vibration of stepped beams. Nuttawit Wattanasakulpong and Arisara Chaikittiratana [17] studied the problems of linear and nonlinear free vibrations of stepped beams using the Adomian modified decomposition method. None of the articles mentioned above have investigated the geometrically nonlinear dynamic response of forced stepped beams.

This study presents an extension of the works previously initiated by Benamar et al. [18–20], with the aim to contribute to the modal analysis of the geometrically nonlinear vibrations of straight thin structures, such as beams, plates, and circular cylindrical shells. Azrar et al. [21] presented an analysis of the geometrically nonlinear forced vibrations of uniform beams based on the multimode model. El Kadiri et al. [22] developed a simplified multimode approach (MMA), which is, according to the authors, easy to use for the study of the geometrically nonlinear forced vibrations of clamped-clamped and simply supported-clamped beams at both sides. El Kadiri and Benamar [23] extended this method to investigate the case of rectangular plates. Merrimi et al. [24] studied a cracked beam excited by a harmonic concentrated force, using both multimode and single-mode approaches. Fakhreddine et al. [25] investigated using the same method the geometrically nonlinear vibrations of beams resting on multiple supports and subjected to concentrated and uniformly distributed harmonic forces. Fakhreddine et al. [26] also extended this method to the case of nonlinear forced vibration of beams carrying concentric masses of different magnitudes and for vibration amplitudes.

The main objective of this work was to present a contribution to the geometrically nonlinear free and forced vibration of beams with multiple steps based on the Euler–Bernoulli and von Karman theories. First, the natural frequencies are calculated by the Newton–Raphson algorithm leading to the determination of linear modes. The discretized expressions for the bending strain energy, the axial strain energy due to the nonlinear stretching forces, and the kinetic energy are then derived. The problem is reduced to a nonlinear algebraic system using Hamilton’s principle and then solved by the so-called second formulation reported in [22]. The numerical method convergence study has been carried out in the free and forced vibrations. The effects of the changes in cross section, inertia, and length ratios are analyzed and then illustrated. Considering a uniformly distributed excitation, the forced vibration case is studied on the basis of a multimode approach. The nonlinear frequency response functions are calculated and presented in the neighborhood of the predominant mode. The effects of the variation in the excitation level, the type of cross section, and the length ratio of stepped beams are investigated and illustrated. The presented method can also be extended for other types of beams including nano-/microbeams, previously investigated in [27] using the modified strain gradient and hyperbolic shear deformation beam theories and in [28] using the Euler–Bernoulli beam theory via the enhanced Eringen differential model. It can also be extended for the case of free vibration nanotubes, earlier studied in [29], using the discrete singular convolution technique.

#### 2. General Formulation

##### 2.1. Linear Formulation

This study concerns nonuniform homogeneous beams of the Euler–Bernoulli type, corresponding to the three configurations shown in Figure 1, characterized by *n* different steps at various points. The first case is a beam with a rectangular cross-section variation in the width (a), the second in the height (b), and the third in both (c). In all cases, the beams are vibrating transversely.

Figure 2 shows that each of the multistepped beams considered is divided into *n* different steps of uniform cross section and is supported by linear and rotational springs at both ends.

According to Taleb and Suppiger [1], the transverse free vibration of any segment *i* of the beam with a uniform cross section is governed by the following equation:where

Each step of the stepped beams has a specified length, thickness, width, Young’s modulus, area of cross section, moment of inertia, mass per unit length, and eigenvalue parameter (, , and ). The variation of each step is presented by the ratios of the cross sections, moments of inertia, step, length, and eigenvalue parameters (, and ).

The ‘active’ dimension denoted by of beams (a), (b), and (c) is, respectively, the width, the height, and both of them, which allow one to write the ratio of steps and moments of inertia as follows.

For beam (a),

For beam (b),

For beam (c),

The general solution of equation (1) leads to the linear mode shapes, which will be used as basic functions in the nonlinear theory detailed in the following. The transverse displacement function of the multistepped beam, shown in Figure 2, can be expressed in each span as given by [25]where , and also,

By varying indices *i* from (1) to (7) and *j* from 1 to *N* + 1, the satisfaction of the boundary conditions together with the compatibility conditions and the transfer matrix leads to a homogeneous system. The boundary conditions for a beam segment supported by linear and rotational springs at both sides are given in [25].

At the left side,

At the right side,

The stiffness values of the transverse and rotary springs at both ends of the beam are denoted as . The compatibility conditions for each beam step are expressed as shown in [13]:

##### 2.2. The Transfer Matrix Method

Passing through the compatibility conditions, it is found that the constants ( , and ) in the step depend on those of the step ( , and ), which can be presented in a matrix form as in [30]:with

From equation (12), one obtainswhere are matrices, which form the transfer matrix expressed by

Equations (17) and (18) lead to

Under successive repetition of equation (12), the four constants of the first segment ( , and ) appear to be related with those of the final segment (, and ), which decreases the number of independent constants to four.

The four constants (, and ) are determined when the boundary conditions are fully satisfied and expressed as given by Chajdi et al. [31], but in this case, the eigenvalue parameter differs from one span to another. The boundary conditions at the right side lead to the following expressions:

Equations (18) and (19) can be expressed in a matrix form as follows:such as

Substitution of (17) into equation (20) leads to

The boundary conditions at the left side can be formulated in a matrix form as given in [31]:

Equations (22) and (23) lead towhere

The homogeneous system is solved by the Newton–Raphson algorithm, where the determinant of matrix *Y* is assigned to be equal to zero in order to allow nontrivial solutions to exist:

After determining the natural frequencies from equation (26), the constants (*A*_{j}, *B*_{j}, *C*_{j} and *D*_{j}) are calculated by the boundary and continuity conditions written above.

#### 3. Nonlinear Formulation

##### 3.1. Free Vibration

The dynamic behaviour of the structure is determined by Hamilton’s principle which is formally expressed as follows:

*T* is the beam kinetic energy. *V* is the total strain energy of the multistep beam, expressed as the sum of the energy due to the bending *V*_{f} and the energy *V*_{a} due to the nonlinear stretching forces induced by the large deflections.

Setting , the kinetic energy can be expressed as in [32]:

Setting , the deformation energy due to bending can be expressed as in [32]:

Equations (28) and (29) can be developed as the following form:

The multistep beam has a variable cross section that leads to a variable moment of inertia, which changes from step to step. The axial deformation energy due to the nonlinear stretching forces is calculated as follows:

The normal load induced by the nonlinear stretching force of the beam is assumed, when the axial inertia forces are negligible, to be constant along the whole beam span as in [33].

Setting , the nonlinear stretching force can be given aswhere

The nonlinear axial strain is expressed as

Substituting equations (33) and (35) in equation (32) gives the following equation:

By assuming a harmonic motion and expanding the transverse displacement as a finite series of basic spatial functions, it is possible to writewhere and are the multistepped beam *i*^{th} linear mode and basic function contribution coefficient, respectively. The replacement of the new form in equations (30), (31), and (36) and the discretization lead to the following expressions for the kinetic energy, the axial strain energy due to nonlinear stretching forces, and the strain energy due to bending as in [22]:in which , , and are the mass matrix due to , the rigidity matrix due to , and the nonlinearity tensor due to , expressed by

Using equations (38)–(40) and applying Hamilton’s principle lead to

According to El Kadiri et al. [22], the above equation can be expressed in a matrix form as follows:where , and are the general terms of , and . The column vector of the contribution coefficients and the frequency are the unknowns. Before they are determined, equation (43) is put in a dimensionless form by replacing the dimensional parameters by the corresponding nondimensional ones:

Replacing equations (44)–(48) in equation (42) or (43) results in the following equation:which can also be written as

The determination of the unknowns is based on the so-called second formulation developed previously in [22], where it was shown that the contribution of the fundamental mode () remains predominant for the full range of the vibration amplitudes considered. The other contributions, which are small in comparison with (), are denoted by . Consequently, following [22], one can write , which allows equation (50) to become, in the vicinity of the *r*^{th} mode,where

This approximate linear system, given in [22], makes it very easy to calculate the contribution coefficients.

##### 3.2. Forced Vibration

The study of nonlinear forced vibrations is made in the three configurations, which differ by the type of the cross section. In all cases, the beam is subjected to a distributed force, which excites the modes of the structure through a series of generalized forces. The distributed harmonic force is denoted by and depends on the expression for , the portion of the beam in which the excitation is applied, and the mode considered. The expressions for the distributed harmonic force are extracted from [20] and written as follows:

To investigate the dynamic behaviour of the conservative system, Hamilton principle, in which the forcing term has been introduced, has been applied as in [26]:

The kinetic energy, the total strain energy, and the work done by the external loads of the multistep beam are denoted by *T*, *V*, and *W*_{F}. After calculations, one gets

The dimensionless generalized distributed forces are denoted by :

Equation (57) is made dimensionless by using the dimensionless parameters in equation (55) to get

According to El Kadiri et al. [22], an approximate linear system has been obtained as in the nonlinear free vibration case developed above, which permits equation (57) to be expressed as follows:

The necessary use of the multimode approach in the neighborhood of the predominant mode makes it very easy to calculate the contribution coefficients, as in [22].

#### 4. Convergence Study

The purpose of this section is to investigate the convergence of the series expansion used to solve the nonlinear vibration problem via the amplitude equation obtained in the numerical approach exposed above. These analyses deal with a two-step fully clamped beam of type (b), where , , and , subjected to a free and forced vibration.

##### 4.1. Free Vibration

In order to solve equation (51), the assumptions stated by El Kadiri et al. were used to simplify the nonlinear algebraic system obtained. The contribution coefficients have been calculated as well as the associated nonlinear frequency. The normalized first nonlinear mode , for a given value of , is obtained as a series of basic functions involving the beam modal parameters depending on and given as follows:

The predominant term, relative to the first linear mode shape, is , while the other terms proportional to the symmetric and antisymmetric higher modes present the corrections due to the geometrical nonlinearity.

Table 1 shows the result of contribution coefficients obtained using equation (51); these contribution coefficients are included in equation (59) in order to calculate the normalized first nonlinear mode . Its maximal value is denoted by *W*_{max} and is located at the middle of the beam. From Table 1, some observations can be clearly made:(i)The contribution of the fundamental mode remains predominant for the full range of the vibration amplitudes considered. The other contributions denoted by remain small with respect to .(ii)The contribution coefficients of the modes decrease as *i* indices increase.(iii)The contribution coefficients increased according to the maximum nondimensional amplitude .

From these observations, it is quite interesting to verify the convergence of the nonlinear beam deflection equation, which is obtained as a series of basic functions. The increase of these basic functions added in equation (59) enhanced the correction due to geometrical nonlinearity. The effect of varying the number of basic function series on the convergence of the numerical method has been investigated by comparing the results of a nonlinear frequency estimated . This latter can be expressed as in [26]:

The nonlinear frequency estimated can be extracted from equation (60), where presents the contribution coefficient vector already obtained from the second formulation, using *N* basic functions. For , the nonlinear frequency is calculated for different cases, in which the number of basic functions *N* increased from 7 to 10. The same test is repeated for different values of , varying from 0.5 to 0.8. The differences between and are calculated, as well as the residual of each test, and summarized hereafter in Table 2.

The numerical results presented in Table 2 show that the differences for the four tests decrease by increasing the number of basic functions, while the residual increases by increasing the *a*_{1} value. The difference of the first test is equal to , and after increasing the number of basic functions to 10, this difference decreases to a low value of , meaning that(i)The use of 6 basic functions induces an error of with respect to (ii)The use of 7 basic functions induces an error of with respect to (iii)The use of 8 basic functions induces an error of with respect to (iv)The use of 9 basic functions induces an error of with respect to

Also, at *a*_{1} = 0.8, the difference for 7 basic functions increases from to , and by adding 3 other terms, the difference decreased to . In other terms,(i)The use of 6 basic functions induces an error of with respect to (ii)The use of 7 basic functions induces an error of with respect to (iii)The use of 8 basic functions induces an error of with respect to (iv)The use of 9 basic functions induces an error of with respect to

When the basic functions’ number increases and exceeds ten, the difference between the final iteration of and its precedent decreases by value less than for the first test and less than for the last one, whereas the residual will be increased by less than these values. These corrections are almost negligible compared to . These results are considered to be tightly convergent in the nonlinear free vibration case using ten or more basic functions.

##### 4.2. Forced Vibration

These analyses deal with a two-step fully clamped beam of type (b), where , , and ), subjected to a uniformly distributed force over an interval of [0.4, 0.6]. Firstly, the contribution coefficients are obtained using multimode approach (58) and then summarized in Table 3.

Afterwards, these contribution coefficients are included in equation (59) in order to calculate the normalized first nonlinear mode . Its maximal value is denoted by *W*_{max} and is located at the middle of the beam. The same observations can be concluded from Table 3. Furthermore, it can be clearly seen that, in the forced case, the contribution coefficients of the even modes have very low values. However, a test has been carried out to predict the effect of neglecting the antisymmetric terms in equation (59) on the normalized first nonlinear mode. Two first normalized nonlinear modes have been calculated: the first one is obtained as shown in equation (59), and the second is also obtained as a series of basic functions involving only the symmetrical terms and denoted by . The results are summarized in Table 4.

The numerical results presented in Table 4 prove that, by neglecting the antisymmetric modes, the use of equation (59) leads to the same results. The normalized first nonlinear mode can also be written as

The effect of varying the basic functions’ number on the convergence of equation (62) has been investigated by comparing the normalized first nonlinear mode obtained for five basic functions for three scenarios.

*Scenario 1. *The nonlinear beam deflection assumes the use of one basic function:

*Scenario 2. *The nonlinear beam deflection assumes the use of two basic functions:

*Scenario 3. *The nonlinear beam deflection assumes the use of three basic functions:The convergence between and has been examined by calculating the integral of these nonlinear modes over an interval of [0, 1] using the Simpson numerical integration method. The results are obtained for different values of in the case where is equal to 100, 200, and 300. The comparison was performed and is summarized in Table 5.

In Table 5, the final normalized first nonlinear mode (62) was compared with the three other scenarios of nonlinear modes as shown in equations (63)–(65). The comparison was performed by varying the frequency ratio for each excitation. The results obtained for the scenario where and show a small difference of for the first scenario, and it decreases to for the last one, in other terms:(i)The use of Scenario 1 induces an error of with respect to (ii)The use of Scenario 2 induces an error of with respect to (iii)The use of Scenario 3 induces an error of with respect to Also, for the same excitation, if , increases until reaching for the first scenario, and increases to for the last one, meaning that(i)The use of Scenario 1 induces an error of with respect to (ii)The use of Scenario 2 induces an error of with respect to (iii)The use of Scenario 3 induces an error of with respect to It can be noticed that the change in the excitation level from 100 to 400 with increases the error involved between and by a small variation. For , this error is equal to , and for , this error increases to . Indeed, the effect of the frequency ratio variation is more remarkable than the force variation effect on the convergence of the normalized first nonlinear mode. The differences between and for the value of the frequency ratio and for the three excitations considered remain small. However, the increase in the number of basic functions exceeding five symmetric terms in equation (62) provides a correction, due to the geometric nonlinearity, less than . This correction is almost negligible compared to ; hence, the numerical method proposed is well converged.

#### 5. Numerical Results and Discussion

##### 5.1. Linear Free Vibration

A computer program has been performed using MATLAB software, which analyses the effects of the number and position of each step of the multiple-stepped beam, as well as the change in the ratios of the step cross sections and moments of inertia on the vibration frequencies and mode shapes of the beam. To verify the validity and precision of this analysis, three cases are examined first, named beam (a), beam (b), and beam (c), in which the beam cross sections vary in width, in depth, and in both. Secondly, the three beams are assumed to be fixed on both sides; for this purpose, an infinite value is adopted for the value of the stiffness of the linear and rotary springs at both ends (, and ). Thirdly, the values of the vibration frequencies for the first four vibration modes are calculated and compared with those available in [6, 13, 14] as shown in Tables 6 and 7.

Finally, the program developed for the case of tapered beams is validated by modelling the tapered beam with a 12-step beam as shown in Figure 3, maintaining the same taper ratio and the same angle of inclination.

The results obtained are compared in Table 8 with those of [34] and show an excellent agreement, with a relative difference not exceeding 0.0094% for the whole range of frequencies considered. The increase in the number of steps is expected to further reduce this difference.

##### 5.2. Nonlinear Free Vibration

The nonlinear vibration of tapered beams similar to that of [33–35] is now investigated in the same way, and a comparison is made between the results obtained here, based on the second formulation and summarized in Table 9, and those of [33–35].

An excellent agreement can be noticed since the relative difference does not exceed 0.2% for beam (a) and 1.8% for beam (b). Increasing the number of steps is expected to further reduce this error.

##### 5.3. The Effect of the Section Ratio

The effect of the cross-section ratio and the vibration amplitude on the beam nonlinear behaviour is illustrated in Figures 4–6 by means of the backbone curves associated to the first nonlinear mode shape of a fully fixed beam with two steps (, , and ), where the value of *m* is increasing for the three beams (a, b, and c).

It is clear that the nonlinearity remains of the hardening type and that the increase in the section ratio () leads to an increase in the frequency ratio.

##### 5.4. The Effect of the Moment of Inertia Ratio

A comparison between the three beam configurations with the same section ratio ( and ) and length ratio () but with a different inertia ratio is performed and illustrated in Figures 7–9.

Figures 7–9 show the stress distributions along the beam for each configuration, derived from the linear theory and the nonlinear model. The effect of geometric nonlinearity is clearly shown in the figures. The calculation of the percentage correction introduced by the nonlinear theory compared to the linear theory shows that it can reach 36% for beam (a), 23% for beam (b), and 59% for beam (c) as summarized in Table 10.

In addition, it can be noticed that the maximum value of the curvature for the linear and nonlinear cases for beam (a) is located near the clamped ends, while it is located in the middle for beam (b). For beam (c), the location of the maximum curvature is not the same in the linear and nonlinear cases: in the linear case, it is located in the middle of the beam, while in the nonlinear case, it is located at its end.

##### 5.5. The Effect of the Step Length Ratio

The effect of the step length ratio on the nonlinear behaviour of beam (b) is illustrated in Figure 10 which gives the backbone curves associated to the first nonlinear mode shape for four cases: , , and .

Figure 10 shows the effect of the length ratios on the geometrically nonlinear beam behaviour. The modification of these ratios is obtained by increasing , resulting in a decrease in and . It should be noted that, by increasing the length ratio , the frequency ratio increases accordingly.

##### 5.6. Nonlinear Forced Vibration

The forced nonlinear dynamic behaviour is investigated for the three types of beams (a, b, and c) with , , and and excited by a uniformly distributed force over an interval of [0.4, 0.6]. In order to examine the forced case, and as mentioned in the theoretical section, it is required to determine the predominant mode corresponding to each level of excitation in order to perform the study in its vicinity. The generalized forces have thereby been calculated and listed in Table 11.

The calculation of the generalized forces led to the conclusion that the first mode is predominant. Figures 11–14 show the nonlinear frequency response functions in the vicinity of the first vibration mode.

Figure 11 shows a comparison between the linear and nonlinear frequency response functions of configuration 2, where the beam is excited by a distributed force uniformly distributed over an interval of [0.4, 0.6].

Figure 12 shows a comparison between the frequency response curves obtained for the three beams considered, having different inertias. The comparison is conducted using a multimodal approach for three beams subjected to a uniformly distributed force over an interval of [0.4, 0.6], with .

Figure 13 shows a comparison between the response curves obtained for the three configurations: , , and . The beam is subjected to a uniformly distributed force in the interval of [0.4, 0.6].

Figure 14 presents a comparison between the frequency response curves of beam (b), excited by different levels of uniformly distributed forces in the range [0.4, 0.6]. The effect of geometrical nonlinearity can be seen with the multivalued region in which the jump phenomenon commonly found in nonlinear vibration systems can occur. The hardening behaviour is also clearly visible, as shown in Figure 11. As the vibration amplitude increases, the frequency ratio increases accordingly. In addition, Figure 12 shows how the nonlinear effect changes with the type of cross section when one moves from one beam to another. It appears from this figure that the hardening effect in beam (b) is stronger than in the other two. Figure 13 shows the effect of the change in length ratio on the geometrically nonlinear behaviour. The increase in the length ratio makes the hardening behaviour more accentuated. Figure 14 illustrates the disproportionate increase in the frequency response amplitude relative to the increase in the intensity of the distributed excitation force.

#### 6. Conclusion

The geometrically nonlinear free and forced vibrations of multistepped beams were examined, based on Euler–Bernoulli’s beam theory and von Karman’s nonlinearity assumptions. The solution of the linear problem was obtained for different configurations of multistepped beams. The calculation of the frequencies was performed using the Newton–Raphson algorithm, and the mode shapes were determined afterwards. Discrete expressions of the strain and kinetic energies were derived, and the problem was reduced using Hamilton’s principle to a nonlinear algebraic system solved by the second formulation. The numerical method convergence study was carried out in the free and forced vibrations, in which the results showed a good convergence. The results obtained were compared with those obtained using the Galerkin method in order to verify the validity and precision of the proposed approach and showed a good agreement since the relative difference does not exceed 1.8% for . The effect of the geometrical nonlinearity on the dynamic behaviour of multistepped beams was illustrated, as well as the effects of the change in the beam cross section, inertia, and length ratios. It was clearly noticed that the nonlinearity remains of the hardening type and that the increase in the section or the length ratios leads to the increase of the frequency ratio. Furthermore, it appears that the location of the maximum curvature can differ from the linear to the nonlinear cases. The forced vibration was analyzed using the multimode approach. The effect of the excitation level, the type of section, and the length ratio on the nonlinear forced behaviour of the beam was illustrated and discussed. Besides, it can be clearly seen that the hardening effect is greatest when the cross section of the beams increases in height or when the second length ratio increases. Furthermore, a convergence study was also presented to show how the results of the numerical method are well convergent.

#### Data Availability

No data were used to support this study.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.