1 Introduction

Bioresorbable polymeric stents (BPS) have been developed as temporary vascular scaffolds, designed to restore patency to obstructed vessels, and provide temporary support before disappearing over time, addressing remaining limitations (i.e., late stent restenosis) [1] and presenting several advantages over permanent devices [24]. BPS have often been manufactured using poly-L-lactide (PLLA), a semi-crystalline polymer which has emerged as a favourable candidate in their design [5, 6]. PLLA degrades in-vivo via chemical hydrolysis, where scission of the polymer backbone chains occurs at vulnerable ester linkages [7] and the active metabolism of degradation products eliminates the material from the body [8]. The total degradation time of PLLA in-vivo is approximately 2 years [9], believed to be a favourable time for superior vessel remodelling to occur [10]. Carboxyl chain ends formed through chain scission lead to an increase in the acidic environment within the polymer matrix, accelerating the rate of hydrolysis in an autocatalytic process and driving heterogeneous degradation [11].

Polyesters such as PLLA also experience changes in crystallinity during hydrolytic degradation [1215], when scissions in the amorphous (tangled) regions provide interwoven polymer chains with extra mobility, allowing them to realign and become packed closely together [16]. As recently reviewed by the authors [17], understanding the dominant physio-chemical processes involved in BPS material degradation can be aided by theoretical frameworks, for example kinetic reaction diffusion equations such as [1821]. For PLLA, microstructure effects due to chain scission, crystallisation and heterogeneous degradation have been experimentally observed to play a significant role in the changing mechanical properties of the material [13, 15, 22]. Currently for the BPS, physio-chemical modelling frameworks do not consider the impacts of such microstructure changes on the mechanical integrity and support capabilities of the device.

Finite element analysis (FEA) presents an attractive means to investigate the mechanical performance of BPS, and a number of such studies have examined the mechanical properties of PLLA for the BPS application [2328]. Geometry and material stiffness were found to significantly impact the radial strength, recoil characteristics and short-term mechanical performance of BPS [23], whereas a similar influence of recoil behaviour on design was noted [28]. Other works have investigated the long-term degradation dependant mechanical performance of PLLA based BPS, including [29] which examined the deformation and degradation of several stent designs using a damage based model established within a thermodynamic framework [30] (expanded on in [31]). A similar model was developed (using long-term experimental radial stiffness and strength data) to examine the mechanical interactions between a degrading BPS and a remodelling vessel [32]. Such studies do not explicitly consider molecular weight, crystallinity, or degradation product diffusion; instead, they typically represent degradation via a scalar parameter that captures the phenomenological response; although this approach provides insight into the changes in mechanical integrity of the BPS with time, it does not elucidate the effect of the degradation reactions on the mechanical behaviour.

Previously, the implementation of a physio-chemical degradation model for amorphous polymers [19] into the FEA software package Abaqus/Standard by the authors [33] allowed for predictions in the degradation rates of various BPS to be made. In the present study, to build on this, a degradation theory applicable to semi-crystalline polymers [34] is incorporated into the user material thermal subroutine (UMATHT) in Abaqus/Standard. The outputs of the degradation models are linked to mechanical behaviour via three different damage models which couple the changes in molecular weight and crystallinity with a hyperelastic constitutive model for PLLA mechanical behaviour. The effects of strut thickness and applied degradation product diffusion boundary conditions on the molecular weight and crystallinity of PLLA BPS under simulated degradation are examined, and the impact of material heterogeneities and mechanical load boundary condition on the scaffolding performance and elastic properties of the degrading stent are investigated.

2 Modelling Framework

The theoretical framework originally developed by Gleadall et al. [34] for semi-crystalline polymers is presented here for clarity and the simpler theory for amorphous polymers [19] (which was previously implemented in Abaqus by the authors [33]) is summarised in Appendix A. The constitutive formulation used for the PLLA material is presented followed by the coupling between the degradation and mechanical components of the framework.

2.1 Degradation and Crystallinity Model

A theory, based on first-order reaction kinetics, for the degradation of semi-crystalline biodegradable polymers ([20, 35]) is summarised here. The rate of polymer chain-scission in semi-crystalline polymers is described as:

$$ \frac{\partial R_{s}}{\partial t} = k_{1} C_{e} + k_{2} C_{e} \left ( \frac{C_{ol}}{1- X_{c}} \right )^{0.5}. $$
(1)

The number of chain-scissions per unit volume, \(\boldsymbol{R}_{\boldsymbol{s}}\) (mol/m3), is the total number of chain scissions that have occurred up to the current time \(\boldsymbol{t}\) and is related to the mole concentration of ester bonds in the amorphous polymer phase, \(C_{e}\) (mol/m3). The other quantities are the volume fraction degree of crystallinity, \(X_{c}\), and the number of ester bonds in oligomers per unit volume, \(C_{ol}\) (mol/m3). Hydrolysis and auto-catalysed degradation are governed by reaction rate constants, \(k_{1}\) (/day) and \(k_{2}\) (\(\sqrt{\mathrm{m}^{3}/ \mathrm{mol}} /\)day), respectively. Equation (2) shows how the time rate behaviour of oligomer degradation products is dependent on the oligomer production rate due to chain scissions, \(\partial R_{ol} / \partial t\), and on the diffusion behaviour:

