Engenharia Mecânica
Application of particle swarm optimization in inverse finite element modeling to determine the cornea ́s mechanical behavior
Aplicação da otimização por enxame de partículas na modelagem inversa de elementos finitos para determinar o comportamento mecânico da córnea
Application of particle swarm optimization in inverse finite element modeling to determine the cornea ́s mechanical behavior
Acta Scientiarum. Technology, vol. 39, no. 3, pp. 325-331, 2017
Universidade Estadual de Maringá

Received: 02 February 2016
Accepted: 20 June 2016
Abstract: Particle Swarm Optimization (PSO) was foregrounded by finite element (FE) modeling to predict the material properties of the human cornea through inverse analysis. Experimental displacements have been obtained for corneas of a donor approximately 50 years old, and loaded by intraocular pressure (IOP). FE inverse analysis based on PSO determined the material parameters of the corneas with reference to first-order, Ogden hyperelastic model. FE analysis was repeated while using the commonly-used commercial optimization software HEEDS, and the rates of the same material parameters were used to validate PSO outcome. In addition, the number of optimization iterations required for PSO and HEEDS were compared to assess the speed of conversion onto a global-optimum solution. Since PSO-based analyses produced similar results with little iteration to HEEDS inverse analyses, PSO capacity in controlling the inverse analysis process to determine the cornea material properties via finite element modeling was demonstrated.
Keywords: inverse analysis, finite element method, swarm intelligence, hyperelastic parameters, human corneas.
Resumo: A otimização por enxame de partículas (PSO) foi usada na modelagem por elementos finitos a fim de se prever as propriedades do material da córnea humana por meio de análise inversa. Deslocamentos experimentais foram obtidos por meio de amostras de córneas provenientes de doadores com aproximadamente 50 anos de idade, as quais sofreram pressão intra-ocular. A análise inversa via Elementos Finitos e baseada no PSO a partir de dados experimentais foi usada para determinar os parâmetros do material de córneas, utilizando como referência o modelo hiperelástico de Ogden de primeira ordem. A análise inversa por Elementos Finitos também foi executada por um programa comercial de otimização e foram encontrados os mesmos valores dos parâmetros do material, o que serviu para validar o resultado do PSO. Além disso, o número de iterações necessárias para a otimização, tanto para o PSO quanto para o programa comercial, foi comparado para avaliar a velocidade de conversão para uma solução global ótima. O estudo mostrou que as análises baseadas no PSO produziram resultados semelhantes com menos iterações para a análise inversa, o que demonstrou o potencial do PSO em controlar o processo de análise inversa para determinar as propriedades mecânicas da córnea por meio de modelagem via Elementos Finitos.
Palavras-chave: análise inversa, método dos elementos finitos, inteligência de enxame, parâmetros hiperelásticos, córneas humanas.
Introduction
Particle Swarm Optimization (PSO), based on socio-psychological principles inspired by swarm intelligence and social behavior, has been widely applied to engineering (Marwala, 2010). PSO method was first forwarded in 1995 (Poli, Kennedy & Blackwell, 2007) and soon found applications in several areas (Liu, Liu, & Duan, 2007; Liu, Liu, & Cartes, 2008; Arumugam & Rao, 2008; Berlinet & Roland, 2008; Thakker, Patil, & Anil, 2009; Santos, Martins, & Santos, 2012).
Marwala, Boulkaibet, & Adhikari (2016) combined PSO with finite element (FE) modeling to get the best fit for numerically-predicted load-deformation behavior with experimental data obtained for a simple beam. More demanding structural mechanics applications were undertaken in later studies which demonstrated PSO reliability. Or rather, it may effectively be combined with FE methods in inverse modeling problems (Tang, Chen, & Peng, 2009; Oliveira eta al., 2011; Chen, Tang, Ge, An, & Guo, 2013; Reutlinger, Bürki, Brandejsky, Ebert, & Büchler, 2014).
Other optimization methods such as Genetic Algorithms (GA) (Shabani & Mazahery, 2013) have been employed with FE (Song , Yang, Yu, & Kang, 2013) or with both (Herath, Natarajan, Prusty, & John, 2014) in optimization problems. PSO has also been compared with the Artificial Bee Colony (ABC) algorithm, demonstrating higher accuracy and better reliability in inverse modeling applications involving FE analysis (Bardsiri, Jawawi, Hashim, & Khatibi, 2013).
The determination of the mechanical properties of human organs, including the eye, using FE inverse modeling, has been quite common, specifically within Biomechanics (Wittek et al., 2013; Perez, Morris, Hart, & Liu, 2013). Research in ocular Biomechanics, focusing mainly on the cornea’s and sclera’s mechanical properties (Elsheikh et al., 2007; Elsheikh, Ross, Alhasso, & Rama, 2009; Elsheikh, Geraghty, Rama, Campanelli, & Meek, 2010; Abyaneh, Wildman, Ashcroft, & Ruiz, 2013) has progressed significantly, although authors used optimization techniques different from PSO.
Current study assessed a simple inverse-modeling implementation of PSO in the application of Biomechanics and compared to the commercial software HEEDS. Commercial software HEEDS has a strong record for the optimization of solutions in mechanics problems by employing a hybrid and adaptive algorithm. HEEDS employs multiple search strategies and adapts to the problem. In fact, it requires significantly fewer model evaluations to identify optimized solution for the first time.
When compared to commercial package HEEDS, PSO is simple to implement and may be coded by research groups with limited resources and fine-tuned to suit their specific needs. The comparisons presented in current paper are related to an issue in which the mechanical properties of the corneas of the 50-year-old donor are estimated in an inverse modeling exercise.
Material and methods
Objective function definition
In previous studies, cornea samples were tested under inflation conditions (Elsheikh et al, 2007). Figure 01 shows the hyperelastic behavior of human corneas under increasing posterior pressure.

