# Non-iterative method of stress projection in beam elements for stent deployment simulation

Arthur Boldron, Benoît Blaysat, Thomas Jailin, Laurent Sarry

► To cite this version:

Arthur Boldron, Benoît Blaysat, Thomas Jailin, Laurent Sarry. Non-iterative method of stress projection in beam elements for stent deployment simulation. 26e Congrès Français de Mécanique, Laboratoire d'Etude des Microstructures et de Mécanique des Matériaux (LEM3 - UMR CNRS 7239), Aug 2025, Metz, France. **hal-05313088**

HAL Id: hal-05313088

https://hal.science/hal-05313088v1

Submitted on 14 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Non-iterative method of stress projection in beam elements for stent deployment simulation

A. BOLDRON<sup>a</sup> @, B. BLAYSAT<sup>b c</sup>, T. JAILIN<sup>b</sup>, L. SARRY<sup>b</sup>

a. Institut Pascal, Université Clermont Auvergne, CHU Clermont-Ferrand, Clermont Auvergne INP, CNRS, F-63000 Clermont-Ferrand, France

b. Institut Pascal, Université Clermont Auvergne, Clermont Auvergne INP, CNRS, F-63000 Clermont-Ferrand, France

c. Institut Universitaire de France, Paris, France

@. arthur.boldron@uca.fr

## Résumé :

Dans ces travaux nous visons à améliorer l'efficacité du déploiement d'un stent (ou endoprothèse vasculaire) à l'intérieur d'une artère pathologique en en développant une simulation numérique de ce déploiement qui soit intégrable dans une routine clinique. Pour cela il a été proposé un maillage du stent par des éléments poutres. Ceux-ci ont permis de réduire significativement les temps de calculs de la simulation mais restent prohibitifs. Afin de diminuer encore les temps de calculs tout en gardant une précision acceptable pour cette application, nous proposons un modèle poutre basé sur le formalisme des torseurs. Nous proposons également d'adapter une nouvelle méthode non-itérative de calcul des contraintes en comportement plastique, développée pour réduire les temps de calculs sur des éléments 2D et 3D, à notre modèle poutre.

## Abstract :

In this work we aim to improve the efficiency of stent deployment inside a pathological artery by developing a numerical simulation of this deployment that can be integrated into a clinical routine. For this purpose, a mesh of the stent by beam elements was proposed. These have allowed to significantly reduce the simulation calculation times but remain prohibitive. In order to further reduce the calculation times while maintaining an acceptable precision for this application, we propose a beam model based on a screw approach. We also propose to adapt a new non-iterative method for stress computations with plastic behavior, developed to reduce calculation times on 2D and 3D elements, to our beam model.

## Mots clefs : Computational plasticity ; Stent ; Beam elements ; Finite element analysis

## 1 Introduction

Coronary artery disease, or ischaemic artery disease, is the formation of lipid plaques hindering blood circulation in the coronary arteries. According to the World Heart Report, the number of deaths from
cardiovascular disease was 18.6 million (a third of all global deaths) including 9.1 million from coronary heart disease. This number has increased over the past 30 years as the population has grown and aged [4].

To treat coronary artery disease, a procedure called percutaneous coronary intervention (PCI) is a popular option consisting of deploying a stent (endovascular metallic prosthesis) inside a pathological artery. This method offers the major advantage of being minimally invasive and non-surgical compared with coronary bypass surgery. To guarantee the efficacy of stenting, it is essential to ensure that the stent is correctly positioned and that the mechanical stress it exerts on the artery after deployment is not excessive. However, this verification procedure can be complex, particularly in view of the tortuous geometry of coronary arteries and the inhomogeneities of wall properties.

To enhance the success rate of stent deployment, we have developed a digital simulation process as part of the ANR VSTENT project (ANR-22-CE19-0014). This process is designed to be integrated into clinical routine, offering the necessary precision for clinical use while being fast enough not to delay PCI. Since the stent has a filiform geometry, a model based on a mesh of beam elements was developed by Krewcun as part of the VSTENT project [1][2]. A Timoshenko based beam approach is used [5], coupled with a von Mises plasticity model with isotropic-kinematic linear hardening. This enables maintaining satisfactory accuracy with limited calculation times. This has then been validated experimentally on in vitro 3D micro-CT scans.