$$ \frac{\partial C_{ol}}{\partial t} = \frac{\partial R_{ol}}{\partial t} + \left ( \nabla \cdot \left ( D\nabla C_{ol} \right ) \right ), $$
(2)

where \(D\) (m2/day) is the effective diffusion coefficient of oligomers given by:

$$ D= D_{a} \left \{ D_{\mathit{matrix}} + \left ( 1.3 p^{2} -0.3 p^{3} \right ) \left ( D_{\mathit{pore}} - D_{\mathit{matrix}} \right ) \right \} . $$
(3)

The porosity, \(p\); temperature dependant material constant, \(D_{a}\); diffusion in liquid-filled pores, \(D_{\mathit{pore}}\); and diffusion within the polymer matrix, \(D_{\mathit{matrix}}\), are further defined in [20, 35]. Diffusion of water is assumed to occur faster than hydrolysis and, therefore, to not substantially influence the kinetics over the long timescale of the stent degradation. The transport of water is purely by diffusion; therefore, advection does not contribute to transport of oligomers. The number of oligomers per unit volume, \(R_{ol}\) (mol/m3), is related to \(R_{s}\) as follows:

$$ \frac{R_{ol}}{C_{e 0}} = \alpha _{2} \left ( \frac{R_{s}}{C_{e 0}} \right )^{\beta }, $$
(4)

in which \(\alpha _{2}\) and \(\beta \) are empirical parameters and \(C_{e0}\) (mol/m3) is the initial concentration of ester bonds. This initial framework for the chain-scission induced crystallinity of biodegradable polymers [20, 35] was developed based on an extension of the classical Avrami theory [3638], and includes differential and integrational equations to represent the time-dependent nucleation and growth of crystallites. In the reduced form of this theory [34], the number of material parameters is greatly decreased by elimination of the integration equations. The assumption that nucleated crystallites grow instantly to a finite size is introduced (i.e., the growth of crystallites is therefore assumed to occur considerably faster than the chain-scission reactions) along with a volume fraction extended degree of crystallinity, \(X_{\mathit{ext}}\), which represents idealised crystal growth without impingement [34]. \(X_{\mathit{ext}}\) is determined by a probability, \(P_{x} \), that each chain scission will nucleate a crystal which grows to a finite volume, \(V_{c}\) (m3), related to Avogadro’s constant, \(\eta _{A}\) (mol−1), through:

$$ X_{\mathit{ext}} = P_{x} V_{c} \eta _{A} R_{s}. $$
(5)

A maximum degree of crystallinity, \(X_{\max }\), is introduced, so that the actual and the proposed extended degree of crystallinity are related by:

$$ \frac{\partial X_{c}}{\partial X_{\mathit{ext}}} = X_{\max } - X_{c}. $$
(6)

Based on this reduced theory, the decrease in \(C_{e}\) in the amorphous phase of the polymer originates from both the hydrolysis of the polymer chains (oligomer production rate due to chain scissions, \(\partial R_{ol} / \partial t\)) and the induced crystallisation of the mobile polymer chains (where \(\omega \) is the inverse molar volume of the crystalline phase (mol/m3)):

$$ \frac{\partial C_{e}}{\partial t} =- \frac{\partial R_{ol}}{\partial t} -\omega \frac{\partial X_{c}}{\partial t}. $$
(7)

The number-averaged molecular weight, \(M_{n}\) (g/mol), is determined from:

$$ M_{n} = \frac{( C_{e} + \omega X_{c} ) M_{0}}{N_{\mathit{chains} 0} + ( R_{s} - \frac{C_{ol}}{m} )}, $$
(8)

where the contribution of the initial molar concentration of the polymer chains, \(N_{\mathit{chains}0}\) (mol/m3), the molar mass of a repeating unit of the polymer, \(M_{0}\) (g/mol), and the average number of repeating units of the oligomers, \(m\), are defined. This theory has previously been shown to capture observed experimental crystallinity and molecular weight trends [34].

2.2 Mechanical Model

A hyperelastic model is used to represent the mechanical behaviour of the PLLA BPS examined in this study. The Knowles hyperelastic model has previously been shown to capture the mechanical response of PLLA under isochoric uniaxial extension [39]. As a thermoplastic polymer, PLLA experiences non-linear deformation behaviour at moderate strains [23, 40, 41]. The strain energy density potential, \(W_{K}\), of a decoupled form of the Knowles model used here takes the form:

$$ W_{K} = \frac{\mu }{2b} \{ \left [ 1+ \frac{b}{n_{h}} \left ( I_{3} -3 \right ) \right ]^{n_{h}} -1\} + \frac{1}{D_{1}} ( J-1)^{2}. $$
(9)