Based on experimental data from human corneas (Elsheikh et al., 2008), the objective function has been obtained and represented by a five-order polynomial equation 1:

where:
x is the intraocular pressure and
y is the displacement.
Equation 1 was used for inverse modeling procedures for PSO and HEEDS methods to compare the rates of material parameters obtained and the speed of the loading process in general.
Finite element model
An FE model of a human cornea involving 12,168 quadratic, triangular, prismatic, hybrid elements, arranged in 3 radial layers, was generated and used in current study (Figure 2). The model was maintained along the edge nodes to simulate the conditions experienced by the donor´s corneas in the experimental program. The model included fluid elements representing the aqueous humor of the anterior chamber, taking an internal pressure of up to 45 mm Hg. The model had an anterior geometry that was sinusoidal, with a central radius of curvature of the anterior surface, Rant = 7.8 mm, a shape factor, p = 0.82, central corneal thickness, CCT = 0.545 mm and a peripheral corneal thickness, PCT = 0.695 mm. The above rates represented the average results arrived at experimentally for the tested corneas and made possible the average values in the numerical modeling conducted in current study.

Inverse modeling
For simplification, corneal tissue was presumed to have an isotropic and hyperelastic behavior, represented by Ogden’s material model (Ogden, Saccomandi, & Sgura, 2004):

where:
W is the strain energy;
µi and ai are material parameters;
i = 1 to equation order N;
are the stretches in the Cartesian system,
j = x, y, z.
The material parameters to be defined from an inverse modeling procedure involve the finite element analysis, in which the input includes the cornea’s geometry and the pressure-deformation behavior, while the expected output is the material behavior represented by µi and ai. The solution comprises an optimization process that identifies the global minimum of the difference between the experimentally-obtained and the numerically-predicted, load-deformation behavior at various points on the corneal surface. This process commonly utilizes commercial tools such as HEEDS or, as in this case, optimization methods such as PSO. In current study, the order of the strain energy material model was assumed as 1 since, in previous studies, this order was adequate to describe experimental behavior of biological tissue (Elsheikh, 2010). Limiting the order to 1 meant the optimizing of rates of only two parameters to reduce cost of analysis.
Particle Swarm Optimization
PSO is an optimization technique based on populations with m particles (m individuals) that evolve within the hyperspace defined by the design´s variable bounds following some random criteria towards the particle with the best performance (usually the particle that is closest to the optimum point) (Barbieri, Barbieri, & Lima, 2015). The PSO algorithm is an alternative which describes each particle (material parameters) by its vector position in the search space as:

where:
w factor is the inertia weight;
v (t) is the particle’s velocity at time t;
x (t) is the particle's current position at time t;
c1 and c2 are weight constants;
p (t) is the particle’s best position;
g (t) is the best known position found by any particle in the swarm;
r1 and r2 are random numbers resulting from an algorithm choice that normally vary between 0 and 1 (Marwala, 2005).
By applying Equation (3), the new velocity, v (t+1), is determined to compute the new particle’s position, x (t+1) as shown in Equation (4).

