Propagation of dissection in a residually-stressed artery model

This paper studies dissection propagation subject to internal pressure in a residually-stressed two-layer arterial model. The artery is assumed to be infinitely long, and the resultant plane strain problem is solved using the extended finite element method. The arterial layers are modelled using the anisotropic hyperelastic Holzapfel–Gasser–Ogden model, and the tissue damage due to tear propagation is described using a linear cohesive traction–separation law. Residual stress in the arterial wall is determined by an opening angle \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\alpha $$\end{document}α in a stress-free configuration. An initial tear is introduced within the artery which is subject to internal pressure. Quasi-static solutions are computed to determine the critical value of the pressure, at which the dissection starts to propagate. Our model shows that the dissection tends to propagate radially outwards. Interestingly, the critical pressure is higher for both very short and very long tears. The simulations also reveal that the inner wall buckles for longer tears, which is supported by clinical CT scans. In all simulated cases, the critical pressure is found to increase with the opening angle. In other words, residual stress acts to protect the artery against tear propagation. The effect of residual stress is more prominent when a tear is of intermediate length (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\simeq $$\end{document}≃90\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^\circ $$\end{document}∘ arc length). There is an intricate balance between tear length, wall buckling, fibre orientation, and residual stress that determines the tear propagation.


Introduction
An arterial dissection is a tear within the wall of a large artery, such as the aorta. The dissection can lead to the creation of a false lumen through which blood flows, and propagation of the tear can quickly lead to death as a result of decreased blood supply to other organs, damage to the aortic valve, and rupture of the artery. The loading conditions on the arterial wall, the geometry of the artery and of the tear, and the material properties of the arterial wall determine whether the tear propagates. A prediction of how the critical condition for tear propagation depends on these factors could help to optimize diagnosis and treatment.
In the absence of loading, many biological soft tissues are not stress-free, but subject to residual stress. At physiological loading, the residual stress in arteries reduces variation in the stress distribution across the arterial wall and decreases the peak stress (Cardamone et al. 2009;Chuong and Fung 1986). Fung (1991) was the first to show that a radial cut along artery can release much of the residual stress. Hence, using an opening angle is a theoretical approach for recovering the stress-free configuration, and the value of the opening angle is often used to quantify the residual stress. For example, Holzapfel et al. (2000) employed this method to obtain the residual stress in a two-layer model of a rabbit carotid arterial wall.
Here we study the effects of residual stress on an arterial dissection or tear. If the artery is subject to hypertension (Golledge and Eagle 2008;Kodolitsch et al. 2000) or the arterial wall becomes weak, e.g. a defect arises in the inner surface of arterial wall, then blood at high pressure may force its way into the wall and propagate longitudinally, creating a dissection, as indicated in Fig. 1.
Tear propagation in the absence of residual stress has been studied previously. For example, Sommer et al. (2008) Tear Fig. 1 A cartoon of arterial dissection. A defect (a tear) that initially begins in the intimal layer of the artery extends into the media where it can grow to create a false lumen. Inflation of the false lumen by blood pressure can restrict the true lumen, restricting blood supply to major organs. In some cases, the tear can rupture the artery. The loading condition on the arterial wall, the geometry of the artery and of the tear, and the material properties of the arterial wall determine whether the tear propagates subjected samples of human aortic media to peeling tests to estimate dissection property. This was followed by a number of computational simulations (Ferrara and Pandolfi 2010;Gasser and Holzapfel 2006). Wang et al. (2014) and Tong et al. (2011) performed similar experiments on human coronary and human carotid arteries, respectively. In these experimental and computational studies, the tear propagation is stably driven by controlling the displacement when peeling apart strips. However, blood pressure drives the propagation of dissections in vivo (Braverman 2010;Rajagopal et al. 2007). There are fewer literatures on studying this pressure-driven propagation. In experiments on a porcine thoracic aorta subject to pressure, Carson and Roach (1990) measured the peak pressure to tear the media and the work per unit area of tissue required to propagate a tear and showed that these values are independent of the tear depth, while Tam et al. (1998) studied the effect of depth of the initial tear on the critical pressure for propagation and showed that the critical pressure decreases as the depth increases. Arterial dissection during balloon angioplasty of an atherosclerotic artery was modelled by Badel et al. (2014), in which the arterial wall is compressed by inflating a balloon controlled by a displacement boundary condition. Recently, we developed a computational scheme to compute the energy release rate, a variable for quantifying the risk of propagation, for pressure-driven dissection propagation using a nonlinear energy argument in a 2D model (Wang et al. 2015), and reported how the critical pressure for arterial dissection can be affected by fibre orientation, tear length, and surrounding tissues.
In this paper, we extend our previous work (Li 2013;Wang et al. 2015) and investigate what role the residual stress plays in the critical condition for propagation of the dissection in a two-layer (media and adventitia) arterial wall. To model the tear propagation, we use the eXtended Finite Ele-ment Method (XFEM) (Moes et al. 1999) implemented in ABAQUS (2014. For the material properties, we use the fibre-reinforced anisotropic hyperelastic Holzapfel-Gasser-Ogden (HGO) constitutive law (Holzapfel et al. 2000). We developed a computational protocol that involves a sequence of novel boundary conditions to close an opening angle, starting from the stress-free configuration, to obtain the unloaded configuration and so introduce the residual stress. Notably, the unloaded configuration is the same for different values of the opening angle, ensuring that the difference between simulations is solely due to the residual stress field. To achieve this, we solve the equilibrium equations analytically to obtain the individual stress-free configuration for different opening angles.
The paper is organized as follows. In Sect. 2, we detail the arterial model including the geometry, the constitutive model, and the cohesive law. In Sect. 3, we show how to calculate analytically the stress-free configurations for different opening angles with a specified unloaded configuration, how to incorporate residual stress using ABAQUS, and how to determine the critical pressure. The results are shown in Sect. 4, followed by discussion and conclusions in Sect. 5.

The model 2.1 Geometry
Our model is based on the widely used two-layer rabbit carotid artery (Holzapfel et al. 2000), but includes the effects of residual stress as shown in Fig. 2. An arterial dissection is modelled as a tear in the wall of thick-walled cylinder in an unloaded and residually-stressed configuration Ω r . The Fig. 2 Cross sections of the arterial wall in the stress-free Ω 0 and unloaded Ω r configurations, subject to the same pressure loading on the inner radial and the tear surfaces. Each layer is a hyperelastic fibrereinforced material, modelled by the HGO strain-energy function. β m and β a define the angles between the two families of fibres in the media and adventitia. The tear subtends an angle η residual stress field in Ω r is calculated from the value of the opening angle α in the zero-stress configuration Ω 0 . The configurations are listed in Table 1. The data for Ω 0 are obtained from (Holzapfel et al. 2000), and the data for Ω r are computed using the analytical approach to be discussed below.
We introduce an idealized dissection (Fig. 2) along an arc of constant radius in Ω r , which is connected to the lumen of the vessel by negligibly small tears, so that the dissection surfaces are subject to the same blood pressure as in the true lumen.

Constitutive law
The description for the two layers of the arterial wall is the same except for different values of parameters in each layer.
For the mechanical response of the arterial wall, we use the incompressible Holzapfel-Gasser-Ogden (HGO) constitutive law (Holzapfel et al. 2000), with strain-energy function where c, k 1 , k 2 are material parameters for each layer of the artery, as listed in Table 2 (Holzapfel et al. 2000), and I 1 , I 4 , and I 6 are the invariants of the right Cauchy-Green strain tensor C and the fibre-structure tensor M n , defined as Fig. 2. The two fibre families are aligned along the two directions, A 1 and A 2 , and only contribute to stress when stretched, i.e.
We assume that the artery is infinitely long (plane-strain problem), and the z-component in these fibre directions is ignored. Therefore, changing β only changes the contribution of fibres in the circumferential direction.

Cohesive law
The initialization and propagation of a tear is modelled using the XFEM in ABAQUS (2014). The existing tear is implemented as an initial condition. The discontinuity of the displacement field at the tear surface is modelled by adding an enrichment term H (r ) u/2 onto the otherwise continuous displacement field, where H (r ) = ±1 when r ≷ r t , r t is the radius of the tear surface in the configuration Ω r , and u is the displacement jump. We assume that the propagation of tear is governed by a linear cohesive traction-separation law (Ferrara and Pandolfi 2010) (Fig. 3), which is specified by the maximum traction T c just before damage, the separation energy G c , and the maximum displacement jump u c , as listed in Table 2. Only two are independent. After each incremental loading step, the maximum principal stress σ mp at the centroid of each element is compared to T c : if σ mp ≥ T c , then the displacement jump u is calculated. When Δu > u c , the tear propagates in the direction perpendicular to the maximum tensile principal stress.
The actual value of T c is material dependent and should be determined by experiments. In the absence of such data, we assume that T c /c = 2 (where c is the value of the media) in this study. This, together with the assumption of the planestrain problem, means that our computed results will be qualitative. However, different values of T c /c are used later to check that the trend we observe is the same.

Methodology
All simulations are based on a residually-stressed configuration Ω r . The residual stress and Ω r can be obtained both analytically and numerically.

Analytical approach
Ω r can be obtained by setting p i = 0 and λ z = 1, where p i is the pressure and λ z is the axial stretch. Let the stress- Table 2 The material parameters for the HGO strainenergy function (Holzapfel et al. 2000) and the cohesive law, for a rabbit carotid artery where R i , R o , and α denote the inner and outer radii, and opening angle, respectively. The corresponding unloaded configuration Ω r in polar coordinates (r, θ) is where r i and r o are the inner and outer radii of Ω r . Incompressibility requires that where where are the principal stretches. The Cauchy stress is given by where p is the Lagrangian multiplier associated with the incompressibility condition. Let σ rr and σ θθ be the radial and circumferential components of the Cauchy (residual) stress tensor in Ω r , which satisfy the momentum balance equation Integration of (9) leads to Using the traction-free boundary condition in (10), we have Substituting (5) and (8) into (12), we obtain a nonlinear integral equation, which is solved for r i using Newton iteration.
Substituting (1) into (8) gives where B = FF T is the left Cauchy-Green tensor, m n = FM n F T is the structure tensor in the residually-stressed unloaded configuration, and p is determined from (10) and (11) p(r ) = −σ rr (r ) + 2cλ 2 r (r ).
To ensure that the difference between various simulations is only due to the residual stress, we can also determine Ω 0 from (5), (7), and (12), given α and Ω r .

Numerical approach
We start from the stress-free configuration Ω 0 with a specified opening angle. The numerical approach is to close the opening angle numerically and obtain the unloaded configuration Ω r . This is achieved in several steps as shown in Fig. 4. During the closing process, the inner radial and tear surface are pressurized to avoid contact. Once the ring is closed, the artificial pressure is removed. The closed configuration is then inflated to simulate the deformation of the arterial wall and dissection propagation subject to the pressure.

Finite element implementation
The computations are performed using the finite element package ABAQUS (6.13). Four-node plane-strain hybrid elements are used to construct the mesh, as shown in Fig. 5. A grid independence test was used to select the optimal number (a) (b) (c) (d) Fig. 4 The closing process from Ω 0 to Ω r , for opening angle α = 160 • and tear length η = 90 • : a the stress-free configuration Ω 0 , b, c move one end and pressurize both the tear and the inner radial surface, d residually-stressed Ω r after the artery is closed and pressure removed The stress-free Ω 0 (grey) and the unloaded Ω r (coloured) configurations, where the colour indicates the magnitude of the maximum principal stress, which is along the circumferential direction, and the opening angle in Ω 0 is 160 • of elements (Table 3), and the computed stresses converge to the analytical solutions as shown in Fig. 6. The relative error in computed residual stress is where σ a rr and σ a θθ are the vectors of all the nodal values of the exact analytical expression for the two components of the residual stress across the wall, and · denotes standard L 2 -norm, i.e. x = i x 2 i 1/2 for a vector x. The intermediate mesh was then used in all the simulations for different opening angles.

The residual stress and critical pressure
The geometries (i.e. the radius R i , thicknesses of the media and adventitia, T m , T a ) of stress-free configurations Ω 0 associated with the same unloaded configuration Ω r (where we specify α, thicknesses of the media and adventitia, t m , t a , and inner radius r i which are equal to T m , T a , and R i when a = 0) are shown in Table 4. The residual stress components computed analytically are shown in Fig. 7 as a function of the opening angle α. The absolute value of σ rr is greatest at the mid-radius of the media and increases with α. The circumferential stress σ θθ is in compression at the inner radius of the media and is in tension at the outer radius of the media and adventitia, and |σ θθ | increases with α. The residual stress is smaller in the adventitia. The simulated configurations for α = 160 • and η = 90 • are shown in Fig. 8. The critical pressure at which the tear starts to propagate is identified in Fig. 8c. Since the value of p c changes with the material properties, we focus on the dimensionless critical pressure p c = p c /c (where c is the value of the media). For the tear of length η = 90 • , the change of p c with the opening angle is plotted in Fig. 9. Notice that p c increases with α in all the cases simulated, suggesting that existence of residual stress makes artery more resistant to the tear propagation.

Inner wall buckling and tear length
Our previous study on a 2D strip (Wang et al. 2015) showed that tear length plays an important role in the dissection; a   Table 4 Geometries of Ω 0 (thicknesses of the media and adventitia T m , T a , and inner radius R i ) corresponding to the specified Ω r and opening angle α longer tear is always more likely to propagate. We now show that this is no longer true for a circular geometry and when residual stress is present. Three additional groups of simulations were carried out with η = 30 • , 150 • , and 210 • , and some interesting results are found when tear length increases. The deformed configurations for α = 160 • , Fig. 10, show that buckling of the inner wall (the material section between the lumen and tear) occurs for η = 150 • and η = 210 • . The maximum principal stress distribution and the critical pressure are also shown in Fig. 10. Beyond the cohesive criterion, the tear propagates. All tear propagations tend to be radially outwards, as shown in Fig. 11. Comparison of the critical pressure for different tear lengths (Fig. 9) shows that the dimensionless critical pressure p c increases with α in all the cases simulated. Notably, the longest tear length studied (η = 210 • ) has a higher value of p c than that of η = 90 • , and 150 • . Hence, the relationship between propagation and tear length is not as simple as that of a 2D strip (Wang et al. 2015). This presumably is due to the buckling of the inner wall, which is more likely to occur for a longer tear. However, for the intermediate case, η = 90 • , inner wall buckling occurs only for α = 0 • , and not when α ≥ 40 • (Fig. 12). This suggests that there is a subtle interplay between dissection length, residual stress, and inner wall buckling.

The effect of fibre orientation
Here we vary the fibre orientation. We refer to the fibre orientation listed in Table 2 as the physiological or 'true'  Fig. 9 The dimensionless critical pressure p c versus the opening angle α, for four different tear lengths η, where letter 'B' indicates that the inner wall is buckled. For discussion of different tear lengths, see Sect. 4.2 case. In the following simulations, we assume that the media and the adventitia have the same fibre orientation, so that β m = β a = β. Simulations were performed with β = 0 • , 10 • , 15 • , 20 • , 30 • , 60 • and 90 • . In addition, a group of simulations were run without fibres (k 1 = 0), referred to as the 'free' case.
To ensure that the differences between the simulations are only due to the fibre orientations, we use the same unloaded configuration Ω r , with Ω 0 calculated for each fibre orientation using the analytical method. The normalized thicknesses of the media and adventitia, T m = T m /t m and T a = T a /t a , are plotted against the opening angle α in Fig. 13, for the different fibre orientations. The thickness of the media T m decreases with the opening angle for smaller values of β (e.g. when β = 0 • , or 10 • ), but the trend changes as β increases. When β ≥ 30 • , it increases monotonically. On the other hand, the thickness of the adventitia T a increases monotonically with the opening angle. For other fibre angles (β = 60 • , 90 • , free, true), the results are identical to that of β = 30 • . This is because fibres beyond this angle are no longer stretched, i.e. I 4 < 1.
The critical pressure also changes with different fibre orientation (Fig. 14). Notice that with the inflation, fibres with β = 30 • also start to bear load, but fibres at greater angles (β = 60 • and 90 • ) still do not take on any load, and hence, the critical pressures for these cases remain the same as in the 'free' case. As expected, the critical pressure when β = 0 • is the highest, since the residual stress is the greatest in this case.

Discussion and conclusions
We have used both analytical and computational approaches to study the effect of residual stress on the propagation of arterial dissection. Our simulations show that the shortest (30 • ) and the longest tear (210 • ) are the most stable (with higher critical pressures), while tears of lengths 150 • and 90 • are most unstable. However, residual stress increases the critical pressure in all cases; the most dramatic improvement is seen for the 90 • tear. This suggests that both the length of tear and the residual stress play an important role in determining the critical pressure for the tear propagation. In particular, we found interesting inner wall buckling associated with the longer tears. Similar buckling was observed in a computed tomography (CT) scans of a patient in the ascending aorta (Fig. 15). Clearly, there is an intricate balance between tear length, buckling of the inner wall, fibre orientation, and residual stress, all of which may affect the likelihood of tear propagation. One plausible explanation is that, although a longer tear propagates more readily, it is also more likely to cause the inner wall to buckle and so increase the value of p c . However, the exact mechanism remains to be explored. We now discuss the limitations of this study. Our model is a plane-strain problem and does not include the effect of axial stretch. In our previous work (Wang et al. 2015), we found that axial stretching of fibres resists the opening of the dissection and significantly decreases the energy release rate for tear propagation. Clinically, dissections may propagate axially and may re-enter the lumen. Our model also cannot predict the absolute value of the critical pressure due to the simplifications mentioned above and lack of data on cohesive parameters, e.g. the value of T c . Nevertheless, this model provides a qualitative description of the variation of the critical pressure with the residual stress. This is further illustrated in Fig. 16 for T c /c = 2, 6, 10 for the rabbit carotid artery, as well for the aged human thoracic aorta, where the shear modulus of the media is much greater, as shown in Table 5. In the human artery model, all the parameters for the HGO model   Fig. 13 Wall thickness T m and T a for the media and adventitia in Ω 0 normalized with respect to their values in Ω r , plotted against the opening angle for different fibre angles. For other fibre angles (β = 60 • , 90 • , free, true), the curves overlap that for β = 30 • (Table 5) are estimated by fitting the experimental data of the cyclic uniaxial tensile tests of 14 human thoracic aortas (60 ± 12 year, mean ± SD) (Fereidoonnezhad et al. 2016).
Our results show that changes of the critical pressure with residual stress are similar for human aorta and rabbit carotid artery.
The other limitation is that we use an isotropic cohesive traction-separation law. The tensile testing of a porcine thoracic aorta performed by MacLean et al. (1999) showed that the stiffness in the radial direction is significantly lower than in the circumferential and longitudinal directions. MacLean et al. (1999)    the behaviour of elastin layers and smooth muscle cells during the aortic dissection. Work is ongoing to develop a three-dimensional and anisotropic arterial dissection model that could make use of the histological data.