Pipe route optimization to avoid undesired vibration by using JuliaFEM

An optimization routine was applied to high pressure fuel pipes to avoid resonance in a heavily vibrating environment. The optimization process and also the natural frequency calculations in every iteration were completely performed with the high-level programming language Julia; the optimization process was performed with the JuMP optimization environment, and the frequencies where calculated with JuliaFEM finite element method solver platform. The benefit of this kind of embedded implementation is a quick response which yields a pleasant development environment to focus on the essential–the choice of the optimization strategy.


Introduction
Designing pipe routes can be utterly tricky especially in a heavily vibrating environment. An optimization routine was applied to avoid possible resonance caused by the vibration. For example, in the high pressure fuel lines of reciprocating medium speed diesel engines and pumps, rotation excitations create undesirable mechanical vibration that is often the main trigger to cause damage to the system [1]. To avoid the need of fixing these issues on the field and to increase reliability, a simulation driven mindset was applied. The optimization procedure was completed with JuMP [2], an optimization environment written in the fast and high-level programming language Julia [3]. JuliaFEM, a Finite Element Method framework also written in the Julia language, was used to create the mesh with 3D structural beam elements and also to solve the natural frequencies of the pipe route alternatives during the optimization [4][5][6][7]. In other words, the chosen approach allows the whole optimization procedure to be performed within the Julia programming environment.
Using only one programming environment has many benefits: the response times are short, and no additional time is spent on launching and closing applications like others have chosen to do [8]. The environment is flexible, and time can be used to focus on testing different optimization strategies and algorithms. Only few researches in the field have considered vibration in the pipe route optimization process. In most cases, the focus has been on pipe length and avoiding obstacles. Vibration examination of pipe systems has mostly been studied independently and not included in the actual optimization process [9,10].
Choosing the most suitable optimization strategy was essential in the pipe route design optimization procedures. It is important to choose a relevant optimization algorithm and equivalent constraints that suite the problem layout. Graph optimization based algorithms has been applied on similar optimization problems [11,12]. Projection and geodesic-based pipe routing algorithms have also been considered [13]. Genetic algorithms are commonly applied since they are most useful when the objective has several local minimums where deterministic methods could easily be stuck in [14][15][16][17][18]. Genetic algorithms are more effective to find the globally best solutions. Particle swarm optimization is also a method that is useful in finding the global minimum for a multi-objective and multi-constrained optimization problem [19]. Pareto optimization provides solutions for optimization problems with multiple goals, such as pipe route optimization [20]. Ant colony algorithms have been applied to classical multi-objective optimization with the weighted sum approach, human-computer cooperation improved optimization and a biobjective automated distributed system design method [21][22][23]. The duct routing has been approached with constraint-based design [24]. Also algorithms for automated generation of pipe routes using simulated annealing optimization schemes have been applied [25]. The augmented Lagrangian algorithm (AUGLAG) [26] was considered in this work since its implementation in JuMP allows nonlinear constraints to be included in the optimization process. Also AUGLAG is pointed out as a good choice if the quality of computed solutions and robustness is important [27].
The objective of this work was to find a suitable optimization strategy to optimize high pressure fuel pipe routes to avoid undesired resonance and to implement the optimization procedure in the Julia programming environment. Also, one goal was to find robust default settings for a working example case.

Optimization strategies
Several optimization strategies were implemented during the study. The objective functions and the emphasis of the objectives in multi-objective cases were varied.
The length of the pipe was included in the objective function in every case in order to minimize costs. The length can be calculated as a sum of all coordinate distances: where x = {p 1 , p 2 , ..., p n } are the current pipe route coordinate points. Figure 1 visualizes the objective function in a case where the pipe route consists of three points and one of them was allowed to move freely in a bounded 2D-coordinate system. The pipe route is expressed by a red line with three black crosses as the coordinate points. In this case, the objective function minima lie on the line connecting the two fixed points.
In addition to the length, also natural frequencies were optimized. Often structural vibrations of machines are quite low. This has lead to design guidelines where the pipelines are to be designed as rigid as possible. Hence frequencies under 125 Hz are undesired in the upcoming examples. Also every 50 Hz between 150 Hz and 1000 Hz are avoided. These are the excitation frequencies caused by the pumps and engine rotation. In the following two of the considered methods, the application of these excitations and frequency limits are shortly introduced.