Here, \(\mu \) (MPa) is the material shear modulus, \(D_{1}\) is related to the material bulk modulus, \(n_{h}\) is a hardening/softening parameter, \(b\) is an empirical constant (\(b\) and \(n_{h} >0\)) [42], \(I_{3}\) is the third principal invariant of the right Cauchy-Green strain tensor (\(\boldsymbol{C}= \boldsymbol{F}^{\boldsymbol{T}} \boldsymbol{F}\)) and \(J\) is the determinant of the deformation gradient, \(\boldsymbol{F}\).

2.3 Damage – Model Coupling

Changes in mechanical material behaviour are coupled to physio-chemical degradation via a scalar degradation (or damage) parameter, \(d\), which is introduced into the model to represent the extent of physio-chemical degradation. Three different couplings are described below, the first of which uses the amorphous degradation theory and the other two use the semi-crystalline theory. The coupling between degradation and mechanical behaviour is captured by evolution of \(\mu \) (\(\mu = \mu _{o} (1-d)\)), where \(d\) is described each model below, with \(d=0\) representative of an undegraded, undamaged material and \(d=1\) for the fully degraded case.

Model 1

Model 1 uses the outputs of the amorphous degradation model (Appendix A), in which the molecular weight is assumed to change in proportion to the changes in the concentration of ester bond in the polymer, \(C_{e}\) (relative to the initial concentration, \(C_{e0}\)). The assumed relationship between \(C_{e}\) and \(M_{w}\) is further discussed in [33]. Therefore, the damage for Model 1 is:

$$ d_{M_{w}} =1- \overline{M}_{w} =1- \frac{C_{e}}{C_{e0}}. $$
(10)

Model 2

For Model 2, the damage parameter is dependent on the normalised number average molecular weight, \(\overline{M}_{n}\):

$$ d_{M_{n}} =1- \overline{M}_{n} = 1- \frac{M_{n}}{M_{n0}}, $$
(11)

where \(M_{n0}\) is the initial value of \(M_{n}\). In this model, \(M_{n}\) is directly computed by the using the degradation/crystallinity theory as described above (equation (8)) for semi-crystalline polymers.

Model 3

Finally, for Model 3, the same degradation and crystallinity equations are used as in Model 2; however, damage is now dependent on \(M_{n}\) and the degree of crystallinity, \(X_{c}\), relative to their original values:

$$ d_{X_{c}} = \left ( 1- \left ( \frac{M_{n}}{M_{n0}} \right ) \left ( \frac{X_{C}}{X_{C0}} \right ) \right ). $$
(12)

2.4 Vessel Behaviour

A three-layered artery, with distinct properties for each layer, is used to apply loading to the BPS, as shown in Fig. 2(b). The thickness of the adventitia, media, and intima layers are 0.34, 0.32, and 0.24 mm respectively [43]. Each layer of the artery is modelled using a sixth order, reduced form of the polynomial model, as described in [44], which captures isotropic deformation of arterial tissue in the simulations. Other anisotropic models of arterial tissue (e.g., [45, 46]) can be used with this framework as no coupling is required to the arterial tissue model; a simpler representation of the artery is used here for efficiency. The specific strain energy density potential applied here, \(W_{P}\), is given by:

$$ W_{P} = \sum _{i=1}^{N} C_{i0} \left ( \overline{I}_{1} -3 \right )^{i} + \sum _{i=1}^{N} \frac{1}{D_{i}} \left ( J-1 \right )^{2i}, $$
(13)

where \(N\) is the number of polynomial terms, \(C_{i0}\) and \(D_{i}\) are experimentally measured material parameters, and \(\overline{I}_{1}\) is the first invariant of \(\overline{\boldsymbol{C}}\), the modified form of \(\boldsymbol{C}\) (\(\overline{\boldsymbol{C}} = J^{-2/3} \boldsymbol{C}\)). Material parameters (Table 1) for each layer are as suggested by [44], based on the experimental testing of arterial tissue reported in [43].

Table 1 Material coefficients for the isotropic strain energy density function used to represent the behaviour of each arterial layer

2.5 Finite Element Model and Simulation Procedure

The first study focuses on the impact of reaction-diffusion boundary conditions on stent degradation. The mechanical behaviour is not included in this phase. Simulated testing of BPS is carried out using Abaqus/Standard (Abaqus v6.17, Dassault Systems Simulia, RI, USA). Stents are modelled using representative unit geometries of a BPS (Fig. 1). BPS degradation and crystallisation are simulated for 52 weeks (approximately 1 year) using the equations in Sect. 2.1 and implemented through a user material subroutine for a thermal analysis (UMATHT), which captures the reaction-diffusion equations by analogy to the heat equation [33]. The resulting molecular weight, \(M_{n}\), and degree of crystallinity, \(X_{c}\), for three different stent thicknesses (in the radial direction) are examined (78 μm, 156 μm and 234 μm). The median stent strut thickness is chosen based on trials of polymeric degradable scaffolds [2] and the other values are half and double this value, respectively, and in line with minimum and maximum reported values [2, 6]. A finite element (FE) mesh of the stent (DC3D20) is created with the mesh density targeted as shown (Fig. 1(b,g)).

Fig. 1
figure 1