Although the use of beam elements enables obtaining limited calculations times, it is not yet possible to achieve real-time numerical simulations under realistic interventional conditions. We therefore propose to modify the proposed beam model based on a screw approach to reduce the number of integration points per element. Furthermore, we implement in our model a new method presented by Yoon and Barlat [3] to solve the plasticity problem without the need of a costly iterative procedure, meaning the common radial return procedure. The aim of this work is therefore to make the combination of these two strategies usable to further reduce calculation time.

## 2 Beam model used

To improve the time efficiency of the model proposed by Krewcun [1], different strategies are used. The first is to lighten the calculation of the stiffness matrix of the beam elements by using the properties of the beam model. Indeed, in Krewcun's model, the integration to obtain the stiffness matrix is performed with 27 Gauss points (3 points per direction) as illustrated by Figure 1, this was to preserve 3-dimensional behavior and avoid neglecting certain terms in the stress tensor as we can have with the assumption of plane stress for instance. We propose to reduce the number of Gauss points for integration to 3, only along the beam axis. To do this, we describe the behavior with a screw approach and generalized section forces and moments.
<figure><img src="image_2.png" /><figcaption>FIGURE 1 – Illustrative diagram of integration points used in each beam element as presented in [1]. $\mathbf{x}_A$ and $\mathbf{x}_B$ are the nodes.</figcaption></figure>

This involves replacing the tensor formulation of stresses and strains with a vector formulation involving the resulting forces and moments after cross-section integration of the stress tensor. The constitutive elasto-plastic tensor law between stress and strain increment ($d\sigma$ and $d\varepsilon$ respectively) on the basis of the associated (normality) flow rule is usually defined by :

$$d\sigma = \mathbf{C}^{\text{ep}} : d\varepsilon \qquad (1)$$

Where $\mathbf{C}^{\text{ep}}$ is the elastoplastic modulus.

For a beam of axial direction $x$, and length $L$, consider $\xi_i = x/L$, the dimensionless axial coordinate of the Gauss point and $\omega_i$ its weight, we can then establish the displacement of an arbitrary point on the beam axis using the shape function matrix $\mathbf{A}$ to obtain $\mathbf{u}(\xi) = \mathbf{A}(\xi) \cdot \mathbf{u}_{\text{nodes}}$. $\mathbf{u}_{\text{nodes}}$ containing the 6 degrees of freedom of the two nodes of the beam element. We can then express the deformation vector as a function of nodal displacements, $\epsilon(\xi) = \mathbf{B}(\xi) \cdot \mathbf{u}_{\text{nodes}}$. Details of shape function expressions are given in [6].

In our model, we have instead a relation between $\mathbf{t}^{\text{beam}}$ the generalized section force and moment vector and $\epsilon_i$, $\chi_i$ the strain components :

$$\mathbf{t}^{\text{beam}} = \begin{pmatrix} f_x \\ f_y \\ f_z \\ M_x \\ M_{fy} \\ M_{fz} \end{pmatrix} = \mathbb{H}^{\text{beam}} \cdot \begin{pmatrix} \epsilon_x = \frac{\partial u_x}{\partial x} \\ \epsilon_y = \frac{\partial u_y}{\partial x} - \theta_z \\ \epsilon_z = \frac{\partial u_z}{\partial x} + \theta_y \\ \chi_x = \frac{\partial \theta_x}{\partial x} \\ \chi_y = \frac{\partial \theta_y}{\partial x} \\ \chi_z = \frac{\partial \theta_z}{\partial x} \end{pmatrix} \qquad (2)$$

$\mathbb{H}^{\text{beam}}$ is the constitutive law matrix, in linear elasticity, is given by :

$$\mathbb{H}^{\text{beam}} = \mathbb{H}^{\text{beam},e} = \text{Diag} \left[ ES \ k_ySG \ k_ySG \ GC \ EI_y \ EI_z \right]$$

Where $\mathbb{H}^{\text{beam},e}$ is the elastic constitutive law matrix. Here $G = \frac{E}{2(1+\nu)}$, where E is the Young's modulus and $\nu$ the Poisson's ratio; S is the cross-sectional area of the beam, and L its length; $I_y$ and $I_z$ are second moments of area; $C = I_y + I_z$; $\Phi_y = \frac{12EI_y}{k_zSGl^2}$ and $\Phi_z = \frac{12EI_z}{k_ySGl^2}$ are shear deformation parameters with
$k_y$ and $k_z S$ the shear area in the the y and z directions. Numerically, the stiffness matrix of the beam element $K_e$ is calculated by integration along the beam axis with 3 Gauss points only.

