About prestretch in homogenized constrained mixture models simulating growth and remodeling in patient-specific aortic geometries

Evolution of mechanical and structural properties in the Ascending Thoracic Aorta (ATA) is the results of complex mechanobiological processes. In this work, we address some numerical challenges in order to elaborate computational models of these processes. For that, we extend the state of the art of homogenized constrained mixture (hCM) models. In these models, prestretches are assigned to the mixed constituents in order to ensure local mechanical equilibrium macroscopically, and to maintain a homeostatic level of tension in collagen fibers microscopically. Although the initial prestretches were assumed as homogeneous in idealized straight tubes, more elaborate prestretch distributions need to be considered for curved geometrical models such as patient-specific ATA. Therefore, we introduce prestretches having a three-dimensional gradient across the ATA geometry in the homeostatic reference state. We test different schemes with the objective to ensure stable growth and remodeling (G&R) simulations on patient-specific curved vessels. In these simulations, aneurysm progression is triggered by tissue changes in the constituents such as mass degradation of intramural elastin. The results show that the initial prestretches are not only critical for the stability of numerical simulations, but they also affect the G&R response. Eventually, we submit that initial conditions required for G&R simulations need to be identified regionally for ensuring realistic patient-specific predictions of aneurysm progression. Supplementary Information The online version contains supplementary material available at 10.1007/s10237-021-01544-3.