The geometry used to represent a 3.3 mm BPS (a) is shown in (b). Two separate regions are created across the BPS strut; an inner core and outer surface, as labelled in (c). The FE mesh and degradation product boundary conditions assigned to the stent are highlighted in (d)-(g), with the different boundary conditions defined as (d) stent free surface (SFS), (e) stent partially embedded (SPE), (f) stent inner surface free (SISF), and (g) stent fully embedded (SFE)

The stent is modelled in a deployed configuration (3.3 mm diameter) and the effects of assumed oligomer boundary condition on degradation and crystallisation are investigated for different boundary conditions (Fig. 1(d)-(g)). In the first of these (Fig. 1(d), stent free surface(SFS)), an oligomer free surface (\(C_{ol} =0\)) is assigned at the inner (lumen) sides of the stent with the remaining surface designated as an impermeable surface (i.e., zero diffusive flux applied, \(\frac{\partial C_{ol}}{\partial n} =0\)) representing a BPS immediately post-implantation into an artery, where any oligomers formed on the assigned surfaces are assumed to be swept away by the physiological actions of blood flow, and oligomers formed on impermeable surfaces are retained. In order to examine partially embedded stents during endothelialisation, the impermeable surface is extended across the stent as shown in Fig. 1, with (e) and (f) signifying stent partially embedded (SPE) and stent inner surface free (SISF) conditions respectively. In the stent fully embedded (SFE) case (Fig. 1(g)), the BPS is assumed to be fully embedded in neo-intimal tissue following a period of implantation in an artery; zero flux assumption assigned to all surfaces.

In the second study, the main focus is the impact of degradation and crystallisation on the long-term mechanical performance of a BPS. In this phase, both the degradation and mechanical components of the framework are used and are coupled using one of the three methods described above. Degradation of a 156 μm BPS is simulated for 2 years (104 weeks) under SFE and SFS conditions using the coupled degradation and mechanical model, again using Abaqus/Standard. Two methods (Fig. 2) are used to apply initial mechanical loading to the BPS based on either (i) a uniform pressure (which idealises physiological load), or (ii) arterial loading, using a three-layered artery (Fig. 2(b)). For each method, the BPS is modelled with the deployed state as the reference configuration and symmetric boundary conditions (Fig. 2(a)) are applied to account for the remainder of the device geometry. An FE mesh of quadratic coupled temperature displacement (C3D20HT) elements is created for the stent, while the artery is meshed using 3D quadratic stress elements (C3D20H) (Fig. 2(d)).

Fig. 2
figure 2

Symmetrical boundary conditions are applied to the BPS for the mechanical simulations (a). A three-layered artery (b) is created to apply loading to the BPS geometry. The steps involved in this method are shown in (c); the initial configuration of the BPS/artery (left) is shown alongside the application of a pressure at the inner surface of the artery (centre), before pressure is gradually reduced and loading is transferred to the BPS (right). The FE mesh for both parts is shown in (d)

In the first loading case, an inward-acting, uniform, external pressure is applied to the outer surface of the BPS (to represent vessel loading) following previous approaches [47]. In the second case (Fig. 2(c)), an approach similar to that detailed in [48, 49] is taken to apply arterial loading to the deployed BPS. To account for the biomechanical regulation of the artery that occurs during stent deployment [43], a 5% axial in-situ pre-stretch is applied to the artery (further described in Zahedmanesh and Lally [50]) and summarised here. In order to simplify the deployment simulation, the artery is expanded radially using a uniform pressure on its internal surface such that it is larger than the deployed stent. The vessel is then moved over the deployed BPS, as depicted in Fig. 2(c). Pressure on the vessel is gradually reduced to a final value of 13.3 kPa (representing a mean blood pressure of 100 mmHg), allowing the artery to relax onto the BPS (which is over-expanded by 10% [23]). Contact is initiated between the artery and BPS during the pressure release step. Contact behaviour in the tangential direction includes a coefficient of friction of 0.25 [51]. Contact behaviour in the normal direction assumes “hard” contact and a linear penalty is used to enforce the pressure-overclosure constraint in Abaqus. This general approach for stent loading has previously been applied to successfully examine the degradation and loading of magnesium stents [48, 52]. A more complete simulation of crimping, deployment, and recoil may be used; however, the focus of the present work is on the computational framework for the coupled degradation-mechanical material behaviour.

As an additional study, simulations are performed on devices sectioned into different regions with varying values of initial crystallinity (Fig. 1) to investigate the impact of heterogeneous material properties. An outer surface and inner core are defined for the deployed stent (Fig. 1(c)) and different initial crystallinity (\(X_{c0}\)) values are assigned to each region.