$$K_e = \frac{L}{2} \sum_{i=1}^{3} \omega_i \mathbb{B} \left( \frac{1 + \xi_i}{2} L \right)^T \cdot H_i^{\text{beam}} \cdot \mathbb{B} \left( \frac{1 + \xi_i}{2} L \right) \quad (3)$$

## 3 Stress computation

The plastic beam model presented by Krewcun has been implemented in the SOFA software. SOFA (Simulation Open Framework Architecture) is routinely used in numerous applications in the medical field. It is an open-source software offering a generic physical simulation platform. It is designed for real-time simulations, as needed when simulation is a decision tool for medical applications such as PCI. Time efficiency of this implementation relies on the single iteration used to solve the non linear part of the mechanical problem. This presents a potentially large approximation when the simulation includes non-linear plasticity models. Indeed, conventional stress calculation methods, such as radial return, require an iterative Newton-Raphson process to handle these non-linearities, which can be computationally expensive.

To address this problem, Yoon and Barlat proposed in [3] a procedure to solve the plasticity problem with a non-iterative stress projection method (NSPM) for anisotropic strain hardening. This method is illustrated by Figure 2 included in [3]. An intermediate stress increment $\Delta\sigma_{n+\frac{1}{2}}$ is first calculated using the elastoplastic constitutive law calculated at the previous increment $n$, $C_n^{\text{ep}}$ and the strain increment $\Delta\varepsilon_{n+1}$. The integrated stress must satisfy the yield condition, generally defined as :

$$F(\mathbf{r}, \boldsymbol{\sigma}, \bar{\varepsilon}) = \bar{\sigma}(\mathbf{r}, \boldsymbol{\sigma}) - H(\bar{\varepsilon}) = 0 \quad (4)$$

$\bar{\sigma}$ is the effective stress which depends on the stress tensor $\boldsymbol{\sigma}$ and internal hardening state variables $\mathbf{r}$. $H$ is the flow stress depending on the effective plastic strain $\bar{\varepsilon}$. After updating $\mathbf{r}$ and $\bar{\varepsilon}$ using $\Delta\sigma_{n+\frac{1}{2}}$, the intermediate stress is projected onto the strain-hardening flow stress $H(\bar{\varepsilon}_{n+1})$, to fulfill the yield criterion (4). More details about this method will be given during the presentation.
<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Schematic view of the non-iterative stress projection method as presented in [3]. The stress increment is determined from the elastoplastic constitutive law. The integrated stress is corrected for the yield condition satisfaction</figcaption></figure>

Simulations used in [3] as a proof-of-concept were performed using ABAQUS shell and brick elements. It showed accurate, robust and efficient stress integration of anisotropic strain-hardening models such as Yoshida Uemori's (YU) nonlinear isotropic-kinematic strain-hardening model with a non-quadratic anisotropic yield function developed by Yoon and Barlat. Computation times were studied with a large-scale metal forming simulation, where NSPM halved them. It enables reducing the computation times by a factor two on large-scale metal forming simulations.

In this work, we propose to combine this method with beam elements and the screw approach for use in stent deployment simulations.

Yoon and Barlat developed the method for 2D and 3D elements and released the code as an Abaqus UMAT (User Material) routine. This UMAT calculates stress tensor $\sigma$, state variables and the constitutive elastoplastic law matrix $C^{ep}$, so the UMAT code also features a tensor formalism. We therefore need to adapt this method to our beam element formalism which reduces the number of nodes and integration points for the element. To do this, we implemented a User Element (UEL) on ABAQUS that calculates the element's stiffness matrix according to the equation (3), as well as the vector of internal forces according to the following equation.

$$F_e = \frac{L}{2} \sum_{i=1}^{3} \omega_i F_i = \frac{L}{2} \sum_{i=1}^{3} \omega_i \mathbb{B}^T \left( L \frac{1 + \xi_i}{2} \right) \cdot t_i^{beam} \quad (5)$$