Since in the Equation (2), µi and ai are the material parameters to be defined, a finite element solver Abaqus/Standard 6.12 (Dassault Systèmes Simulia Corp., Rhode Island, USA) was used with PSO. For the numerical analysis, Abaqus was used with a PSO code written in Visual Basic (VB) to provide the best fitness between FE results and experimental data. A Python code was written to extract the displacement curve in z direction of the apical node of the cornea, using FE. The PSO code used this target curve to calculate the cost function as the sum of squared error (SSE) between FE displacements from the cornea apex and experimental data (Equation 1). The algorithm takes the discrete data points of a target curve, fit the curve to the data and compare the equation of the curve with FE results. A summary of the PSO process, specific for this kind of problem, is shown in a flow chart (Figure 3).

Boundary rates were defined to constrain the search space between maximum and minimum rates for each PSO iteration. This limitation is required since rates that exceed maximum and minimum boundaries lead to unrealistic material parameters and convergence issue for the FE solver.
Ten particles, parameters c1 and c2 equal to 1.494, w equal to 0.729 and population size equal to 30 were set based on previous studies. Rates were recommended by Clerc (2013) and tested by Tuppadung & Kurutach (2011). Trelea (2003) used them in a large number of simulation experiments for convergence and sensitivity analysis, with good results. In case of non-converged simulation in the proposed algorithm, a specific line was added in the PSO code to discard the iteration in which errors occurred.
Results and discussion
Material parameters obtained from PSO have been compared with those by HEEDS. Analysis was performed with six random boundaries situations (Table 1), justified by the variation of the limits (maximum and minimum) of the parameters to certify the stability of the proposed methodology.
Figure 4 shows graphs from PSO and HEEDS results for the first boundary condition. Similarly, the other boundary conditions were also evaluated (Table 1).

Figure 4 shows fewer iterations for PSO when compared with those for HEEDS. It should be underscored that each iteration displayed on the graphs represents the position of ten particles.
| Boundary | µi (min-max) | ai (min-max) | µi results (PSO) | µi results (HEEDS) | ai results (PSO) | ai results (HEEDS) | SSE results (PSO) | SSE results (HEEDS) |
| 1 | 0.00001-1 | 20-400 | 0.054 | 0.055 | 110.45 | 109.3 | 0.0016 | 0.0017 |
| 2 | 0.00001-1 | 10-200 | 0.054 | 0.055 | 110.43 | 109.75 | 0.0016 | 0.0016 |
| 3 | 0.02-0.2 | 20-150 | 0.054 | 0.054 | 110.46 | 110.35 | 0.0016 | 0.0016 |
| 4 | 0.01-0.1 | 60-240 | 0.054 | 0.054 | 110.45 | 110.40 | 0.0016 | 0.0016 |
| 5 | 0.025-0.075 | 60-150 | 0.054 | 0.054 | 110.45 | 110.40 | 0.0016 | 0.0016 |
| 6 | 0.001-0.8 | 40-300 | 0.054 | 0.053 | 110.45 | 111.50 | 0.0016 | 0.0017 |
This means that the number of FE simulations tested was the same for both methods. Taking the above into consideration coupled to HEEDS rates in ascending order, PSO achieved the optimized rate first. Figure 5 compares SSE rates for the two methods with boundary 2 (Table 1) as example.

Figure 5 revealed that, by SSE values, PSO converged after 400 simulations, or rather, prior to HEEDS which converged after 500 simulations. However, HEEDS ran simulations in parallel using a variety of optimization algorithms such as Multi-objective SHERPA, Genetic algorithm, Quadratic programming, Simulated annealing, Response surface, Multi-start local search and Nelder-Mead simplex. In other words, each method took advantage of the best attributes and solutions found from other parallel searches (Abedrabboa, Worswicka, Mayer, & Riemsdijkc, 2009).
Regarding to PSO convergence, it may be observed that running the same boundary condition more than once, all rates after each iteration remained the same even in different PC. Consequently, six different boundaries were implemented to test the convergence. Figure 6 represents the convergence of parameters and for the six boundaries tested (from Table 1).
Even though the benefits of PSO extend to other cases from a specific case (healthy 53-year-old human corneas modeled with Ogden constitutive law and loaded to 45 mm Hg), it may be noted that PSO is more stable than HEEDS. On the other hand, HEEDS software provided a graphical user-friendly interface easy to use and did not need any programming knowledge, or rather, assets benefitted HEEDS.
Ogden model is highly useful in the cornea for material characterization (Elsheikh et al., 2007; Lago et al., 2015), even though it has certain disadvantages, such as isotropy. One option is the use of PSO optimization scheme for other models (e.g., Holzapfel model, Fung model) not used here. It is common knowledge that the hyperelastic material model is not the only choice to represent the corneal mechanical response, whilst viscoelasticity (Lombardo, Serrao, Rosati, & Lombardo, 2014) and anisotropy (Whitford, Studerb, Booteb, Meekc, & Elsheikh, 2015) should be considered in future research works with PSO.