The FE implementation of the semi-crystalline degradation equations (Sect. 2.1) in Abaqus, as presented here, is verified for PLLA against previous, simpler, implementations of the degradation-crystallinity equations (which neglected diffusion) [34], against semi-analytical (MATLAB) solutions, and against experimental data [53, 54]. This verification uses model parameters for PLLA at three different degradation temperatures (reported in the experimental study [54]) as given by Gleadall et al. [34]. An Arrhenius relation [35] is used to capture the temperature dependent changes in properties. For the verification, Abaqus/Standard simulations are run on a single element cube, with an assumed oligomer flux of zero acting on all surfaces (i.e., no diffusion or spatial derivatives). The system of ordinary differential equations is solved in Abaqus using a direct Euler scheme for the model variables (initial conditions of \(C_{e} = C_{e0} \left ( 1- X_{c0} \right );\ X_{c} = X_{c0};\ R_{s}, C_{ol}, R_{ol}, X_{\mathit{ext}} =0\)). In parallel to this, semi analytical solutions to the equations are obtained in Matlab using a Runge-Kutta solution method (ODE45). We have confirmed that the Abaqus and Matlab solutions are in agreement, verifying the implementation and demonstrating temporal convergence of the Abaqus time incrementation (further details and data are shown in Appendix B).

3 Material Parameters

The parameters used above in the verification are previously published values [34] which used an assumption of no diffusion during the calibration. Here, we recalibrate to the experimental data while including diffusion to more correctly represent the experimental conditions. The updated values of \(D_{a}\), \(k_{1}\), and \(k_{2}\) are shown below for PLLA at 37 °C. For the recalibration, oligomer diffusion is assumed to occur at a much faster rate in the liquid filled pores than in the polymer matrix (ratio of \(\frac{D_{\mathit{pore}}}{D_{a}} =1000\) is assumed) [35]. All other parameters are as given by Gleadall et al. (key parameters are in Table 2) [34].

Table 2 Calibrated parameters used in the fitting of the reaction rate and diffusion constants for the degradation-crystallinity and diffusion theory implementation in Abaqus

The BPS are assigned mechanical properties representative of PLLA at 37 °C based on experimental observations of PLLA [41]; properties, including the hardening/softening parameter, \(n_{h}\), were previously calibrated to capture the experimentally observed behaviour. For all simulations \(\mu = 612.5\text{ MPa}\), \(D_{1} =0.001\text{ MPa}^{-1}\), \(b = 0.00249\) and \(n_{h} = 8.1 \times 10^{-5}\).

4 Results

4.1 Degradation, Crystallinity and Diffusion in BPS

Using the semi-crystalline degradation model, predictions of \(M_{n}\) and \(X_{C}\) in BPS after 52 weeks are shown in the contour plots in Fig. 3 for each of the stent boundary conditions investigated. Fully embedded stents (SFE) (Fig. 3(d)) exhibit significantly faster degradation (\(M_{n}\) of 13536 g/mol) than is predicted for the other SFS, SPS, and SISF conditions (Fig. 3(a), (b) and (c) respectively) which range from 79055 to 80525 g/mol. Predictions of \(X_{c}\) show a rapid increase in SFE conditions, where a maximum of 0.64 (64%) is predicted after 1 year (Fig. 3(d)), compared to the initial value of 0.48. Molecular weight trends are further examined in Fig. 4, where BPS predictions after 104 weeks are compared using a volume averaged form, \(\overline{M}_{n}^{V}\), for each case.

Fig. 3
figure 3

Contour plots show the predicted \(\boldsymbol{X}_{\boldsymbol{C}}\) and \(\boldsymbol{M}_{\boldsymbol{n}}\) (g/mol) in the stent after one year of simulated degradation under different degradation boundary conditions; (a) stent inner surface free (SISF), (b) stent partially embedded (SPE), (c) stent free surface (SFS) and (d) stent fully embedded (SFE)

Fig. 4
figure 4

Predictions of \(\overline{\boldsymbol{M}}_{\boldsymbol{n}}^{\boldsymbol{V}}\) (b), (c) and \(\boldsymbol{X}_{\boldsymbol{c}}\) (d), (b) are shown for the degradation models and various boundary conditions after two years of degradation. SFE predictions (full lines) are shown alongside SFS, SPE and SISF trends (dashed lines). Note that SFS, SPE, and SIFS predictions are similar to the extent that the lines are coincident as shown. In the case of stents with initial heterogeneity in \(\boldsymbol{X}_{\boldsymbol{c}}\), volume average quantities are shown in (b) and (d) and the values in the inner core (triangles) and outer surface (squares) regions are shown in (c) and (e), respectively

Figure 4 shows the evolution with time of \(\overline{M}_{n}^{V}\) and \(X_{c}\) for the semi-crystalline model; the change in \(\overline{M}_{n}^{V}\) for the amorphous degradation model are included for comparison. A more gradual decrease in \(\overline{M}_{n}^{V}\) is predicted for amorphous polymers in SFS, SPE and SISF cases compared to predictions for semi-crystalline polymers under similar degradation boundary conditions. In contrast for the SFE condition, the changes in \(\overline{M}_{n}^{V}\) for both theories are similar (Fig. 4(b)), but the drop in molecular weight is substantially faster. At the one-year point, the semi-crystalline model shows a reduction of ∼40% for SFS, SPE, and SISF, whereas the SFE condition has in effect completely degraded (>90% reduction).

