Numerical Approach and Expert Estimations of Multi-Criteria Optimization of Precision Constructions Sergey Doronin1 and Alexey Rogalev2 1 Institute of Computational Technologies SB RAS, Krasnoyarsk Branch Office, Krasnoyarsk, Russia 2 Institute of Computational Modeling SB RAS, Krasnoyarsk, Russia rogalyov@icm.krasn.ru Abstract. Structurally complex constructions with high accuracy are considered as hierarchically complex systems. For these systems, the re- quirements for the manufacture of elements, for the formation of control actions, and for instrumental monitoring of system parameters are for- mulated. The peculiarities of applying multilevel optimization methods to technical object designs are that multiparametric optimization is nec- essary at all hierarchical levels. At the same time, rigorous statements and solutions of extremal problems are difficult. The problems of opti- mization of precision constructions in the literature are considered for an extremely narrow class of technical objects. A practical (numerical) approach to multilevel optimization of structures is proposed. In the framework of this approach, multi-parameter weakly formalizable prob- lems are solved at various hierarchy levels with additional requirements for high accuracy of parameter monitoring. Multilevel optimization is applied to a wide class of technical objects. Keywords: Precision designs · Multi-level optimization · Pareto- optimal technical solutions 1 Introduction In most cases precision constructions are structurally complex technical objects. To one or several elements of these designs, special requirements are imposed on the accuracy of the manufacture of elements, the requirements are imposed on the formation of control actions, the requirements are imposed on the instru- ment control or on the estimated evaluation of the parameters of the system. The methodology of design, analysis, optimization of structurally complex designs of technical objects presupposes the consideration of them as complex hierarchical Copyright c by the paper’s authors. Copying permitted for private and academic purposes. In: S. Belim et al. (eds.): OPTA-SCL 2018, Omsk, Russia, published at http://ceur-ws.org 324 S. Doronin, A. Rogalev systems. The development of methods for optimizing hierarchical structures is one of the important problems of developing the theory of hierarchical struc- tures [5], [12]. A significant part of the results in this area is based on multi-level approaches to finding optimal solutions [11]. We note that due to the complex- ity of constructive-technological and circuit solutions, rigorous statements and solutions of extremal problems are difficult. These difficulties require the use of multicriteria approaches to optimization often [13], [21]. There is an extensive positive experience of multi-level optimization in various subject areas. In [17], hierarchical stiffened shells are examined, using various optimization methods at each hierarchical level. Multi-level optimization of the contour of the aero- dynamic profile of the wing was made in [10]. The problem of optimizing of the cost of a set of goods using a multi-level nested logit model was solved [7]. The results of the rationale for the optimal design of a geothermal power plant using a multilevel artificial neural network are presented in [1]. Multilevel op- timization of the method of field development in the context of uncertainty of geological information is reduced to justifying the distribution of wells along the field of occurrence and pressure in the wells [2]. In [15] the experience of two-level optimization of the intelligent transport system is considered. The specificity of optimization of precision designs lies in the fact that, in addition to methods of proper optimization, additional tools are used to improve the accuracy of the results. These tools are technical solutions of structures, spe- cialized software and algorithms, additional mathematical models in the subject area, as well as control systems. In [16], a precision moving platform is considered. Precision is provided by a control system including sensors and actuators. The sensor characteristics as well as the control actions produced are included in the optimization procedure. Optimization of the quantization accuracy of the second harmonic of the digital gyroscope is carried out with the help of feedback organization and development of compensating effects by the control system [20]. Providing high accuracy of photogrammetric measurements is achieved by applying additional procedures to optimize the camera’s temperature, the positioning settings of the measurement object, and illuminations [4]. In [19], an ultra-precision ion beam control device is considered for polishing the reflecting surface of mirror antennas, providing the primary optical accuracy of the mirror. Optimization of the dynamic char- acteristics of the device is realized using known dependencies between dynamic characteristics and polishing accuracy. To find the exact position of the wheel of the cutting tool for grooving the grooves [8], a global optimization procedure based on the improved particle swarm method was used. An additional model of the influence of the wheel position on the geometry of the groove is applied. Obtaining numerical estimates of the deformed state of precision structures is one of the tools used in the design calculations of precision structures In this case, the accuracy, low level of errors of such estimates is one of the most important factors for ensuring precision. In this paper, we propose an approach to multi-level optimization of precision designs based on the application of their Multi-Criteria Optimization of Precision Constructions 325 finite element models in conjunction with an additional procedure for inverse analysis of errors in numerical solutions. 2 The Concept of Multi-Level Optimization of Precision Constructions The following concept of multi-level optimization of precision constructions is proposed. This treatment of the process is developed on the assumption that the parameters of precision (a quantity, to the accuracy of determining which special requirements are imposed) are the displacements (absolute deformations) of individual elements of the system expressed in natural terms (units of length) or the integral indices (dispersion, mean-square deviation of displacements). The principal characteristic of the concept is that a special procedure for the analysis of a computational error is included in the decision contour when searching for the optimal variant, the evaluation of which is taken into account in the inequalities expressing the constraints. Then in the general case the constraints take the form [p] ≤ pcalc − perr or (1) pcalc + perr ≤ [p], (2) where p− is the precision parameter; pcalc is the calculated value of the preci- sion parameter; perr - evaluation of the computational error of pcalc ; [p]− the maximum (allowed) value of the parameter. The concept assumes 1) multi-level decomposition of the design with the selection of nodes (groups of elements) according to functional and technological features: the unit includes elements that share some function of the system; 2) selection of the parameter (parameters) of precision and the assignment of its (their) limit values; 3) formulation of local optimization problems (selection of design variables, objective functions, constraints) and obtaining Pareto-optimal solutions for each node on the basis of finite element analysis of its deformed state; 4) global optimization of the entire system in the collection; 5) evaluation and accounting for equation (1) of computational errors in determining the precision parameter at all optimization levels. Within the framework of the concept of various hierarchical levels, the so- lution of multiparameter weakly formalizable optimization problems is real- ized. When comparing constructive options from the positions of the theory of decision-making, we are dealing with the problem of ordering objects by pref- erence in the case of multicriteria choice on a finite set of alternatives. This problem is solved using the hierarchy analysis method, which is based on the use of comparative object scores by criteria using paired comparisons. 326 S. Doronin, A. Rogalev 3 The Method of Inverse Analysis of Numerical Solutions Errors Let us consider the resolving system of linear algebraic equations (SLAE) arising from finite element analysis Ax = b, (3) where A is the square sparse matrix of coefficients (stiffness matrix); b is the vector of the right side (load vector). In this paper, we propose modified methods for a posteriori error analysis [18], [3], which control the effect of rounding errors on finite element analysis of power structures. Without loss of generality, we assume that all calculations are carried out on one computer within the frame- work of the implementation of one floating-point arithmetic, and errors in the matrix of coefficients of the system arise when the matrix is represented in the computer, that is, the error of the coefficient matrix is imposed by the condition A−1 · k∆Ak < 1 Accumulation of the error in the numerical solution of algebraic equations is the total effect of roundings made at individual steps of the computational process. Therefore, for a priori estimation of the total effect of rounding errors in numerical methods of linear algebra, one can use a scheme of so-called inverse analysis. If the inverse analysis scheme is applied to evaluate the solution of SLAE, then this scheme is as follows. The uh SLA solution, calculated by the direct method M, does not satisfy the original system (3), but can be represented as an exact solution of the perturbed system (A + ∆A) u = b + ∆b (4) The quality of the direct method is estimated by the best a priori estimate [6], which can be given for the norms of the matrix ∆A and the vector ∆b. Such ”best” ∆A and ∆b are called respectively the matrix and the vector of the equivalent perturbation for the method M. If the estimates for ∆A and ∆b were obtained, then theoretically the error of the approximate solution uh can be estimated by inequality u − uh   cond(A) k∆Ak k∆bk ≤ + . kuk 1 − kA−1 k k∆Ak kAk kbk As already noted, the estimate of the norm of the inverse matrix A−1 is rarely known, and the main reason for this inequality is the ability to compare the quality of different methods. Below is a view of some typical estimates for the matrix ∆A. For methods that use orthogonal transforms and floating-point arithmetic (the coefficients of the system A, b are assumed to be real numbers) k∆AkE ≤ f (n) · kAkE · ε. (5) In the estimate (5) ε− the relative accuracy of arithmetic operations in a com- P 1 2 2 puter, kAkE = a i,j ij is the Euclidean matrix norm, f (n)− a function of Multi-Criteria Optimization of Precision Constructions 327 the form Cnk , where n− is the order of the system. The exact values of the con- stant C and the exponent k are determined by such details of the computational process as the rounding method used, and the type of operation of accumulation of scalar products. In the case of methods of Gauss type, the right-hand side of the estimate (5) also contains the factor g(A), which reflects the possibility of growth of elements of the matrix A at intermediate steps of the method com- pared to the original level (such growth is absent in orthogonal methods ). To reduce the value of g(A), different methods of selecting the leading element are used, which prevent the elements of the matrix from increasing. Let us describe the methods of a posteriori analysis of SLAE errors based on the inverse analysis of errors. Let there be given two SLAEs having one matrix of coefficients A (rigidity matrix), the first system Au = binit (6) has the vector binit of the right-hand sides, for the second system Az = bnew , (7) we form the vector of right-hand sides bnew in such a way that the resulting system has a known solution. Let us single out two variants of the formation of systems for which inverse error analysis is applicable. In the first variant, solving the system (6) by the numerical method M , we find the solution unum , then substituting this solution vector in the left-hand side of the system (6), calculate the vector bnew = Aunum . Thus, we construct a system of linear algebraic equa- tions of the form (7), the exact solution of which is the vector of the numerical solution of the system (6) unum , as can be seen by simple substitution. The norm of the difference vector kunum − znum k ≈ kunum − uk, that is, the vector of the difference unum − znum approximates the error vector of the numerical method K). Further, the second variant of the reverse error analysis changes the way of constructing the right part of the SLAE. Let us formulate the basis of this method of construction. Using estimates of the norm of the solution vector of the SLAE, it is easy to obtain the inequality kbk kuk ≥ . kAk We construct a solution z = (zk )k=1,n whose components kbk rk (A) zk = , kAk or kφ(b)k rk (A) zk = . kψ(A)k The form of the functions φ, ψ, rk is chosen so that after substituting these functions into the SLAE the norm of the right-hand side was either equal to the 328 S. Doronin, A. Rogalev norm of the right-hand side or majorized it. For example, if the components of the solution vector are put equal   max |bi | i zk =  P n , max j=1 aij i then after substitution in the system (7) the components of the vector of the right-hand side of this system take the form  Pn  a j=1 ij bnew,k = max |binit,i | z  P n . (8) i max j=1 aij i Then for the systems (6), (7), the matrices of the coefficients coincide and the norms of the right parts are equal. All test cases have confirmed estimates for errors in numerical solutions of these systems ku − unum k ≤ kz − znum k . (9) It is possible to use other types of formulas for exact solutions of SLAE with the same coefficient matrix and the right-hand side close. The proposed method is approved for the estimation of errors of nodal displacements of discrete models of a number of highly relevant technical objects. Table 3. Variants of implementation of various meth- ods for estimating the errors of a numerical solution Reverse Numerical solution unum of the sys- Error estimation of the nu- analysis tem (A + ∆A) u = binit merical solution of Wilkin- is represented as the exact solution ku − unum k ≤ son’s of the perturbed system kxk error (A + ∆A) u = b + ∆b = bnew cond(A) ≤ −1 · 1 − kA k k∆Ak  k∆Ak k∆bk · + kbk kAk Multi-Criteria Optimization of Precision Constructions 329 Inverse Numerical solution unum of the Numerically solved by one nu- analysis system merical method M of the sys- newline Au = b tem of equations Au = b and method I Az = bnew . Az = bnew , where bnew = Aunum Estimation of the error of the numerical solution of the orig- inal system k u − unum k≈k unum − znum k Reverse er- A system of linear equations is con- Solved numerically by the ror analysis structed  Az = bnew  , method M systems of equa- method II tions Au = b Az = bnew .  max |bl |  l Estimation of the error of zk =   max n aij  , P  j=1 the numerical solution of the i k original system bnew,k =   ku − unum k∞ ≤ Pn kz − znum k∞ .  j=1 aij  max |binit,l |  P n  = l  max j=1 aij  i   k = max |binit,l | · αi , l i=1,..,n The norms of the right parts of the constructed and initial systems are equal to The table was constructed (Table 3) comparing our two approaches to the inverse estimation of the error in the numerical solution of a system of equa- tions with a stiffness matrix and the inverse error analysis proposed earlier in Wilkinson’s papers. The semantic characteristic of the estimates for methods I and II is as follows. The more accurately we solved the original SLAE, the less the discrepancy, the closer the exact and approximate solution. In this case, the evaluation by method I will be better and is optimistic. The estimation by method II is the upper limit of all possible errors and, accordingly, is pessimistic. 4 Multi-Criteria Optimization of the Precision Design of a Large Mirror Antenna The practical application of the approach is considered on the example of search- ing for the optimal constructive variant of a large-sized precision reflector of a parabolic mirror antenna of terrestrial satellite communication systems. Param- eters to which high requirements for accuracy of determination and provision are made, are the maximum values of absolute deformations and the standard deviation of the reflecting surface of the mirror from the ideal paraboloid. The 330 S. Doronin, A. Rogalev Fig. 1. The distribution of total translational displacements (absolute deformations) in the reflecting segments of the first form factor connected to the frame by means of seven brackets under the action of wind pressure (40 m/s) and its own weigth Fig. 2. The distribution of total translational displacements (absolute deformations) in the reflecting segments of the second form factor connected to the frame by means of eight brackets under the action of wind pressure (40 m/s) and its own weigth Multi-Criteria Optimization of Precision Constructions 331 Fig. 3. The distribution of normal stresses (MPa) in the rod elements of the frame under the action of wind pressure (40m/s) and the own weight of the frame and mirror 332 S. Doronin, A. Rogalev Fig. 4. The first three waveforms with frequencies of 14.7, 14.7, 15.9 Hz, respectively proposed procedure for the search for a rational constructive shape of a large- sized precision reflector is multilevel (in accordance with the structure of the object: mirror optimization – level 1, frame – level 2, hubs – level 3), multi- criteria and based on obtaining Pareto-optimal technical solutions. When comparing constructive options from the positions of the theory of decision-making, we are dealing with the problem of ordering objects by pref- erence in the case of multicriteria choice on a finite set of alternatives. This problem is solved using the hierarchy analysis method, which is based on the use of comparative object scores by criteria using paired comparisons [14]. At levels 1-3, the effects of air flow, solar radiation, low and high temperature of the medium, and its own weight are considered. At level 3, the effects of sinusoidal vibration and acoustic noise are additionally taken into account. Level 1. Optimization of the mirror is to determine the best - the shapes and sizes of segments (the method for cutting a paraboloid) and their cross sections; - variants of mechanical connections of segments. We represent the space of project parameters of the mirror in the form of a set Hm = {{Vm }; {Gm }; {Lm }}, where{Vm }− is the set of variable characteristics of the mirror; {Gm }− a set of objective functions for the mirror;{Lm }− a lot of restrictions for the mirror: Hm = {{Vm }; {Gm }; {Lm }} {Vm } = {{Tm }; {Jm }; {Em }{A}} {Gm } = {fm → max, σm → min, m → min, um → min}, {Lm } = {fm ≥ [fm ]; σm ≤ [σm ]; mm ≤ [mm ]} where Tm − a set of parameters describing the topology of reflecting mirror segments; Jm − set of parameters of cross-sections of mirror segments; Em − set Multi-Criteria Optimization of Precision Constructions 333 Fig. 5. Distribution of total translational movements (mm) under the action of wind pressure (40 m / s) and the own weight of the reflector assembly Fig. 6. Distribution of stresses equivalent to Mises (MPa) under the action of wind pressure (40 m / s) and its own weight 334 S. Doronin, A. Rogalev of parameters of rigidity of constructional materials of mirror segments; A− the set of parameters describing the location of the support brackets of mirror segments. The quantities fm , σm , mm , um denote, respectively, the lowest eigenfre- quency of free oscillations of the mirror, the maximum stresses in the cross sections of the segments, the mass of the mirror, and the maximum deforma- tion (deflection) of the mirror. The values in square brackets correspond to the allowed values of the quantities. As a result of multivariate numerical analysis, the features of the mechanical behavior of reflecting mirror segments are established (Fig.1), which is the basis for optimizing their geometric and rigidity characteristics, the conditions for their mechanical interaction with the framework (the required coordinates of the support brackets). Level 2. Optimization of the frame is reduced to justifying the best - the spatial structure of the frame (number, spatial orientation of the rods and con- nections between them); - the shape and dimensions of the cross-sections of the rods. The space of design parameters of the frame is represented by a set Hs = {{Vs }, {Gs }, {Ls }}, where Vs − is the set of variable characteristics of the frame; Gs − a set of target functions for the frame; Ls − many restrictions for the skeleton: {Vs } = {{Ts }; {Js }; {Es }} {Gs } = {fs + m → max, σs → min, ms=m → min, us → min, Fb → max}, {Ls } = {fs+m ≥ [fs1 ]; σs ≤ [σs ]; ms+m ≤ [ms+m ]; Fb ≥ [Fb ]}, where {Ts }− is the set of parameters describing the topology of the framework; {Js −} a set of parameters of cross-sections of rod elements;{Es −} set of param- eters of rigidity of used structural materials of a skeleton; fs+m , σs , ms+m , us , Fb − respectively, the lowest eigenfrequency of the free vibrations of the mirror and the carcass assembly, the total mass of the carcass and mirror, the maximum deformation of the carcass and the mirror assembly, the critical strength of the buckling resistance (Euler force). The problem of optimal frame design is solved in two stages. At the first stage, the substantiation of the placement of the rods in the frame volume is carried out, i.e. ”A picture of a lattice” (”pattern of a structure”) of a skeleton. This is a weakly formalized procedure, based on the coordinates of the location of the support brackets of mirror segments obtained at the first level, Mazhid’s theorem on structural changes [9], rational reasoning and engineering intuition. A purposeful search for ways to rationally shape the frame was carried out by means of sequential introduction into the design scheme and exclusion of certain types of structural elements from it. Further, in the second stage, the geometric characteristics of the cross sections of the rod elements are justified. The ratio of objective functions and constraints by different criteria is such that the requirements for ensuring geometric stability and stability come to the fore, while the strength conditions are fulfilled with large reserves (Fig.3). Level 3. Optimization of the hub is to justify the best - the shape and di- mensions of the panels and their cross sections; - the spatial structure of the hub Multi-Criteria Optimization of Precision Constructions 335 (number, spatial orientation of the panels). The space of design parameters of the hub is represented by a set Hp = {{Vp }, {Gp }, {Lp }}, where {Vp }− a set of variable characteristics of the hub; {Gp }− a set of target functions for the hub; {Lp }− many restrictions for the hub: {Vp } = {{Tp }; {Jp }; {Ep }}; {Gp } = {fp+s+m → max, σp → min, δ → min}, {Ls } = {fp+s+m ≥ [fp+s+m ]; σp ≤ [σp ]; mp+s+m ≤ [mp+s+m ]; δ ≤ [δ]}, where {Ts }− a set of parameters that describe the hub topology; {Jp }− set of parameters of cross-sections of hub panels; {Ep }− set of parameters of rigidity of applied constructional materials; fp+s+m , σp , mp+s+m , δ− are respectively the lowest eigenfrequency of the free oscillations of the reflector assembly, the max- imum stresses in the hub elements, the total mass of the reflector, the standard deviation of the surface of the deformed mirror from the ideal paraboloid. In the analysis at the third level, the construction of the mirror and the frame is known. The main requirement for the projected hub is to ensure the maximum realization of the strength and rigidity potential inherent in the design of the mirror and the frame.As a result of optimization of the hub design, the compliance (validation) of the characteristics of the reflector in the assembly is checked by the requirements imposed on the lower frequencies of free oscillations (Fig.4), the maximum displacements (Fig.5), and the stresses (Fig.6), and other parameters, which is part of the limitations when performing the optimization procedure. 5 Conclusion A multi-level technology for the practical optimization of a precision reflector has been developed. This technology consists in the sequential substantiation of the constructive forms of the mirror, frame and hub, and is characterized by the continuity and inheritance of the system of constraints and objective functions applied at different levels. Numerical substantiation and practical optimization of rational constructive forms of the main functional subsystems of large-sized precision reflectors are performed using multivariate multivariate computational experiments. The results are satisfied by the application of the procedure of reverse error analysis. In this case, it consists in the development of several alternative versions of finite-element models and the estimation for each variant of the errors in numerical solutions. Perr . The displacement of the structures (Fig. 5) was considered as computed values of the precision parameter Pcalc .The substitution of Pcalc and Perr in (2) ensured that the permissible values of [P ] were not exceeded by moves and guaranteed achievement of structural precision. Acknowledgement. This work was done during the complex project and was financially supported by the Russian Federation Government (Ministry of Edu- cation and Science of the Russian Federation). Contract 02.G25.31.0147. 336 S. Doronin, A. Rogalev References 1. Agat, H., Arslan, O.: Optimization of district heating system aided by geothermal heat pump: a novel multistage with multilevel. Applied Thermal Engineering, ANN modelling 111, 608–623 (2017) 2. Aliyev, E., Durlofsky, L.J.: Multilevel field development optimization under un- certainty using a sequence of upscaled models. Mathematical Geosciences 49(3), 307–339 (2016) 3. Babuska, I., Rheinboldt, W.: A posteriori error analysis of finite element solutions for one dimensional problems. SIAM J. Numer. Anal. 18, 565–589 (1981) 4. Dauvin, L., Drass, H., Vanzi, L., et al.: Optimization of temperature, targets and il- lumination for high precision photogrammetric measurements. IEEE Sensors Jour- nal 18 (4), 1449–1456 (2018) 5. Dementiev, V.T., Erzin, A.I., Larin, R.M., Shamardin, Yu.V.: Problems of Opti- mization of Hierarchical Structures. Publishing house of Novosibirsk University, Novosibirsk (1996) 6. Higham, N.: Accuracy and Stability of Numerical Algorithms. SIAM, USA, Philadelphia (2002) 7. Jiang, H., Chen, R., Sun, H:. Multiproduct price optimization under the multilevel nested logit model. Annals of Operations Research 254(1-2), 131–164 (2017) 8. Li, G., Zhou, H., Jing, X., Tian, G., Li,. L. An intelligent wheel position searching algorithm for cutting tool grooves with diverse machining precision requirements. International Journal of Machine Tools and Manufacture, 149–160 (2017) 9. Majid, K.I.: Optimum Design of Structures. Newnes-Butterworth, London (1974) 10. Masters, D. A., Taylor, N. J, Rendall, C. S., Allen, C. B.: Multilevel subdivision pa- rameterization scheme for aerodynamic shape optimization. AIAA Journal 55(10), 3288–3303 (2017) 11. Migdalas, A., Pardalos, P., Varbrand, P., et al.: Multilevel Optimization: Al- gorithms and Applications. Kluwer Academic Publishers, Springer US (1998). https://doi.org/10.1007/978-1-4613-0307-7 12. Moiseyev, N.N., et al.: The Current State of the Theory of Operations Research. Science, Moscow (1979) 13. Pardalos, P., Siskos, Y., Zopounidis, C., et al.: Advances in Multicriteria Analysis. Springer-Verlag US, Springer US (1995). https://doi.org/10.1007/978-1-4757-2383- 0 14. Saaty, T.L.: The Analytic Hierarchy Process. NcGraw-Hill, New York (1980) 15. Stoilova, K., Stoilov, T., Ivamov, V.: Bi-level optimization as a tool for implemen- tation of intelligent transportation systems. Cybernetics and Information technolo- gies, 17(2), 97–105 (2017) 16. Van der Veen, G., Langelaar, M., van der Meulen, S., Laro, D., Aangenent, W., van Keulen, F.: Integrating topology optimization in precision motion system design for optimal closed-loop control performance. Mechatronics 47, 1-13 (2017) 17. Wang, B., Tian, K., Zhao, H., Hao, P., et al.: Multilevel optimization framework for hierarchical stiffened ahells accelerated by adaptive equivalent strategy. Applied Composite Materials 24(3), 575–592 (2016) 18. Wilkinson, J.: Rounding Errors in Algebraic Processes. Her Majesty’s Stationary Office, London (1963) 19. Yang, B., Xie, X., Zhou, L., Hu, H.: Design of a large five-axis ultra-precision ion beam figuring machine: structure optimization and dynamic performance analysis. International Journal of Advanced Manufactiring Technology 92, 3413–3424 (2017) Multi-Criteria Optimization of Precision Constructions 337 20. Ying, K., Wang, Z., Mao, J., Jin, Z.: Optimization of second-harmonics quantiza- tion precision for intensity modulation noise suppressing in a digital RFOG. Optics Communications 405, 114–119 (2017) 21. Zopounidis, C., Pardalos, P., et al.: Handbook of Multicriteria Analysis. Springer-Verlag Berlin Heidelberg, Springer-Verlag, Berlin Heidelberg (2010). https://doi.org/10.1007/978-3-540-92828-7