where this total deformation incorporate an elastic and inelastic deformation for each jth constituent in the mixture, equation 1. It is assumed that the material mixture is hyperelastic and depends only on the elastic deformation of each constituent, i.e. on the elastic Cauchy-Green tensor where C C C = F F F T F F F is the Cauchy-Green tensor related to the deformation gradient defined in Eq. 1. Thus, from the equations 3 and S.2 is possible to define stress in the jth constituent as, where the sum of the all individual stresses give the second Piola-Kirchhoff stress of the whole mixture (∑ n j φ j S S S j ), with the constituent mass fraction φ j = ρ j R ρ R and the density of the mixture (ρ R = ∑ n j ρ j R During G&R the tissue is continously adapting and consequently is changing its structure with depositon/removal of mass of the constituents, for instance, elastin degradation, SMCs apoptosis/proliferation or collagen production by fibroblasts/SMCs. The mass changes in the mixture can be written like, whereṁ is the depositon/removal of mass in the mixture, ρ R V is the density times the volume in reference configuration, ρv is the density times the volume in spatial configuration. The reference volume V does nont change, the spatial density is constant in time and homogeneous, and with the relation between the volumes v = JV (ρ R = ρJ), we can write the mass changes in reference and spatial configurations as it follows Rodriguez et al. (1994), then, the reference density (per unit reference volume) change over time ρ R = ρ R (t) Braeu et al. (2017). Furthermore, G&R occurs at slow time scales, then it can be assumed that it is quasi-static and the linear momentum balance equals zero, where v v v is the velocity of the system, divσ σ σ is the diveregence of the mixture Cauchy stress σ σ σ and the body force b b b, in the spatial configuration. At the surface of the body the conditions can be, given deformations (Dirichlet), external loads applied on the body surface (Neumann) or deformation-dependent forces (Robin, Moireau et al. (2012); Nama et al. (2020)). The Robin boundary conditions are introduced by the following expression: σ σ σ ·n n n = pn n n + ku u u, (S.8) where p denotes the pressure, n n n the normal to the surface, the elastic forces with stiffness k and displacement u u u. This forces are appropriate for the modelization of the lumen pressure in blood vessels due to the blood flow, and the elastic forces are usefull to give flexible and stable boundary displacements.

B Constitutive Models
At the equation 2 is assumed that the mixture strain energy W per unit reference volume is the sum of the constituents strain energies Ψ j per unit mass. Therefore, we consider our material to be vascular tissue composed by three constituents such as, elastin, SMCs and collagen. The elastin is assumed to be isotropic; the smooth muscle have an active, (•) act , and a passive ,(•) pas , behavior; and the collagen is composed by four fibre families aligned in circumferential, longitudinal and two diagonal directions, respectively. The total strain energy may be written such as for each component we assumed a strain energy to represent the corresponding hyperelastic behavior. We modeled the elastin with a Neo-Hookean hyperelastic model as in Mousavi and Avril (2017); Mousavi et al.
where µ l and κ l are stress-like material parameters (shear and bulk modulus), whileC C C l e and J l e are the elastic isochoric right Cauchy-Green tensor and the elastic jacobian of elastin, respectively. The isochoric right Cauchy-Green tensor is related to the right Cauchy-Green tensor byC C C l e = (J l e ) 2/3 C C C l e and the elastin elastic jacobian is J l e = det(F F F l e ).
We modeled the collagen fibre families by an anisotropic Fung-type exponential function such as where k c i 1 and k c i 2 are a stress-like and dimensionless material parameters, respectively, while λ c i e is the elastic stretch contribution of the collagen fibre obtained as where λ c i and λ c i r are the total and remodeling stretch (cf. 1) of the fibre, respectively. We also modeled the passive behavior of SMCs by an anisotropic Fung-type exponential function such as where k m 1 and k m 2 are stress-like and dimensionless material parameters, respectively, λ m e is the elastic contribution of SMCs calculated such as where λ m and λ m r are the total and remodeling stretch (cf. 1) of the fibre, respectively. While we modeled its active behavior according to Braeu et al Braeu et al. (2017), with σ actmax the maximal active Cauchy stress, λ act is the active stretch in the fibre direction, λ m 0 and λ m max are the zero and maximum active stretches and ρ R0 denotes the total mixture density in the homeostatic reference configuration.

C Growth and Remodeling
The idea of the hCM models is to pool all the sequential mass additions within one single change using temporal homogenization (Fig. 1). To do so, three assumptions are made: (i) the mechanical properties are changed by G&R, (ii) survival mass (mass turnover) functions are exponential and (iii) inelastic deformations, F F F j gr = F F F j g F F F j r , are in turn decomposed into growth-related, F F F j g , and remodeling-related (turnover-related), F F F j r , contributions. In this model a single local average inelastic deformation gradient F F F j gr is defined by constituent. The model can handle isotropic or anisotropic growth, the latter being more relevant for arteries and manifesting with thickening or thinning effects Matsumoto and Hayashi (1996). We assumed that G&R is a stress mediated process which tends to minimize deviations between the current stress and a reference stress metrics named homeostatic stress (σ σ σ h ). Therefore, the rate of mass degradation and deposition at time t for the jth constituent is expressed aṡ where ρ j R is the mass density (per unit reference volume) of the jth constituent at time t, k j σ denotes a mass-gain parameter, σ j is the spatial stress along the fibre (σ j = (a a a Then, the growth is measured from the changes of the reference mass density (ρ R ) respect to the initial (t 0 = 0) reference mass density (ρ R0 ) through noting that the mixture volume changes are measured by the determinant of the deformation gradient J = det(F F F) and assuming the elastic and remodeling processes are isochoric (det(F F F j e ) = det(F F F j r ) = 1), so, the mixture volume changes remains equal to the growth J = det(F F F g ). If the growth is assumed to be anisotropic and along the thickness direction (a a a ⊥ 0 ), can be expressed in tensorial form as where I I I is the identity second order tensor. Therefore, due to the continuous mass deposition and removal, the traction-free configuration change during G&R (Fig. 1), even if there is a balance between mass deposition and removal (ρ j R = 0), this occur with a prestress which is different from the current stress at which mass is removed or deposited. Altogether leads to changes of tissue microstructure referred as remodeling. Therefore, assuming that remodeling occurs at a constant volume and along a fibre in the direction a a a j 0 , the evolution of the remodeling of the jth constituent at time t is expressed such as where subscript "pre" indicates prestress,λ j r denotes the remodeling velocity and T j is the average turnover time during which old mass increment is degraded and replaced by a new mass increment. According to proposition 1 from Cyron and Humphrey , the prestress σ i pre is equal to the homeostatic stress σ i h . The remodeling of the fibre can be represented in tensorial form as To calculate the G&R deformation gradient over time we solved the system composed of Eqs. S.16, S.17 and S.19 by performing temporal integration, this is applied in a precedure to carry-out G&R simulation within a FE code, figure 1.

D Material properties
The mechanical properties of the patient-specific model are fitted from data available in the literature Mousavi et al. (2018) (bulge test). For this issue we assume that the material in the bulge have bi-tangential deformation (F F F = diag[λ θ 1 λ θ 2 λ r ]) and it is incompressible (detF F F = 1). The computation of the experimental measures as stretch and stress are made in base of the formulation presented in Rossi et al. (2020), from where we distribute the tangential stress σ b within the specific-layer stresses σ M and σ A , media and adventitia, respectively, with α as the media thickness ratio, Fig. 2. The definition of the constituent densities is made in base of the histological observation showed in Humphrey and Holzapfel (2012), where the thoracic aorta is composed by 35% of SMCs, 35% of elastin and 30% of collagen, according to our approach of three constituent in the arterial wall. After we place the constituents ratios in a unit square (left square in Fig. 3), and from the new constituent rectangles is possible to get their vertical and horizontal dimensions. The next step is to draw two new unit squares subjected to the media/adventitia (50%/50%) ratio; so the media gets 70% of muscle, 16.2% of elastin and 13.8% of collagen; while the adventitia gets 53.8% of elastin and 46.2% of collagen. Finally, the constituent areas (or proportion) are multiplied by the mixture density, for instance, ρ R0 = 1050[kg/m 3 ], which either correspond to the whole arterial wall or to each layer, Fig. 3.
' & $ % Initialize state variables ρ j R (0) and F F F j gr (0) Set boundary conditions Assemble the internal stiffness matrix K K K int Let the time flow n years • If t = 0. Update state variables (G&R) and internal stiffness matrix x x x(t + 1) = x x x(t) +∆ ∆ ∆ * Compute internal stiffness matrix K K K int and internal forces T T T int * Compute external stiffness matrix K K K ext and external forces F F F ext * Update stiffness matrix K K K and residual R R R  Fig. 3 The first unit square at the left represents the whole arterial wall divided in the mixture constituents, which is split into two new unit squares according to the media/adventitia ratio, each one with its corresponding mixture constituents.