Although Fig. 3 shows negligible gradients across the stent strut in both crystallinity and molecular weight, the effect of diffusion is clear when the SFS, SPE, and SIFS cases (which all have at least one surface on which oligomers are removed) are compared to the SFE case (which has no oligomer flux on any surface). The degradation rate in the SFE case is much faster due to the accumulation of acidic oligomer ends which lead to autocatalysis and faster degradation. Changes in crystallinity are similarly different as the different evolution in molecular weight also affects the evolution of crystallinity with the SFE case increasing in crystallinity much faster. The impacts of these changes on the elastic modulus are discussed below, but it is also worth noting that a large increase in crystallinity is likely to cause a more brittle response. Consistent with these observations that no significant spatial gradient in molecular weight or crystallinity is computed, stent strut size does not have an effect on the degradation rates (data not shown).

In the results discussed above, we have focused on the homogeneous case, i.e., where the initial conditions for the polymer, including crystallinity, are homogenous throughout the stent strut albeit with small heterogeneity occurring during degradation (Fig. 3). Motivated by the observations of Wang et al. [55], we investigate the effect of initial heterogeneity in crystallinity on the degradation response. The crystallinity of the outer core is increased to 0.61, which leads to an initial volume averaged crystallinity of 0.54 as seen at \(t=0\) in Fig. 4 (d) and (e). Generally, the degradation proceeds similar to the homogenous initial condition. The molecular weight changes at the same rate for each boundary condition (SFE, SFS, SPE, SISF) and although the initial values of volume averaged crystallinity are different, the SFE case plateaus at the same maximum possible value and the SFS, SPE and SIFS cases increase in crystallinity at a similar rate (Fig. 4 (d) and (e)).

4.2 Mechanical Performance of BPS

Following the above simulations for a subsection of a stent strut, simulations of a stent repeating unit subject to mechanical loads during degradation are performed. As no substantial difference was seen between the SFS, SPE and SIFS boundary conditions, only SFS and SFE conditions are now considered. The changes in stent diameter during degradation are used as a measure of lumen loss, \(\phi _{L}\):

$$ \mathit{Loss}\ \mathit{of}\ \phi _{L} = \frac{D_{I} - D_{t}}{D_{I}}, $$
(14)

where \(D_{I}\) is the initial internal diameter of the deployed stent geometry at time zero and \(D_{t}\) is the current internal diameter at each time point. As described above two loading regimes are considered: (i) an idealised uniform pressure load and (ii) a physiological load from a pseudo-deployment into an artery. The three damage models (and associated degradation models) are compared (Fig. 5).

Fig. 5
figure 5

Normalised radial force (N/mm) predictions for BPS under arterial loading (b) are compared over 2 years for the various models (summarised in (a)) under SFS (dashed lines) and SFE (full lines) boundary conditions. Predicted loss of \(\boldsymbol{\phi }_{\boldsymbol{L}}\) for BPS under an applied uniform pressure (c) and arterial load (d) are examined alongside predictions of material damage ((e) and (f) for uniform pressure and arterial loading cases respectively)

In the case of arterial loading, the stent after deployment exerts a radial outward force per unit length of the stent of ∼1.4 N/mm. This radial force is calculated by integration of the contact pressure over the contact area. In the case of the SFS conditions, the radial force decreases linearly over the two years considered for all versions of the degradation theories and damage models; however, the drop in force is not substantial. In contrast, the radial force decreases by half for all models after one year and the stent provides no scaffolding support at two years. Examining the arterial loading case further, Fig. 5(d) shows the loss of lumen patency is substantially greater at one year for the SFE condition than for SFS. This in turn is explained by the almost complete loss of mechanical stiffness due to the damage approaching unity for each model. In the previous degradation analysis, molecular weight is seen to decrease faster for the SFE case which consequently leads to much faster damage accumulation. In Models 2 and 3, the degradation proceeds at the same rate (both use the semicrystalline theory), however lumen loss for Model 3 is slower due the compensatory effect of increasing crystallinity (Fig. 4(d)) which, according to equation (14), would lead to an increase in modulus.

Initial heterogeneity in crystallinity does not have a substantial effect on the mechanical behaviour during degradation for the SFE case (Fig. 5(b,d,e)) as both the inner and outer regions increase in crystallinity at similar rates until the maximum level is reached (Fig. 4(d)). The increased crystallinity of the outer surface reduces the amount of lumen loss in the SFS case, however, the total amount of lumen loss in both the homogenous and heterogeneous cases is not large for SFS conditions and furthermore, SFS conditions are not representative of long degradation times, where endothelialisation would lead to the fully embedded boundary condition.

In the arterial loading case, as lumen loss proceeded, the load transferred to the stent is also reduced, therefore the uniform pressure loading regime is also considered (Fig. 5(c,e)). Here, the rapid accumulation of damage and loss of diameter is apparent for the SFE case; a more gradual lumen loss is evident for SFS, however, the lumen loss in this case if larger than observed for the arterial loading. Furthermore, the small differences in degradation and changes in mechanical behaviour are amplified for the SFS case. The semi-crystalline theory coupled with the Model 2 damage behaviour (i.e., damage only dependent on molecular weight) shows the fastest degradation and would therefore form the basis of the most conservative estimate.