Based on the UMAT proposed by Yoon and Barlat, the screw beam approach was implemented to compute non-iteratively $t^{beam}$ and $\mathbb{H}^{beam}$. The yield condition (6) is similar to (4), replacing the Cauchy stress tensor $\sigma$ by $t^{beam}$ and using adapted formulas for $\bar{\sigma}$ and $H$. We use the formulas presented in the work by Verwaerde [6].
$$F(\mathbf{r}, \mathbf{t}^{\text{beam}}, \bar{\varepsilon}) = \bar{\sigma}(\mathbf{r}, \mathbf{t}^{\text{beam}}) - H(\bar{\varepsilon}) = 0 \quad (6)$$

For instance, the yield condition used by Verwaerde [6] used state variables to take into account the fact that the beam cross-section gradually transits from full elastic to full plastic behavior under bending stress. For the sake of conciseness, the full equations of this yield criterion are not detailed here but can be found in [6].

# 4 Implementation of the model on ABAQUS, first feedback

To develop and validate the method we modelled a bending test of a cantilever beam as shown below Figure 3. The implementation was first validated with a linear elastic case using the analytical solution of the problem.

<figure><img src="image_4.png" /><figcaption>FIGURE 3 – Simulation used to test our UEL, the AB beam can be divided into several elements to have significant computation times to quantify</figcaption></figure>

The Yoon-Barlat's NSPM method in plasticity was then implemented for our beam elements. As in their paper, the results are compared with a fully-implicit iterative stress integration method : the closest point projection method (CPPM). It shows a good implementation of the beam approach.

These first results showed that the non-iterative method reduces the user time on ABAQUS by a third, where the resolution of a load increment required up to 6 iterations with the CPPM method.

# 5 Conclusion

In order to increase the success of percutaneous coronary intervention, a popular method of treating coronary artery disease, a simulation of stent deployment was developed with the aim of making it part of a clinical routine. For this purpose, a model based on a beam element mesh was developed, which significantly reduced the computation times that were still prohibitive for integration into clinical routine. In order to further reduce calculation times, while maintaining the necessary accuracy, a new formulation of beam elements with a screw approach was proposed to reduce the number of integration points per element. This approach was combined with a non-iterative method for calculating stresses. Initial tests show that this non-iterative method should be applicable to our beam elements and may reduce calculation times.

However, we have so far have difficulty characterizing the accuracy, in fact the simulation converges with the assumption of small displacements, the addition of the possibility of large displacements with in particular the updating of the local reference of each element as in [7] poses convergence problems
still under study. Once this beam model has been correctly produced on ABAQUS, we will implement
this method on SOFA on the stent simulations produced by Krewcun in order to improve them.
At the same time, we will implement the non-iterative stress calculation method on SOFA with Krew-
cun's model, but which can be used for other types of elements present in this software. We plan to carry
out this implementation, test it on SOFA and present the results at the conference.

Références

[1] C. Krewcun, Simulation numérique de déploiement de stent coronaire en vue d'un traitement per-
sonnalisé - validation sur fantôme acquis en imagerie microdensitométrique. Thèse, Université Cler-
mont Auvergne, 2020.

[2] C. Krewcun, É. Péry, N. Combaret, P. Motreff, L.Sarry, Simulation of Balloon-Expandable Coronary Stent Apposition with Plastic Beam Elements, in : D. Shen, T. Liu, T. M. Peters, L. H. Staib, C. Essert, S. Zhou, P-T. Yap, A. Khan (ed.), Medical Image Computing and Computer Assisted Intervention (MICCAI 2019), Lecture Notes in Computer Science, 11768 (2019) 191-199.

[3] S. Yoon, F. Barlat, Non-iterative stress projection method for anisotropic hardening, Mechanics of Materials, 183 (2023) 104683.

[4] World Heart Report 2023 : Confronting the World's Number One Killer, Geneva, Switzerland, World Heart Federation, 2023

[5] J.S. Przemieniecki, Theory of Matrix Structural Analysis, McGraw-Hill, New York, 1968.

[6] R. Verwaerde, Stratégie d'élaboration, d'identification et optimisation d'un modèle de connecteur pour la simulation d'assemblages boulonnés. Thèse, Université Paris-Saclay, 2020.

[7] C. Duriez, S. Cotin, J. Lenoir, P. Neumann, New approaches to catheter navigation for interventional radiology simulation. Computer Aided Surgery, 11(2006), 300-308.