Strategy 1: Minimizing distances to excitations
In the first method, all of the natural frequencies under 1000 Hz are compared to the closest frequency to be avoided. 1000 Hz is considered as the upper limit since it is the highest undesired frequency in the examples. Then the inversion of the sum of all the shortest frequency distances presents the frequency emphasis of the objective function: for which the exponential penalty function is given by and where f i are the natural frequencies of the current design, f j are the frequencies to be avoided that are closest to the current natural frequencies of the pipe route, f 1 is the currently lowest natural frequency, f min is the lowest allowed frequency, and n is the number of natural frequencies. The objective function gives an exponential penalty when In the examples presented in this paper f min = 125 Hz. Number one is added to the nominator in eq. 2 because the sum of minimum distances might become zero if all the current frequency distances are zero. But this is highly unlikely, since there are several frequencies to be compared and the current pipe route frequencies are highly unlikely the same as the frequencies to be avoided. The length is considered as in eq. 1. The common multi-objective optimization strategy of weighted sum is where a and b are factors which depend on how the optimization features are weighted in the objective, f r norm is the normalized frequency emphasis eq. 2 and L norm is the normalized length emphasis eq. 1. The normalization is accomplished so that the current results of these equations are divided by the results calculated with the original pipe route properties. After the normalization both emphases of the objective function have the same unit, and their relationship can be altered with a and b. For example if a = 0.8 and b = 0.2, the objective function focuses 80 % to frequency optimization and 20 % to length optimization. Hence the objective normalization, the initial guess is f (x 1 ) = 1. Figure 2 demonstrates the objective function as a contour plot for a 2D case when only the frequencies are considered in the objective function, and then, in figure 3, the length parameter is also included so that a = b = 0.5. The figures show that when the length is included the plot becomes smoother.

Strategy 2: Minimizing vibration energy
In the second method, the objective function is physically inspired by the diagonalization of the vibrating system by the eigenvalue decomposition [28] where the natural frequencies remain the same. In this transformed system, the system is described by detached, onedegree of freedom vibration systems. The response of a single system is proportional to the unknown magnitude of the excitations. Velocity response was chosen as its root mean square value correlates with the vibration energy. However, the function correctly penalizes for a natural frequency being close to an excitation frequency. The damping factor also ensures that no null division occurs when a natural frequency coincides with an excitation frequency. The frequency response of a one-degree of freedom system can be expressed as where ξ is the damping factor, f are the excitation frequencies, f i is the current natural frequency for which the response is calculated and w is is a weight factor for defining the most critical excitation frequencies, if the system's natural frequencies are in resonance.
In the example analysis w = 1. Then vibration speed calculated as the root mean square (RMS) is applied to all the responses: where n is the amount of natural frequencies calculated. A multiplicative approach was chosen for strategy: If the vibration speed approaches zero, the objective function reduces to a pipe length minimization problem. Figure 4 demonstrates the case when only RMS is considered in the objective, and then, in figure 5, also the length is considered in the same configuration as earlier. It can be seen that RMS has spread minima while, in the RMS and length combination, the minima are closer to the line connecting the end points. In both cases, when x = −100 and y = 101, the natural frequencies are f = [368.5, 371.4, 1015.6, 1015. 6, 1991.5, 1992.0] Hz.

Optimization algorithm and constraints
The considered optimization algorithm in the test cases presented in this paper is the augmented Lagrangian method (AUGLAG) [26]. AUGLAG combines the objective function with constraints by adding penalty for any violated constraints, resulting in one single function. The updated objective function is then passed to another-a so-called subsidiary optimization algorithm-which has no nonlinear constraints. The subsidiary algorithm will then perform the minimization process. If the constraints are violated by the solution the subsidiary algorithm gives, constraint penalties are increased, and the process is repeated until an optimal solution is reached.
The JuMP environment in Julia allows any number of linear/nonlinear constraints to be added to the optimization process depending on the optimization algorithm. The AUGLAG algorithm allows nonlinear equality and inequality constraints to be included. All constraints that where included in the pipe route optimization were nonlinear inequality constraints.
The first constraint forces the pipe route to stay inside the feasible domain: where dist 1 (x) calculates the minimum distance between the coordinate points and the feasible domain. The second constraint sets the minimum distance between adjacent coordinate points: where p i and p i−1 are two adjacent coordinate points and n is the total of coordinate points. The third constraint prevents the pipe colliding with itself: where dist 3 (x) calculates the minimum of the distances from a coordinate point to anywhere in the previous pipe route. The fourth constraint prevents the pipe roundings to overlap: where dist 4 (x) is the minimum distance between two adjacent roundings, and it is calculated with the equation where p 1 and p 2 are two adjacent coordinate points of the pipe route, r is the rounding radius, and α and β are the angles shown in Figure 7. Both strategies presented before apply all these constraints, but strategy 2 also has a fifth constraint to set the lowest allowed frequency: where f 1 (x) is the first natural frequency for the current pipe route. Also other rules that could be baked into the design parameters were included in the optimization process. For example, the pipe must have at least a certain straight length in both ends to make it easier to attach.

Bounding box, supports and roundings
The pipe route is supposed to stay inside a feasible spatial domain, a so-called bounding box, which is a user-defined 3D space. The designer defines the geometry for the bounding box using 3D shapes so that the pipe can avoid obstacles, and the pipe stays inside the desired space in the optimization process.
Supports or fasteners support and stiffen the pipe and attach it to the desired location. One support is defined as a beam element. One end is located on a user-defined plane, Figure 6. Mesh of a pipe that has one support. The support is a beam element the first node of which lies on a user-defined plane, and the second node is the closest pipe node in the current iteration. and the other one is the pipe node nearest to the plane in the current iteration. This kind of implementation helps the optimization algorithm by smoothing the problem. The length of the support, the attachment location on the pipe and also the starting point on the support plane will automatically vary in every iteration. The shorter the support beam is, the more it will stiffen the pipe. Thus the supports have a great impact on the optimization process. Figure 6 illustrates the discretization of the pipe route for a finite element analysis with one support plane.
In real life, the pipe route has no harsh corners. A mold with a certain size is used to create a rounding at the corners in the manufacturing process. In the optimization process, the roundings are considered in an inequality constraint that determines the minimum distance between adjacent roundings. This is accomplished since the roundings must not overlap. The rounding radius is defined by the user. The distance between two roundings is calculted in equation eq. 12.