The stress distributions in the stent are shown in Fig. 6 after 20 weeks of simulated degradation under uniform pressure loading for each of the three damage models and also for Model 3 with the heterogenous initial crystallinity. The initial stress distribution can be seen on the inside of the curved regions and is approximately 70 MPa. During degradation under the uniform pressure load, the stress distribution remains similar, although the maximum stress has a larger reduction in SFE (values of 57.4 MPa, 59.6 MPa and 61.6 MPa for Model 1, 2 and 3 respectively) than in SFS (69.7 MPa, 65.6 MPa and 66.2 MPa). Heterogeneous BPS (in this case Model 3 applied only) display a marginally higher initial stress than the uniform stents, with the maximum \(\boldsymbol{\sigma }_{\boldsymbol{VM}}\) in heterogenous BPS predicted to reduce by about 4% and 1.4% for SFE and SFS cases respectively (Fig. 6(d)) at 20 weeks.

Fig. 6
figure 6

Contour plots of \(\boldsymbol{\sigma }_{\boldsymbol{VM}}\) (MPa) in BPS with uniform (homogeneous) ((a)-(c)) and heterogeneous (d) material properties following initial application of uniform pressure and 20 weeks of simulated degradation under stent fully embedded (SFE) and stent free surface (SFS) boundary conditions. The stress predictions in each model (Model 1 (a); Model 2 (b); Model 3 (c), (d)) are also compared

The maximum principal stress, \(\boldsymbol{\sigma }_{\boldsymbol{MP}}\), in the artery with heterogeneous BPS degradation under SFS and SFE boundary conditions is examined in Fig. 7. Maximum values of stress are generated in the artery in both cases after contact is initiated with the BPS (day 0 in degradation simulations). \(\boldsymbol{\sigma }_{\boldsymbol{MP}}\) in the artery shows a gradual decrease for the SFS case in comparison to the SFE simulation, where a 26% decrease is predicted after 189 days, or 27 weeks (compared to 8% for the SFS case at the same time point).

Fig. 7
figure 7

The \(\boldsymbol{\sigma }_{\boldsymbol{MP}}\) in the artery following heterogeneous BPS degradation under (a) SFS and (b) SFE conditions is investigated at a number of time points throughout the degradation

5 Discussion

In this study, a framework for simulating polymer degradation and crystallisation and the impact on the scaffolding performance of bioresorbable stents is developed. The approach combines equations of a degradation model for semi-crystalline polymers [34] with a mechanical model for PLLA; three different coupling frameworks are introduced here each of which uses a scalar damage term to alter the mechanical response as a function of the polymer properties (molecular weight and crystallinity). Predictions of BPS degradation and crystallisation are examined for various diffusion boundary condition assumptions, before the predicted damage and scaffolding ability (\(\phi _{L}\) lumen loss) of several BPS configurations is evaluated.

The BPS simulations suggest a considerable dependency of predictions on the assumed degradation product boundary conditions. As is evident for \(M_{n}\) and \(X_{c}\) trends, significant increases in the predicted rates of chain scission (hydrolysis leading to \(M_{n}\) decrease) and the formation of crystalline regions (increasing \(X_{c}\)) are apparent for BPS where degradation products are considered to be confined within neo-intimal tissue (i.e., SFE), compared to when oligomer flux out of the BPS occurs. The inclusion of an oligomer free surface condition at any surface of the stent (representing the eradicating actions of in-vivo blood flow) lessens the autocatalytic effects observed for the PLLA material. The autocatalytic nature of polyester degradation has been widely investigated in experimental studies (for example, in Siparsky et al. and Yuan et al. [18, 22]), while few computational studies of BPS have considered this phenomenon [21, 56]. Degradation driven changes in \(X_{c}\) are believed to contribute to the increased stiffness and higher tensile properties often observed for degrading PLLA specimens [13, 54, 57, 58] and here were shown to alter the rate of lumen loss for BPS under uniform pressure loading. Additionally, device processing and manufacturing methods can lead to heterogenous material properties [55], which as examined here for \(X_{c}\), can impact the degradation rates. In the current work, the stress state does not influence the reaction kinetics that drive degradation; inclusion of stress-dependant coefficients similar to Gatica et al. [59] would allow for model enhancement. Experimental studies which have examined the impact of in-vivo stress and loading on the degradation behaviour of BPS have shown the importance of static [60] and dynamic loading [41, 57] on the lifetime prediction of BPS. The framework presented here may be further enhanced through knowledge of BPS processing, manufacturing, and deployment techniques (that can alter molecular weight, crystallinity, material heterogeneity, and the strain history) which have recently been investigated as a key influence on BPS performance [55].