Conclusion
Current research focused on a novel application for Particle swarm optimization coupled to finite element modeling to predict the material properties of the human cornea through inverse analysis. One advantage of the PSO code is that it does not need the initial guess since it chooses random rates between the boundary values established by the user. However, HEEDS runs simulations in parallel with several different optimization algorithms. At the end of all simulations, it chooses the best solution from the parallel searches. PSO depends on the previous particle’s position, generating a limitation in the way that the proposed application has been implemented in terms of time consuming.
It is known that results´ stability with inverse analysis used is also an advantage. Stability is mostly a function of the change in results from different initial guesses. In current case, results demonstrated no concern with PSO with regard to stability since six different boundary conditions have been used with the same results after optimization.
PSO has also the advantage that it converged first, considering the same number of iterations as HEEDS. Further, PSO is an open source code whereas HEEDS is a commercial software. Taking all the above aspects into consideration, the capacity of PSO in controlling the inverse analysis process to determine the cornea´s material properties via finite element modeling should be underscored.
Acknowledgements
Current research was funded CAPES under Grant 2623-13-7.
References
Abedrabboa, N., Worswicka, M., Mayer R., & Riemsdijkc, I.V. (2009). Optimization methods for the tube hydroforming process applied to advanced high-strength steels with experimental verification. Journal of Materials Processing Technology, 209(1), 110-123.
Abyaneh, M. H., Wildman, R. D., Ashcroft, I. A., & Ruiz, P. D. (2013). A hybrid approach to determining cornea mechanical properties in vivo using a combination of nano-indentation and inverse finite element analysis. Journal of the Mechanical Behavior of Biomedical Materials, 27(1), 239-248.
Arumugam, M. S., & Rao, M. V. C. (2008). Molecular docking with multi-objective particle swarm optimization. Applied Soft Computing, 8(1), 666-675.
Barbieri, R., Barbieri, N., & Lima, K. F. (2015). Some applications of the PSO for optimization of acoustic filters. Applied Acoustics, 89(1), 62-70.
Bardsiri, V. K., Jawawi, D. N. A., Hashim, S. Z. M., & Khatibi, E. (2013). A PSO-based model to increase the accuracy of software development effort estimation. Software Quality Journal, 21(3), 501-526.
Berlinet A, & Roland C. (2008). A novel particle swarm optimization algorithm for permutation flow-shop scheduling to minimize makespan. Chaos, Solitons & Fractals, 35(5), 851-861.
Chen, J., Tang, Y., Ge, R., An, Q., & Guo, X. (2013). Reliability design optimization of composite structures based on PSO together with FEA. Chinese Journal of Aeronautics, 26(2), 343-349.
Clerc, M. (2013). Particle Swarm Optimization. New York, NY: John Wiley & Sons.
Elsheikh, A. (2010). Understanding corneal biomechanics through experimental assessment and numerical simulation. In J. H. Levy (Ed.). Biomechanics: principles, trends and applications (p. 57-110). New York, NY: Nova Science Publishers Incorporated.
Elsheikh, A., Alhasso, D., & Rama, P. (2008). Biomechanical properties of human and porcine corneas. Experimental Eye Research, 86(5), 783-790.
Elsheikh, A., Geraghty, B., Rama, P., Campanelli, M., & Meek, K.M. (2010). Characterization of age-related variation in corneal biomechanical properties. Journal of the Royal Society Interface, 7(51), 1475-1485.
Elsheikh, A., Ross, S., Alhasso, D., & Rama, P. (2009). Numerical study of the effect of corneal layered structure on ocular biomechanics. Current Eye Research, 34(1), 26-35.
Elsheikh, A., Wang, D., Brown, M., Rama, P., Campanelli, M., & Pye, D. (2007). Assessment of corneal biomechanical properties and their variation with age. Current Eye Research, 32(1), 11-19.
Herath, M. T., Natarajan, S., Prusty, B. G., & John, N. S. (2014). Smoothed finite element and genetic algorithm based optimization for shape adaptive composite marine propellers. Composite Structures, 109(1), 189-197.
Lago, M. A., Rupérez, M. J., Martínez-Martínez, F., Monserrat, C., Larra, E., Güell, J. L., & Peris-Martínez, C. (2015). A new methodology for the in vivo estimation of the elastic constants that characterize the patient-specific biomechanical behavior of the human cornea. Journal of Biomechanics, 48(1), 38-43.
Liu, L., Liu, W., & Cartes, D. A. (2008). Particle swarm optimization-based parameter identification applied to permanent magnet synchronous motors. Engineering Applications of Artificial Intelligence, 88(21-22), 1092-1100.
Liu, X., Liu, H., & Duan, H. (2007). Particle swarm optimization based on dynamic niche technology with applications to conceptual design. Advances in Engineering Software, 38(10), 668-676.
Lombardo, G., Serrao, S., Rosati, M., & Lombardo, M. (2014). Analysis of the Viscoelastic Properties of the Human Cornea Using Scheimpflug Imaging in Inflation Experiment of Eye Globes. Plos One, 9(11), e112169.
Marwala, T. (2005). Finite element model updating using particle swarm optimization. International Journal of Engineering Simulation, 6(2), 25-30.
Marwala, T. (2010). Finite-element-model updating using computional intelligence techniques: applications to structural dynamics. London, UK: Springer-Verlag London.
Marwala, T., Boulkaibet, I., & Adhikari, S. (2016) Model selection in finite element model updating, in probabilistic finite element model updating using bayesian statistics: applications to aeronautical and mechanical engineering. Chichester, UK: John Wiley & Sons, Ltd.
Ogden, R. W., Saccomandi, G., & Sgura, I. (2004). Fitting hyperelastic models to experimental data. Computational Mechanics, 34(6), 484-502.
Oliveira, M. E., Hasler, C. C., Zheng, G., Studer, D., Schneider, J., & Büchler, P. (2011). A multi-criteria decision support for optimal instrumentation in scoliosis spine surgery. Structural and Multidisciplinary Optimization, 45(6), 917-929.
Perez, B., Morris, H., Hart, R., & Liu, J. (2013). Finite element modeling of the viscoelastic responses of the eye during microvolumetricchanges. Journal of Biomedical Science and Engineering, 6(12A), 29-37.
Poli, R., Kennedy, J., & Blackwell, T. (2007). Particle swarm optimization: An overview. Swarm Intelligence, 1(1), 33-57.
Reutlinger, C., Bürki, A., Brandejsky, V., Ebert, L., & Büchler, P. (2014). Specimen specific parameter identification of ovine lumbar intervertebral discs: On the influence of fibre-matrix and fibre-fibre shear interactions. Journal of the Mechanical Behavior of Biomedical Materials, 30(1), 279-289.
Santos, R. P. B., Martins, C. H., & Santos, F. L. (2012). Simplified particle swarm optimization algorithm. Acta Scientiarum. Technology, 34(1), 21-25.
Shabani, M. O, & Mazahery. A. (2013). Application of GA to optimize the process conditions of Al Matrix nano-composites. Composites Part B: Engineering, 45(1), 185-191.
Song ,J. H., Yang, S. M., Yu, H. S., & Kang, H. Y. (2013). Fatigue strength under optimal conditions of spot weldment using GA and FEM. Applied Mechanics and Materials, 284-287(1), 123-126.
Tang, Y., Chen, J., & Peng, W. (2009). Probabilistic optimization of laminated composites considering both ply failure and delamination based on PSO and FEM. Tsinghua Science and Technology, 14(S2), 89-93.
Thakker, R. A., Patil, M. B., & Anil, K. G. (2009). Parameter extraction for PSP MOSFET model using hierarchical particle swarm optimization. Engineering Applications of Artificial Intelligence, 22(2), 317-328.
Trelea, I.C. (2003). The particle swarm optimization algorithm: convergence analysis and parameter selection. Information Processing Letters, 85(6), 317-325.
Tuppadung, Y., & Kurutach, W. (2011). Comparing nonlinear inertia weights and constriction factors in particle swarm optimization, International Journal of Knowledge-based and Intelligent Engineering Systems, 15(2), 65-70.
Whitford, C., Studerb, H., Booteb, C., Meekc, K. M., & Elsheikh, A. (2015). Biomechanical model of the human cornea: considering shear stiffness and regional variation of collagen anisotropy and density. Journal of the Mechanical Behavior of Biomedical Materials, 42(1), 76-87.
Wittek, A., Karatolios, K., Bihari, P., Schmitz-Rixen, T., Moosdorf, R., Vogt, S., & Blasé, C. (2013). In vivo determination of elastic properties of the human aorta based on 4D ultrasound data. Journal of the Mechanical Behavior of Biomedical Materials, 27(1), 167-183.
Author notes
ricardorm@deg.ufla.br