Design process
The idea is that the designer defines a few parameters before the optimization process: start and end points for the pipe, bounding box geometry, support planes, frequencies to be avoided, and some other parameters, such as inner and outer radius of the pipe and material parameters. The optimization would then run on default settings, although customizing the optimization process itself, such as choosing the optimization algorithm and the constraints depending on the case, are not completely ruled out.
The initial guess is a route that goes in the middle of the feasible domain. The route Figure 7. Distance between two roundings. The orange line presents the real pipe route with rounded corners. Distance L must be above zero so that the pipe can be manufactured. L can be computed with trigonometry since radius r and the coordinate points p 1 and p 2 are known.
is divided into a pretermined number of coordinate points that will be optimized until the most desirable shape for the pipe route is found. The objective function is either eq. 4 or eq. 7, depending on which optimization strategy is used. Also constraints are included in the optimization process at this stage. The necessary constraints depend on the chosen strategy.
The optimization stops when an optimal solution is reached. The solution might become unfeasible if the bounding box or constraints are too restrictive. In cases where a feasible solution cannot be reached, the optimization stops when a certain stopping criteria (tolerance, time limit, etc.) is reached.

Example and results
The example pipe route optimization and the results are presented next. Both strategies presented before were tested in the example case, and the results were compared.
In the example applying strategy 1, the factors a and b of the objective function eq. 4 were both 0.5, which means the frequencies and length optimization were equally weighted. Table 1 presents the pipe properties together with the undesired frequencies in the example case. Table 2 presents the bounding box settings in the example optimization process where the bounding box is made of four cylindrical shapes. Settings were the same for both the example optimization processes. Figures 8 and 9 present the convergence of the objective functions applying both strategies. In both strategies the objective function peaks at first very high, but later it converges back. Also in strategy 2, the objective has a few lower peaks after the first peak before the curve settles. Figures 10 and 11 show the natural frequency convergence applying both strategies. In both strategies, the frequencies drop close to zero at first but then climb back to a feasible level. Figures 12 and 13 present the convergence of the constraint functions. In the figures the functions g1 and g2 refer to eq. 8 and 9, and functions g4 -g6 refer to eq. 10, 11 and 13. Constraint function g6, the constraint that defines the lowest allowed frequency, is only applied in strategy 2, since in strategy 1 the lowest allowed frequency is included in the exponential penalty function eq. 3.
All constraints should be negative in a feasible optimization result. In both strategies, two constraints (the constraint limiting the minimum distance between the pipe and the feasible domain and the constraint limiting the minimum distance between adjacent roundings) peak at the beginning into different directions, but later they converge back to the feasible region and stay below zero until the optimal solution for pipe length and natural frequencies is reached.
An optimization result, where the final pipe route design proceeds inside the given bounds and the natural frequencies are within the permissible range and all constraints are satisfied, is reached with both strategies. Figure 14 presents the optimization results in 2D with both strategies. The results are compared in table 3, showing that strategy 1 performs better on length optimization and the objective converges with less iterations than with strategy 2 which, in turn, performs better in optimizing frequencies. In strategy 2 the result pipe is more bent than in strategy 1, since it is more tender to create bends to the pipe to reach optimal frequencies at the expense of length.

Conclusions and future work
The goal of this work was to find a suitable optimization strategy to optimize high pressure fuel pipe routes in order to avoid undesired resonance. The objective was to implement the optimization procedure completely in the Julia programming environment.
JuMP was applied to the optimization process, and the frequency analysis in each iteration was performed using JuliaFEM. Both these platforms are written in Julia. Applying Julia, the whole optimization process could be completed within one programming environment.
Two strategies were presented and tested on an example case. Both strategies took into account both length and natural frequencies in the objective function. The length was optimized by minimizing the sum of distances between adjacent coordinate points of the pipe route. Frequency optimization was carried out differently in both strategies. In the first strategy, the sum of minimum distances between pipe frequencies and the excitations were minimized. In the second strategy, the vibration energy was minimized.
The example optimization case proofed that both optimization strategies work, while the 1:st strategy was faster and better on length optimization and the 2:nd strategy performed better on optimizing natural frequencies. In future other strategies and also other algorithms than AUGLAG could be tested.   Table 2. Bounding box settings in the example optimization.         Example pipe optimization result comparison presented in a 2D coordinate system. In strategy 1, the distances to excitations were minimized, while strategy 2 was based on minimizing vibration energy.