The importance of stent endothelialisation in predicting device degradation rates is highlighted through the simulations in this work. The considerable impact of boundary condition on the predicted rates of damage and \(\phi _{L}\) loss (Fig. 5) draws attention to the need for accurate knowledge of tissue response in examining stent performance. BPS degradation in-vivo occurs alongside several key tissue responses (e.g., inflammation, endothelial and smooth muscle cell proliferation, neo-intima formation and endothelialisation during the early stages of implantation, followed by neo-intimal remodelling (NIR) [61]). The majority of this active response of arterial tissue post stent implantation is generally believed to occur within the first 3-6 months [62, 63], suggesting that the boundary conditions in reality transform and fluctuate throughout the lifetime of the implanted device. Alongside this, incomplete vessel healing leading to uncovered struts [64] would further impact the diffusion of degradation products. In-vivo animal and clinical studies designed to evaluate tissue response alongside BPS performance would aid enhancement of BPS degradation models. In computational studies of bioresorbable magnesium stents (BRMS), NIR significantly impacts on the mechanics of the stent, along with the predicted in-vivo degradation rates [48, 52, 65].

The framework presented here allows for integration with other computational techniques that were not included in the model, such as the NIR models discussed above, or material models which also incorporate the effect of plasticity in the stent material. The role of plasticity during degradation is not well described in the literature. Experimental studies have suggested that plasticity during crimping and deployment can influence the degradation rates [55, 66] and some computational models have proposed linking the degradation rate to plastic strain [60, 67]; however, these models do not account for other effects such as hydrolysis or autocatalysis. While some studies have included plasticity in computational investigation of BPS performance [27, 28], these approaches have focused on fitting of empirical degradation data to elastic-plastic formulations; the authors are not aware of a computational modelling framework that considers changes to the polymer (i.e., molecular weight and crystallinity) and plastic mechanical properties (changes in yield criteria or hardening laws). For other non-polymeric bioresorbables such as magnesium, plasticity is more readily included as the degradation mechanism is primarily surface based and changes to mechanical properties are not the dominant mechanism [48]. In the present study, the focus is on the post recoil behaviour; future developments of this framework could include the full plastic loading history and changes to the kinetic formulation to capture changes in plastic behaviour during degradation.

In the three damage models described here, the damage variable \(d\) is written in terms of the physical properties (molecular weight and crystallinity). These forms are chosen to allow for comparison to literature reports of these properties during degradation, for example [13, 15]. In model 1, the term \(C_{e}\) is used as a proxy for molecular weight. This approach allows the framework to be readily adapted to include other degradation process, for example enzymatic degradation in poly(4-hydroxybutyrate) [68]. The dependence of modulus on damage can also be modified for other material constitutive laws.

The simulations presented here have focused on a small representative section of a vessel and used a simple material model for the vessel mechanics for computational efficiency and to limit the number of variables considered. However, the framework is readily adapted to use more recent models for arterial tissue which include anisotropy [45, 69], patient specific vessel geometries [70], vessel material failure [71] or muscle contractility in the vessel wall [72]. The main conclusions presented here in terms of the importance of the transport boundary conditions and the use of this framework as a means to assess and compare device performance are not influenced by these assumptions; although the more detailed arterial models would facilitate the investigation of tissue responses discussed above.

Radial force testing is often applied to evaluate stent long-term performance; values of radial strength and radial stiffness [26, 73] allow for insight into a stent’s ability to provide scaffolding support to the vessel wall. Such properties are of particular importance in the polymeric stents, due to the reduced mechanical properties of polymers compared to metallic materials [5]. The impact of BPS degradation on radial stiffness is evaluated here; the value of initial radial force (1.41 N/mm) is similar to that experimentally reported by Wang et al. [26] and predicted by Qui et al. [32], with predictions for decreasing radial properties for the two examined mechanical load boundary conditions falling within the range of those observed in Qui et al. [32].

Clinical guidelines for use of BPS include recommendations for device deployment, and in general a pre-dilation, sizing, post-dilation (PSP) strategy is observed for implanting BPS in diseased arteries [74]. Consideration of such guidelines in the computational modelling of BPS deployment is important in order to simulate realistic stresses in both the device and the artery. BPS examined here are modelled in a deployed state, with 10% over expansion in the artery [23]. BPS deployment and implantation produce high stresses in the arterial wall, which can often trigger stent restenosis [61, 75] or stent thrombosis [76] due to vessel injury (which may occur as a result of over-expansion, incomplete BPS expansion or strut fracture during deployment). The present work allows for the investigation of vessel stresses and the predictions are in agreement with that reported in Qiu et al. [32], which investigated the mechanical interactions between degrading BPS in a combined artery and plaque model. Vessel stresses during vessel implantation are linked to arterial remodelling [61, 77] and inclusion of such processes in a computational framework (e.g., [65]), would allow for prediction of clinical outcomes. Endothelialisation and neointimal growth can be combined with the lumen loss predicted here to form an overall measure of vessel patency, which is often reported as a primary endpoint in clinical studies of BPS in-vivo performance [63].

In conclusion, the capability of a computational framework for simulating physio-chemical and damage-based degradation of BPS and predicting in-vivo scaffolding performance is demonstrated. Using this framework, the loss of \(\phi _{L}\) in stents degraded under various conditions is predicted, which provides insight into their long-term scaffolding ability. Although the test cases presented here make a number of simplifications for computational efficiency, the framework presented here can be used as presented to investigate performance of full BPS geometries in patient specific vessel models with more complex arterial tissue constitutive relations or include the effects of vessel health, lesions or plaques.