A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy

ARTICLE IN PRESS

JID: JJBE

[m5G;July 6, 2017;16:0]

Medical Engineering and Physics 0 0 0 (2017) 1–13

Contents lists available at ScienceDirect

Medical Engineering and Physics journal homepage: www.elsevier.com/locate/medengphy

A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy Daniele Bianchi a,∗, Elisabetta Monaldo b, Alessio Gizzi c, Michele Marino d, Simonetta Filippi c, Giuseppe Vairo a a

Department of Civil Engineering and Computer Science (DICII), Universitá degli Studi di Roma “Tor Vergata”, Via del Politecnico 1, Rome 00133, Italy Department of Engineering, Universitá degli Studi “Niccoló Cusano” - Telematica, Roma, Via Don C. Gnocchi 3, Rome 00166, Italy Department of Engineering, Unit of Nonlinear Physics and Mathematical Modeling, University Campus Bio-Medico of Rome, Via A. del Portillo 21, Rome 00128, Italy d Institute of Continuum Mechanics, Leibniz Universität Hannover, Appelstr. 11, Hannover 30167, Germany b c

a r t i c l e

i n f o

Article history: Received 21 February 2017 Revised 5 May 2017 Accepted 1 June 2017 Available online xxx Keywords: Multiscale tissue mechanics Nonlinear material modeling Fluid-structure interaction Abdominal aortic aneurysm Wall shear stress Hemodynamic risk indexes

a b s t r a c t A novel fluid-structure computational framework for vascular applications is herein presented. It is developed by combining the double multi-scale nature of vascular physiopathology in terms of both tissue properties and blood flow. Addressing arterial tissues, they are modelled via a nonlinear multiscale constitutive rationale, based only on parameters having a clear histological and biochemical meaning. Moreover, blood flow is described by coupling a three-dimensional fluid domain (undergoing physiological inflow conditions) with a zero-dimensional model, which allows to reproduce the influence of the downstream vasculature, furnishing a realistic description of the outflow proximal pressure. The fluidstructure interaction is managed through an explicit time-marching approach, able to accurately describe tissue nonlinearities within each computational step for the fluid problem. A case study associated to a patient-specific aortic abdominal aneurysmatic geometry is numerically investigated, highlighting advantages gained from the proposed multiscale strategy, as well as showing soundness and effectiveness of the established framework for assessing useful clinical quantities and risk indexes. © 2017 IPEM. Published by Elsevier Ltd. All rights reserved.

1. Introduction Vascular pathologies are triggered by hemodynamical [1,2] and biochemical factors [3,4]. As a matter of fact, in the context of Abdominal Aortic Aneurysm (AAA), irregular hemodynamics have been shown to accelerate the progression of such a disease [5,6] with a common basis provided by diabetes, male gender, smoke and hypertension [7]. Moreover, the imbalance of tissue biochemical pathways drives the onset of pathological remodeling, and thereby of histological alterations affecting vascular mechanics [8]. For instance, pathological remodeling in aneurysms has been associated with a reduction of elastin content, as well as with the alteration of collagen microstructural and biochemical properties, such as the reduction of the radius in collagen fibers or the density increase of cross-links among collagen molecules [9–11]. In this context, the development of advanced computational models aims at capturing the complexity of the biological tissue physiopathology [12,13] and at providing clinically-useful indica-

∗

Corresponding author. E-mail address: [email protected] (D. Bianchi).

tions [14]. Nevertheless, reliable tools capable of effective predictive diagnostics are still in their infancy stage [15]. Advances have been obtained thanks to the rigorous merging of biomedical imaging and numerical modeling [16,17], though the translation of these approaches into the clinical practice is still ongoing and novel tools for a real-time decision-making protocol are strongly demanded. As a matter of fact, effective computational tools would allow to obtain a reliable assessment of clinical quantities, such as hemodynamic risk indexes based on the measure of the wall shear stress (WSS), for supporting clinicians on decision process towards elective therapeutic approaches [18–21]. Moreover, a suitable in silico platform might be employed for carrying out extensive parametric analyses, able to reproduce different clinical scenarios and to furnish useful indications for accelerating some clinical outcomes. This occurrence might further shade a light on a better understanding of the imbalance of mechanobiological pathways driving the etiology of some diseases. Main challenges in the development of an effective computational framework for vascular physiopathology are inherited by the coupling of mechanical, biochemical and fluid phenomena [15,22]. Indeed, accurate representations of both physiological and patho-

http://dx.doi.org/10.1016/j.medengphy.2017.06.028 1350-4533/© 2017 IPEM. Published by Elsevier Ltd. All rights reserved.

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE 2

ARTICLE IN PRESS

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

logical flow conditions require a fine complete description of vessel wall mechanics and hemondynamics [23]. Among others, specific issues on the tissue constitutive description and on the modeling of the downstream vasculature for a vessel segments are herein recalled. Addressing arterial tissue mechanics, isotropic and/or linearlyelastic assumptions have been employed [24,25], although they bring to inaccurate results [26,27]. More refined constitutive models for cardiovascular tissues are based on hyperelastic nonlinear and anisotropic formulations [28,29]. Nevertheless, these approaches generally exhibit a weak relationship between model parameters and tissue histological/biochemical features, since they are essentially based on phenomenological descriptions of collagen-fiber strain energy, with fiber mean orientation being the only structural information [30,31]. Accordingly, some important microstructural and biochemical alterations experienced in aneurysms could be accounted for only via a suitable calibration of a set of phenomenological parameters, though to be validated and targeted towards patient-specific clinical cases [32,33]. As an alternative, a multiscale constitutive approach has been recently proposed for the description of the mechanical response of soft tissues [16,34–37]. In detail, thanks to the explicit description of the hierarchical collagen organization, this approach introduces model parameters describing well-defined histological and biochemical features only (e.g., radius and crimp amplitude of collagen fibers, density of cross-links among collagen molecules), thereby opening to a straightforward model calibration based on clinical/experimental evidence. As regards, the fluid problem pressure and blood flow within a vessel segment are driven by the tight interplay between wall structural features and hemodynamic boundaries that mimic the rest of the circulation outside the limited computational threedimensional domain [40,41]. For instance, addressing the downstream vasculature, appropriate boundary conditions at the fluid outlet are mandatory in patient-specific preoperative simulations for virtual surgical options [38,39]. In this regard, a number of reduced-order models and different numerical techniques exist [42,43]. One- and zero-dimensional treatments of the outflow boundary conditions have been compared and a different predictive accuracy has been observed [44,45]. On such a basis, a suitable and computationally-effective choice consists in coupling, via a lumped-parameter strategy, the three-dimensional fluid domain with a zero-dimensional model, the latter reproducing the effects of the distal pressure, and accounting for the downstream vascular resistance and capacitance. In other words, following this approach, a multiscale rationale is gained on fluid mechanics, generally referred to as 3D-0D coupling. The main aim of the present work is the development of a novel computational tool that (i) encompasses the state-of-the-art of multiscale cardiovascular modeling, (ii) is amenable for large scale numerical simulations, and (iii) is well suited for direct evaluation of clinical quantities, such as hemodynamic risk indexes. To reach this goal, a fluid-structure interaction (FSI) computational framework, based on a novel flow-tissue multiscale strategy, is presented. As a result, the proposed simulation tool combines the double multiscale nature of vascular physiopathology in terms of both material properties and blood flow. The constitutive model for the arterial tissue is based on the afore-introduced multiscale constitutive rationale proposed and described in Bianchi et al. [16], Maceri et al. [34], Maceri et al. [35], Marino and Vairo [36], Marino and Wriggers [37]. This modeling approach has been proved to be effective for describing anisotropic and non-linear features of arterial tissues. Nevertheless, the employed multiscale constitutive description has been validated on simple case studies [34–36] and applied to patient-specific geometries only by addressing static cases not involving blood flow [16].

Therefore, it has never been employed in a dynamic FSI computational framework, resulting in a significant lack for the applicability of such a refined multiscale constitutive approach in clinics. The fluid-structure interaction is based on an explicit timemarching approach, consisting in a discretized physical coupling between the structural problem and the fluid one in terms of compatibility and equilibrium relationships. As a result, a time discretization of the fluid problem is conceived as a series of subsequent computational steps each of them faced via an Euler formulation. Outflow boundary conditions are obtained from a 3D-0D coupling, based on a three element Windkessel model [46] and accounting for the downstream vasculature [41]. The post-processing phase is conceived in such a way that an accurate assessment of WSS, and hence of some WSS-based clinical indexes (e.g., the Oscillatory Shear Index -OSI- or the ones based on the novel ThreeBand Decomposition approach -TBD- [20]) can be performed. Hence, the proposed computational framework allows to extend some recent contributions on accurate hemodynamic risk indexes [19–21] accounting for a refined tissue modeling. The paper is organized as follows. In Section 2, the multiscale modeling of the blood flow is presented (Section 2.1), together with the multiscale description of the arterial structural mechanics (Section 2.2). The FSI procedure that allows to couple the present double-multiscale strategy, as well as the clinical quantities addressed in the proposed illustrative applications, are introduced in Sections 2.3 and 2.4, respectively. In Section 3, an AAA exemplary computational case study, reproducing a patient-specific geometry obtained from computed tomography (CT) images, is presented and discussed showing the effectiveness of the entire procedure, as well as the advantages gained from the proposed multiscale rationale and the assessment capability of useful clinical quantities. Some future perspectives and concluding remarks are traced in Section 4. 2. Materials and methods The domain of an aortic segment, reconstructed via CT images, is regarded as a continuum body undergoing deformation process during the cardiac cycle, the latter being described by the time variable t ∈ [0, T], where T denotes the duration of a full cardiac cycle and t = 0 identifies the diastolic phase. The current configuration of the aortic segment is denoted by  = (t ), representative of a cylinder-like structure (see Fig. 1). The domain is composed of the arterial tissue, represented by the solid domain s = s (t ), and of blood flow, represented by the fluid domain  f =  f (t ), such that (0 ) = s (0 ) ∪  f (0 ). Denoting by x the position vector of a material point in the current configuration, the geometry of the solid domain is described by the mid-surface s = s (t ), parametrized in x = x (t ), and by the thickness function δ = δ (x , t ) (assumed to be small and slightly variable with x ). Therefore, under the assumption of small local curvature, and with a little abuse of notation, it results s = s × [−δ /2, δ /2]. The boundary of s is denoted by s = ∂ s , with outward normal unit vector ns = ns (x, t ). The surface  s is regarded as s = si ∪ se ∪ s+ ∪ s− , where si (resp., se ) is the lateral internal (resp., external) vessel surface and s± are the end boundaries. The mid-surface  s is endowed with the orthonormal reference system (k, t, n), where n = n(x , t ) is the outward unit normal to  s , t = t(x , t ) is orthogonal to n and it reflects the centerline direction of the aortic segment, k = k(x , t ) is such that k = t × n [16]. With a further little abuse of notation, k, t and n can be considered to denote respectively the circumferential, axial and radial directions (see Fig. 1). Addressing the fluid domain, the boundary of f is denoted by  f = ∂  f , with outward normal n f = n f (x, t ). The surface  f is regarded as  f =  lf ∪  + ∪ − , where  lf = si is the lateral f f

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE

ARTICLE IN PRESS

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

3

Fig. 1. Schematic representation of the adopted multiscale methodology. Starting from biomedical CT images (left), solid and fluid computational domains are reconstructed (center). The solid domain is characterized by multiscale homogenised constitutive features accounting for nano-micro-macro mechanisms. The fluid domain is defined via a 3D-0D multiscale coupling, that allows to adapt the outlet pressure pˆ (t ) in order to account for the influence of the downstream vasculature. The simulation of the flow-tissue interaction allows to perform a reliable risk assessment (right), deriving clinical insights from both microstructural features and wall shear stress.

surface of the fluid domain and  ± are the end (namely, inflow f

and outflow) boundaries. Furthermore, it results ns · n f = −1 on

with initial and boundary conditions respectively given by



 lf .

Vascular mechanics along a cardiac cycle is investigated by coupling a fluid problem for blood flow and a structural problem for the aortic tissue. Addressing the fluid problem (see Section 2.1), the unknown quantities are the velocity v f = v f (x, t ) and the pressure p f = p f (x, t ) fields of the blood flow (i.e., for x ∈ f ). Moreover, introducing xo as the position vector of a material point in the aortic solid configuration at t = 0, the displacement field us = us (xo, t ) of the solid domain (i.e., for xo ∈ s (0)) is the unknown quantity for the structural problem (Section 2.2), with us (xo, 0 ) = 0. By adopting a multiscale strategy, the coupling between the fluid and the structural problem is accounted for by subdividing the overall time-path [0, T] in B equal intervals t long, such that t = T /B, with B a positive integer much greater than one. In particular, at each discrete time value tb+1 = tb + t (with b ∈ {1, . . . , B}), compatibility and equilibrium conditions are prescribed between solid and fluid problems (see Section 2.3). To this aim, at the surface si (t ) =  lf (t ) the fluid pressure field pf (x, t) is employed for prescribing the loading conditions of the structural problem. For the sake of notation, let t¯ = tb (for a certain b ∈ {1, . . . , B}) denote a reference time value and, when necessary, the bar superscript denote known quantities at the reference time t¯.

vf (x¯ , t¯ ) = v¯ f (x¯ ) p f (x¯ , t¯ ) = p¯ f (x¯ )

in  f (t¯ ) , in  f (t¯ )

⎧ QIN (t ) ⎪ ¯ ¯ ¯ ⎪ ⎨vf (x, t ) = − A+ n f (x, t ) v f (x¯ , t ) = 0 ⎪ ⎪ ⎩ p f (x¯ , t ) = pˆ (t )

on  +f (t¯ )

on  lf (t¯ ) . on  −f (t¯ )

(1b)

Here, ∇ (•) and ∇ · (•) are the gradient and divergence operators, Sym(A) is the symmetric part of the second-order tensor A, ρ f is the blood density, μf is the blood viscosity, σ f is the second-order stress tensor, and I the second-order identity tensor. Moreover, v¯ f and p¯ f are velocity and pressure fields at the reference time t¯ (see Section 2.3). Furthermore, QIN is the inlet blood flow, A+ is the area measure of the surface  + (t¯ ), and pˆ (t ) is a f given proximal blood pressure. The latter is determined by taking into account the effects of the downstream vasculature via a zero-dimensional (0D) model (see Fig. 1), built as a three-elements Windkessel circuit [46]. It consists of a proximal resistance (Rp ), a capacitance (C), and a distal resistance (Rd ), and is governed by the following relationship:



pc (t ) dpc (t ) +C Rd dt pˆ (t ) = QOUT (t )R p + pc (t )

QOUT (t ) =

(1c)

where pc (t) is the distal blood pressure and QOUT (t ) =  ¯ − v f · n f dA represents the blood flow through the reference f

¯ − =  − (t¯ ). outflow boundary  f f

2.1. Fluid flow model The fluid problem is governed by the Navier–Stokes equations, herein formulated by assuming blood as an incompressible Newtonian fluid and by neglecting gravity load. Starting from the reference time t¯ and within a given time interval [t¯, t¯ + t], velocity v f = v f (x¯ , t ) and pressure p f = p f (x¯ , t ) fields are obtained, in an Eulerian framework [40], from

⎧ ∂v ⎨ρ f f + ρ f (v f · ∇ )v f − ∇ · σ f = 0 ∂t ⎩∇ · v f = 0 σ f = 2μ f Sym(∇ v f ) − p f I

in  f (t¯ ),

(1a)

2.2. Structural model The deformation of the aortic solid domain for time range t ∈ [t¯, t¯ + t] is described by the displacement field us (xo , t). In the framework of an updated-Lagrangian approach, each time interval t is endowed for the solid problem with an additional internal time variable τ s ∈ [0, t] and is subdivided in equal time-steps dτs = t/Ks , Ks being a positive integer much greater than one. Within a given time interval t, a succession of linearly elastic problems are solved, neglecting any inertial effect. To this aim, let tsk be the discrete time value along the discretized variable τ s and

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

ARTICLE IN PRESS

JID: JJBE 4

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

Fig. 2. Nano-to-micro multiscale rationale employed for modeling mechanics of collagen fibers. Notation.

let xks be the position vector of a material point in the configuration s (tsk ), such that

tsk = t¯ + (k − 1 )dτs ,

xks = xo + us (xo, tsk )

(2a)

with k ∈ {1, . . . , Ks + 1}. Starting from the known displacement field us at the reference time t¯, the displacement field us (xo, t¯ + τs ) along the discretized variable τ s is obtained as

us (xo, tsk ) = us (xo, t¯ ) +

k−1

dus (xsj , tsj ) ,

(2b)

j=1

with k ∈ {1, . . . , Ks + 1} and dus = dus (xks , tsk ) being the incremental displacement from s (tsk ) to s (tsk+1 ). Introducing the corresponding incremental stress and strain tensors, dσ s = dσ s (x, t ) and dD = dD(x, t ), respectively, and starting from an equilibrated and compatible state at t = tsk , the incremental structural problem formulated in terms of dus (x, t) reads as

 ∇ · dσ s (xks , tsk ) = 0 in s (tsk ) , dD(xks , tsk ) = Sym(dH(xks , tsk )) k k k k k k dσ s (xs , ts ) = C(xs , ts ) dD(xs , ts )

(3a)

with boundary conditions expressed by



dσ s (xks , tsk ) n s (xks , tsk ) = dps (xks , tsk ) on si (tsk ) k k k k e k k dσ s (xs , ts ) n s (xs , ts ) = K dus (xs , ts ) on se (tsk ) . dσ s (xks , tsk ) n s (xks , tsk ) = K± dus (xks , tsk ) on s± (tsk )

(3b)

Here, dH(x, t ) = ∇ (dus (x, t )) is the displacement gradient increment, Ke and K± are the tangent stiffness tensors describing the coupling of the arterial segment with the surrounding environment and the outer arterial portions connected to the end boundaries, respectively. Moreover, dps (xks , tsk ) is the incremental surface traction field on the boundary portion si (tsk ), obtained from the coupling with the fluid problem and associated to the fluid pressure field p f (x¯ , t¯ + t ) (see Section 2.3). Finally, the (fourth-order) tangent stiffness tensor C describes the mechanical properties of the aortic tissue and it is representative of a nonlinearly elastic anisotropic behavior. Tensor C is assumed to be constant across the thickness and variable along the mid-surface  s , that is

˜ (x , D (x , t ), S (x , t )) , C ( x, t ) = C ( x + s n ( x , t ), t ) = C

(3c)

with s ∈ [−δ (x , t )/2, δ (x , t )/2] and x ∈  s (t). The equivalent ˜ is obtained as the remid-surface fourth-order stiffness tensor C sult of a multiscale homogenisation approach that explicitly considers the hierarchical arrangement of collagenous constituents. The multiscale formulation allows to obtain the explicit functional ˜ (and hence, of C) on a measure of the current dependency of C deformation state, described by the set D = D (x , t ), and on structural properties in terms of histological, biochemical and biophysical features, described by the set of parameters S = S (x , t ). These

latter may generally vary in time as a result of damage and remodeling effects [47]. For the sake of clarity, the definition of D and S will be given after the presentation of the multiscale rationale. The proposed approach is based on a refined description of collagen fiber mechanics (presented in the following Section 2.2.1) coupled with a lamellar description of the aortic tissue microstructure (described in the following Section 2.2.2). Although wellestablished from previous works [16,34–37,48], main details are herein briefly recalled for the sake of completeness. 2.2.1. Collagen fiber mechanics In the current configuration, the geometry of a crimped collagen fiber is described by means of its centerline shape function f = f (x, εF ), where x is the along-the-chord coordinate and ε F is the fiber nominal along-the-chord strain. Function f(x, ε F ) is assumed to be periodic of period L = L(εF ) = Lo (1 + εF ), where Lo is the length period of collagen fibers in the reference configuration. Let amplitude H = H (εF ) be introduced as a measure of the fiber crimp (see Fig. 2)

H ( εF ) =

max

{ f (x, εF )} .

(4)

x∈[0,L(εF )]

The mechanics of a crimped collagen fiber is described by introducing a curvilinear beam model [37]. At each incremental step from tks (with k ∈ {1, . . . , Ks }), the shape of collagen fibers is updated from the incremental displacement field computed by applying the Principle of Virtual Works [37]. In the reference configuration, fibers are assumed sinusoidal in shape with amplitude Ho , that is

f (x, 0 ) = Ho sin

2π x 

,

Lo

(5)

Therefore, fiber shape function f(x, ε F ) in the current configuration can be fully defined starting from reference period Lo and reference amplitude Ho . Moreover, a further application of the Principle of Virtual Works gives the equivalent along-the-chord tangent stiffness EF = EF (εF ) of the fibrous constituents as [37]



EF (εF ) = E f (ε f ) g1 (εF ) + g2 (εF )/rF2

−1

,

(6)

where rF is the radius of collagen fibers, g1 and g2 are functions accounting for geometric nonlinearities related to collagen crimp [37]

g1 ( εF ) =

g2 ( εF ) =

1 L ( εF )

1 L ( εF )

 L ( εF ) 0

 L ( εF ) 0





∂ 1+ f (x, εF ) ∂x 



2 −1/2

dx ,

∂ f 2 (x, εF ) 1 + f (x, εF ) ∂x

(7)

2 1/2 dx ,

(8)

and Ef represents the tangent modulus of collagen fibrils, which are structures at the mesoscale between micro- and nano-scale and

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

ARTICLE IN PRESS

JID: JJBE

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

whose bundles make up fibers. Fibril tangent modulus depends on the nominal strain ε f , that is E f = E f (ε f ), obtained from the interscale (from micro- to mesoscale) condition

ε f = ε f ( εF ) = L ( εF )/L ( 0 ) − 1 ,

(9)

where L(εF ) is the period curvilinear length of collagen fibers

L ( εF ) =

 L ( εF )



0



∂ 1+ f (x, εF ) ∂x

2 1/2

dx .

(10)

Fibril tangent modulus Ef is obtained by considering the elasticity series associated to collagen elongation and sliding mechanisms, resulting in Marino and Vairo [36]



E f (ε f ) =

1 Em ( εm )

+

Am λc kc m,o

−1

,

(11)

with λc and kc being respectively the density and the stiffness of intermolecular cross-links, Am , m, o and Em the cross-sectional area, the reference length and the tangent modulus of collagen molecules. In particular, m, o is obtained from the contour length

c of collagen molecules and from the length k of molecular kinks, which are labile domains of collagen triple-helices activated by thermal fluctuations [48], resulting in m,o = c − k . The tangent modulus Em = Em (εm ) of collagen molecules depends on their nominal strain ε m . Along the deformation path, ε m is updated starting from the fibril elongation by means of the interscale (from meso to nano-scale) incremental rule [36]

d εm =

E f (ε f ) dε , Em ( εm ) f

(12)

resulting in εm = εm (ε f ). Based on theoretical results from the Worm-Like Chain model and atomistic computations [48], the molecular tangent modulus Em reads as:

Em ( εm ) =

s Em (εms )Emh (εmh )

,

(13)

(εms ) + Emh (εmh ) s and ε h represent molecular nominal strains associated where εm m s Em

with entropy-related (i.e., thermal fluctuations of molecular kinks) and energy-related (i.e., stretching of covalent bonds within triplehelix polypeptides) elongation mechanisms, respectively. Analos and E h are the entropy-related and energy-related tangously, Em m gent moduli of collagen molecules, respectively, given by [48] s Em =

h Em

kB ϑ m,o

p c Am

m,o =

c





3c

s 2[ c − m,o (1 + εm )] 3



+1 ,

(14a)



Eˆ h 1 + exp[−η ( m,oεm / c − εoh )]

+ Eˆo ,

(14b)

with kB and ϑ being the Boltzmann constant and the absolute temperature, respectively, p the persistence length of collagen molecules, Eˆo and Eˆ the low-strain and high-strain collagen tangent moduli, εoh and η the uncoiling strain and resistance of collas and gen triple-helices [48]. Along the deformation path, strains εm h are updated starting from the molecular elongation by means of εm the interscale (from nano to atomistic scale) relationships [48] s d εm =

Em ( εm ) dε , s Em (εms ) m

h d εm =

Em ( εm ) h Em (εmh )

d εm ,

(15)

s = ε s (ε ) and ε h = ε h (ε ). resulting in εm m m m m m In conclusion, by adopting Eqs. (5), (6), (11), (13) and (14), collagen fiber mechanics is fully defined in terms of the set SF of physically-meaningful parameters (see Fig. 2)

SF = {Ho, Lo, rF , λc , kc , M} ,

(16)

with the set M = { k , c , p , Am , Eˆo, Eˆ , η, εoh } collecting molecular biophysical properties.

5

2.2.2. Aortic tangent stiffness Well-established histological evidence allows to regard the aortic tissue as a lamellar structure. In particular, addressing media layer only [16], tissue is described by means of N lamellar units distributed along the tissue thickness [49]. From inner to outer, each media lamellar unit (MLU) is made up by the concentric arrangement of a layer describing the elastin network and NIL layers describing the interlamellar substance. Hence, the total number Na of layers describing the aortic microstructure is given by Na = N (NIL + 1 ). The (fourth-order) tangent stiffness tensor of each layer is denoted by Ci with i = 1, . . . , Na . Let E and I be the sets collecting the indexes associated with the layers representing the elastin network and the interlamellar substance, respectively. Addressing indexes i ∈ E (i.e., the elastin network), the tangent stiffness tensor Ci corresponds to the one of an isotropic material with tangent modulus Ee and Poisson’s ratio νe. On the other hand, addressing indexes i ∈ I (i.e., the interlamellar substance), each layer is modeled as a unidirectional fiberreinforced material, with fibers embedded in a non-collagenous matrix. The reinforcement phase consists in crimped collagen fibers. In the current configuration, the main direction of collagen fibers is helically wrapped along the vessel axis and it is described by the layer-dependent wrapping angle θ i with respect to the circumferential direction. In agreement with arguments recalled in Section 2.2.1, the deformation state of crimped collagen fibers is expressed in terms of their nominal along-the-chord strain ε F, i , where subscript i highlights the layer-dependency. The realignment of collagen fibers along the deformation path is taken into account by updating angle θ i on the basis of the incremental tensor dD = dD (x , t ), which represents the strain tensor increment dD (obtained from the solution of Eq. (3a)) on the mid-surface  s . For the sake of simplicity, dD is here evaluated on the basis of the circumferential strain only, that is dD (x , t ) = dεϕ (x , t )[k(x , t )  k(x , t )] with dεϕ (x , t ) = dD(x , t ) : [k(x , t )  k(x , t )]. Therefore, current fiber orientation θi = θi (x , t ) can be straight obtained starting from the corresponding value θ o, i in the reference configuration. Moreover, the incremental fiber strain dεF,i = dεF,i (x , t ) is given by the interscale (from macro- to micro-scale) relationship [16]:

dεF,i (x , t ) = dD (x , t ) : [ci (x , t )  ci (x , t )] ,

(17)

with ci (x , t ) = cos (θi (x , t ) )k(x , t ) + sin (θi (x , t ) )t(x , t ) being the layer-dependent mean fiber direction. It is worth pointing out that, the structural properties collected in S may vary both layer-by-layer (along the thickness) and with x (along the mid-surface  s ), resulting defined as the collection of layerwise subsets Si = Si (x , t ). Thereby, the tangent stiffness tensor Ci for the ith layer of interlamellar substance is the one of a transversely isotropic material, described by means of the along-the-fiber collagen tangent stiffness EF (ε F, i ) (see Eq. (6)) and the non-collagenous matrix constants, namely Ee and ν e [16]. Therefore, the mechanical properties of the ith composite layer depend on the fiber strain ε F, i and on the set of parameters Si = {θo,i , SF,i , Ee , νe }, where SF,i collects the layer-dependent values associated to the set introduced in Eq. (16) and describing collagen fiber mechanics. Therefore, the dependence Ci = Ci (εF,i , Si ) is clearly highlighted. The mechanics of the aortic tissue is described via the ho˜ , given by mogenised fourth-order stiffness tensor C

˜ (x , D (x , t ), S (x , t )) C =

Na 1 δi (x , t ) Ci (εF,i (x , t ), Si (x , t )) δ ( x , t )

(18)

i=1

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

ARTICLE IN PRESS

JID: JJBE 6

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

where δi = δi (x , t ) is the thickness of the ith layer, and all tensors Ci are considered expressed in a global reference system. Finally, it is worth highlighting that, following the afore-introduced multiscale rationale, it results

D (x , t ) = {εF,i (x , t ) | i ∈ I } ,

S ( x  , t ) = {S i ( x  , t ) | i ∈ I } .

(3) Return to step 2 by enforcing that t¯ = tb+1 and thereby considering that (a) the reference velocity and pressure fields are updated as

v¯ f (x ) = v f (x, tb + t ) ,

(19) 2.3. Fluid-structure coupling The procedure employed for the fluid-structure interaction problem within the proposed flow-tissue multiscale strategy is summarized as follows: (1) Let t¯ = tb be considered for a certain b ∈ {1, . . . , B}. (2) Within a given time interval [t¯, t¯ + t]: (a) on the basis of the known fields v¯ f and p¯ f (whose values are given from the solution of the previous computational step), the fluid problem in Eqs. (1) is solved, by adopting the multiscale description based on the 3D-0D coupling given by Eq. (1c) and allowing to account for downstream vasculature; (b) on the basis of the fluid pressure field p f (x¯ , t¯ + t ) computed from the step 2(a), the displacement field us (xo, t¯ + τs ) in Eq. (2b) is obtained as the solution of the incremental structural problem introduced in Section 2.2. In detail, for any k ∈ {1, . . . , Ks }, Eq. (3a) are solved by employing the incremental surface traction field dps , defined as

p¯ f (x ) = p f (x, tb + t ) , (21)

where fields vf (x, t) and pf (x, t) are computed from the step 2(a); ¯ f = (b) the reference fluid configuration is updated as   f (tb + t ) by using the displacement fields us (x, tb + t ) computed from the step 2(b). To this aim, the updated fluid lateral surface is defined as:

 lf (tb + t ) =

{x|x=xo + us (xo, tb +t ), ∀ xo∈ lf (0 )≡si (0 )}.(22)

An analogous procedure is employed for updating surfaces  ± (tb + t ). Thereby, the updated fluid domain is f simply defined as the volume region within the updated closed boundary  lf ∪  + ∪ − referred to the time tb + f f

t. In this way, the compatibility between solid and fluid regions is naturally verified, resulting in  lf (tb + t ) ≡ si (tb + t ).

Along the time-path described via the variable τ s within [t¯, t¯ + t], the definition of the tangent stiffness tensor C in Eq. (3c) allows to account for microscale and nanoscale mechanisms of collagenous constituents, and hence for the multiscale description of aortic tissue mechanics. At each incremental sub-step dτ s , the solid domain is updated on the basis of the computed incremental displacement field dus , that is (see Eq. (2a))

It is worth observing that, the fluid problem within each time interval t is treated via an Eulerian approach and without considering a strong coupling with the structural problem. Accordingly, the fluid-structure interaction is managed via an explicit time-marching approach, not involving an Arbitrary–Lagrangian– Eulerian formulation and neglecting advective effects induced by the wall moving within the fluid domain at each time interval t. This assumption, formally corresponding to an approximation, can be surely retained as acceptable since local changes in structural configuration induced within the time interval t are expected sufficiently small. Moreover, in order to properly account for tissue nonlinearities, the subincrements dτ s have been introduced for treating the structural problem, allowing to incrementally describe the fluid-pressure loading associated to each time interval t. As initial conditions for the fluid problem, nonuniform pressure and velocity fields are enforced, determining v¯ f (x¯ ) and p¯ f (x¯ ) in Eq. (1b) at t¯ = 0 as the result of a preliminary filling stage. By referring to the boundary conditions introduced in Eq. (1b) and (1c), the latter consists in the simulation of nf cardiac cycles on the CT-based model of the aortic segment considered as rigid (i.e., disregarding any fluid-structure interaction), and starting from an empty domain (namely, starting from null values for blood pressure and velocity in the fluid region). Along the lines of the step 2(b), the pressure field p¯ f (x¯ ) at t¯ = 0 obtained from such a preliminary filling stage is employed for computing the initial deformation state in arterial walls (namely, the pre-strain field in s (0)), and thereby for determining D (x , 0 ) in Eq. (19). As matter of fact, due to the multiscale relationship between deformation state and arterial stiffness (see Eq. (3c)), an heterogeneous distribution of stiffness tensor C is therefore attained at t = 0, even if a constant wall thickness δ and a homogeneous distribution of the structural features S are considered.

s (tsk + dτs ) = {x | x = xks + dus (xks , tsk ), ∀ xks ∈ s (tsk )} ,

2.4. Clinical quantities

dps (xks , tsk ) = ps (xks , tsk + dτs ) − ps (xks , tsk ),

(20a)

where

p s ( x , t ) = p f o ( x ) + [ p f e ( x ) − p f o ( x )]

(t − t¯ ) . t

(20b)

Functions pfo and pfe result from

p f o (xks ) = p f (x¯ (xks ), t¯ ) n f (x¯ (xks ), t¯ ) ,

(20c)

p f e (xks ) = p f (x¯ (xks ), t¯ + t ) n f (x¯ (xks ), t¯ ) ,

(20d)

where, for each position xks at time tsk the corresponding reference postion x¯ at time t¯ is computed as

x¯ (xks ) = xks −

k−1

dus (xsj , tsj ) ,

(20e)

j=1

with k ∈ {1, . . . , Ks + 1}. This is consistent with the fact that, within the time interval [t¯, t¯ + t], pressure field and normal direction n f are referred to the configura¯. tion 

(20f)

s (ts1 )

where the domain (i.e., for k = 1) corresponds to ¯ at t = t¯. Accordthe known reference configuration  ingly, Eq. (20f) allows to compute for k = Ks the updated structural domain s at time t¯ + t = t¯ + Ks dτs , as well as the corresponding boundary s (t¯ + t ).

Starting from the results of FSI simulations, a post-processing phase has been conceived and implemented in order to compute and to monitor several clinical quantities. As a matter of fact, the adopted multiscale constitutive framework allows to predict the local evolution of microstructural tissue features (e.g., fiber crimp and orientation) providing important

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

ARTICLE IN PRESS

JID: JJBE

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

outcomes for properly characterizing local mechanical stimuli affecting tissue cellular environment in physiopathology. In the following exemplary case study, results are shown in terms of the evolution of collagen fiber amplitude along the cardiac cycle. Indeed, alterations in this microstructural feature have been recently recognized as strictly related to the onset of aneurysms [10,11]. It is worth highlighting that, on the basis of Eq. (4) and since εF,i = εF,i (x , t ), fiber amplitude in the current configuration is layer-dependent even if all layers are assumed to be endowed with fibers having identical amplitude in the reference unloaded configuration. Denoting with Hi the layer-specific fiber amplitude, it results Hi = Hi (x , t ) and it is convenient to introduce the average along-the-thickness amplitude H˜ = H˜ (x , t ): Na 1 Hi (x , t ) . Na

H˜ (x , t ) =

(23)

i=1

Moreover, reference will be made to a measure of the circumferential strain ε ϕ , defined as

ε ϕ ( x  , t ) =

1 δ ( x , t )

 t  δ (x ,ζ )/2 0

−δ (x ,ζ )/2

dD(x + s n(x , ζ ), ζ )

: [k ( x  , ζ )  k ( x  , ζ ) ] d s d ζ ,

(24)

representing the average (across-the-thickness) differential (with respect to the diastolic configuration at t = 0) circumferential strain. Results will be also presented in term of the circumferential tangent stiffness Cϕ = Cϕ (x , t )

 S ( x, t ) = 0

S ( x, t ) 0

 S − ( x, t ) =

7

for − Sˆ < S(x, t ) < Sˆ , for S(x, t ) ≤ −Sˆ and S(x, t ) ≥ Sˆ

(29)

for S(x, t ) ≤ −Sˆ . for S(x, t ) > −Sˆ

(30)

S ( x, t ) 0

The closed time support of function S+ (x, t ) is made up by N+ time intervals, where N+ denotes the numerosity, and is computed as the number of time intervals [to+α , t1+α ] ⊆ [0, T ] where S(x, t ) ≥ Sˆ. The union of these N+ time intervals gives the total time duration  + + + + 0 T+ = N α =1 (t1α − t0α ) when S = 0. With obvious notation, let N 0 and T (resp., N− and T − ) be the numerosity and the total duration associated to S0 (x, t) = 0 (resp., S− (x, t ) = 0). It is worth highlighting that, the afore-introduced quantities are functions of the scalar risk factor Sˆ, namely N±,0 = N±,0 (Sˆ) and T ±,0 = T ±,0 (Sˆ). The rich information contained in the dynamic evolution of the WSS signal may be then condensed in two novel synthetic risk predictors, TBDN and TBDT , based on the numerosities N ±,0 (Sˆ) and on the total durations T ±,0 (Sˆ). In particular, considering J equispaced values Sˆ j for the risk factor Sˆ in the range [0, Smax (x)], where Smax (x) is the maximum of S(x, t) over the interval [0, T], the TBDN and TBDT are defined as

J



TBDN (x ) = J j=0

j=0

Sˆ j N− (Sˆ j )

Sˆ j N− (Sˆ j ) + N+ (Sˆ j ) + N0 (Sˆ j )

˜ (x , D (x , t ), S (x )) Cϕ (x , t ) = C

J

: [ k ( x  , t )  k ( x  , t )  k ( x  , t )  k ( x  , t )] ,

(25)

TBDT (x ) = J j=0



j=0

,

(31a)

.

(31b)

Sˆ j T − (Sˆ j )

Sˆ j T − (Sˆ j ) + T + (Sˆ j ) + T 0 (Sˆ j )

˜ resulting from Eq. (18). with C Post-processing phase allows also to compute hemodynamic risk indexes based on the WSS, the latter represented by the timevarying local signal S = S(x, t ), exerted by the blood flow on the arterial wall. Function S(x, t) is computed as the tangential component of the traction vector on  lf in the mean direction of the fluid [40,50], assumed to coincide with the one adopted for the vessel centerline (see Fig. 1), that is:

It is noted that high values of WSS in the direction of the fluid (i.e., high values of functions N+ (Sˆ) and T + (Sˆ)) induces low values of TBDN and TBDT . On the other hand, high values of the indexes in Eqs. (31) are related to the presence of low and retrograde WSS signals (i.e., high values of functions N− (Sˆ) and T − (Sˆ)), in turn linked with the possible deposition and growth of intraluminal thrombi.

S ( x, t ) = −

3. Applications



σ f (x, t )n f (x, t ) · t(x˜  (x, t ), t ) ,

(26)

with σ f obtained from Eq. (1a), x ∈ and function x˜  (x, t ) obtained as the solution of x˜  = x + [δ (x˜  , t )/2]n(x˜  , t ). In the following reference is made to the oscillatory shear index (OSI), which is one of the most used risk indicators, and is defined along the period of the cardiac cycle [0, T] as

 lf ,

OSI(x ) =

   | T S(x, t )dt | 1 1 −  T0 2 |S(x, t )|dt

(27)

0

representing a measure of the time variability of WSS and ranging in the interval [0, 0.5], with the value 0 associated to a unidirectional (healty) shear rate and 0.5 to a highly oscillating (pathological) one. Alternative indexes, have been also proposed in the recent literature, trying to extend the predictive power of OSI, e.g. the residency time [51], the shear index [52], or those based on the three-band decomposition (TBD) approach [19]. The TBD-analysis methodology is characterized by the concept of running threshold that leads to a multi-value segmentation of WSS signal, highlighting directionality effects. Given the WSS signal S(x, t), and a scalar risk factor Sˆ, the TBD analysis defines a triplet of functions [19,20]



S + ( x, t ) =

S ( x, t ) 0

for S(x, t ) ≥ Sˆ , for S(x, t ) < Sˆ

(28)

A computational case study based on the proposed FSI multiscale strategy is presented in what follows. The analysis is conducted on computational domains reconstructed through patientspecific CT images. A finite element discretization is employed for solving fluid and solid problems formulated in Sections 2.1 and 2.2, respectively, exploiting the libraries provided by Comsol Multiphysics (COMSOL Inc., Burlington, MA) [53]. The FSI coupling procedure detailed in Section 2.3 is implemented by means of a home-made Matlab code (MathWorks Inc., Natick, MA) which drives also the enforcement of both the 3D-0D coupling for the fluid problem (see Eq. (1c)) and the multiscale tissue constitutive description for the structural one (see Eq. (3c) and Sections 2.2.1 and 2.2.2). 3.1. Computational case study A portion of a human abdominal aorta has been reconstructed identifying the initial three-dimensional shell configuration  s (0). By employing the Geodesic Active Contours algorithm [54], an image segmentation process is applied, providing a description of the aortic geometry via a surface mesh (see [16] for details). In order to define the solid domain s (0), a surface extrution procedure is then performed by assigning a constant thickness δ = 1 mm over the whole  s (0). Finally, the fluid domain f (0) is obtained via

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE 8

ARTICLE IN PRESS

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

Fig. 3. Computational case study. Frontal (left) and lateral (right) views of the spatial mesh discretization with indication of three control regions (R1: ascending healthy aorta, R2: aneurysm sac, R3: iliac bifurcation) and four control points (P1 to P4). The adopted fluid inflow rate profile (QIN ) is shown at the inlet boundary for a representative beat [21]. The three-element Windkessel circuits (defining the 0D fluid outflow boundary conditions) are also reported at both the LCIA (left common iliac artery) and the RCIA (right common iliac artery), with the corresponding computed outflow pressures pˆ L (t ) and pˆ R (t ), respectively. A sketch of the spring-like mechanical boundary conditions is provided for the inlet, outlet and lateral tissue surfaces, highlighting the corresponding stiffnesses (i.e., K+ , Ke and K− ). Referring to the fluid inflow profile, the control time value t = t ∗ near the systolic peak and adopted for presenting some numerical results is also indicated.

boolean operations and defined as the internal volume delimited by s (0) (Fig. 3). With reference to Fig. 3, the computational domain consists of a healthy descending segment (R1), about 4 cm long and 2 cm in diameter, followed by an aneurysmatic sac (R2), about 4 cm long and 5 cm in diameter. The latter ends with the iliac bifurcation (R3) characterized by left (LCIA) and right (RCIA) portions, about 2 cm long and with an outflow cross section of about 4.7 mm2 (resp., 3.2 mm2 ) for RCIA (resp., LCIA). Fig. 3 shows also a schematic representation of the boundary conditions for fluid and solid problems. The mixed boundary conditions on external tissue surfaces (see Eq. (3b)) are applied by considering stiffness tensors Ke = k¯ e I and K± = k¯ ± I, where k¯ e and k¯ ± are stiffness constant parameters accounting for surrounding tissues and for the arterial portions connected to the end boundaries, respectively. For numerical applications, the values for stiffness parameters k¯ e and k¯ ± have been fine tuned in order to reproduce physiological displacements in the aortic segment, in agreement with the evidence reported by Hallock and Benson [55] in 1937 and discussed by Maceri et al. [34]. Fluid inflow profile QIN is obtained on the basis of in vivo measures, referring to a cardiac cycle having a duration of T = 1 s, and fully in agreement with data adopted by Nestola et al. [21]. Two different 0D Windkessel-type circuits are considered at the fluid outlet surfaces  LCIA and  RCIA , such that  − =  RCIA ∪ f f f f

 LCIA , thus implying the separate numerical solution of Eq. (1c) on f

the two outflow-boundary portions for computing the proximal outflow blood pressures pˆ L and pˆ R , respectively associated with LCIA and RCIA. Following the modeling strategy adopted in Gizzi et al. [20], Nestola et al. [21], the parameters of the 0D Windkessel circuits are fine tuned in order to reproduce physiological pressure levels within the aortic segment. Four selected points along the different regions and at the blood-vessel interface are also identified as control locations for the post-processing analysis: P1 along the healthy segment (i.e., within R1), P2 and P3 within the aneurysm sac (i.e., within R2), and P4 at the iliac bifurcation (i.e., within R3).

Table 1 Values of model parameters for fluid problem employed in numerical simulations [20,21]. Parameter

Symbol

Value

Blood viscosity Blood density Capacitance left common iliac artery Distal resistance left common iliac artery Proximal resistance left common iliac artery Capacitance right common iliac artery Distal resistance right common iliac artery Proximal resistance right common iliac artery

μf [Pa s] ρ f [Kg/m3 ]

0.0035 1050 1 × 10−4 1.35 × 104 600 1 × 10−4 1.65 × 104 800

Cl [cm5 /dyne] Rld [dyne s/cm5 ] Rlp [dyne s/cm5 ] Cr [cm5 /dyne] Rrd [dyne s/cm5 ] Rrp [dyne s/cm5 ]

The aortic tissue is described by structural, biophysical and biochemical parameters Si , chosen in agreement with available data [35,36,48]. In particular, initial fiber orientation pattern is assumed to be the same among different MLUs, resulting univocally described by the NIL -components vector θ F defining fiber orientation within a single MLU. Moreover, structural features collected within the set S (x , t ) (see Section 2.2) are assumed to be uniformly distributed on  s (0) and constant in time (thereby, no damage and remodeling effects are accounted for). It is worth pointing out that structural model parameters have a clear histological meaning and they can be straight deduced from experimental evidence. In this study structural features defining the tissue model have been chosen by referring to histological mean values available in literature, as reported by Bianchi et al. [16]. A sensitivity analysis of the tissue mechanical description to model parameters has been performed in previous studies by some of the authors [35–37], clearly highlighting also how model parameters affect tissue mechanical response when physiopathological conditions are considered [16]. Values of fluid and structural model parameters employed for numerical simulations are summarized respectively in Tables 1 and 2. The preliminary filling stage has been simulated, as detailed in Section 2.3, and by considering n f = 4 cardiac cycles. Addressing TBDN and TBDT risk indices in Eqs. (31), J = 25 discretized values of the scalar risk factor are employed. Some details of the adopted finite element discretization are provided in Fig. 3. In particular, tetrahedral Lagrange quadratic el-

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

ARTICLE IN PRESS

JID: JJBE

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

9

Table 2 Values of model parameters for structural problem employed in numerical simulations [16]. Parameter

Symbol

Value

Surrounding tissues stiffness constant End boundaries stiffness constant MLUs number Number of IL sub-layers per MLU Aortic thickness Collagen volume fraction MLU collagen-fiber orientation vector Collagen fiber period Collagen fiber radius Collagen fiber crimp amplitude Collagen cross-link stiffness density Non-collageneous matrix Young modulus Non-collageneous matrix Poisson ratio Body temperature Collagen molecular cross-sectional area Collagen molecular persistence length Collagen molecular contour length Collagen molecular kinks length Collagen molecular high-strain modulus Collagen molecular low-strain modulus Collagen molecular uncoiling resistance Collagen molecular uncoiling strain

k¯ e [N/cm3 ] k¯ ± [N/cm3 ]

0.1 1 60 6 1 20.0 {0, 20, 40, 0, −40, −20} 5.0 0.04 Lo 0.3 Lo 10.0 1.0 0.49 310 1.41 14.5 287.0 14.0 80.0 1.0 22.5 0.1

N NIL δ [mm] VF [%] θ F [°] Lo [μm] rF Ho λc kc [pN/nm] EM [MPa]

νM

ϑ [K] Am [nm2 ]

p [nm]

c [nm]

k [nm] Eˆ [GPa] Eˆo [GPa]

η εoh

Fig. 4. Time evolution along a cardiac cycle of: (a,b) wall shear stress signal S(x, t); (c,d) average (through-the-thickness) differential circumferential strain ε ϕ (x , t). With reference to the notation introduced in Fig. 3, values at control points P1 to P4 are presented in (a,c), while spatially-averaged quantities within regions R1 to R3 in (b,d). As a notation, y(x) R denotes the average of the function y(x) over the region R ⊆ .

ements based on a pure-displacement formulation have been used for the solid problem, while P1-P1 elements (where both velocity and pressure fields are described by linear shape functions), stabilized by means of an artificial diffusion based on a consistent streamline and crossflow approach, have been adopted for the fluid problem [53]. As a result of a preliminary convergence analysis and as a suitable compromise between accuracy and computational cost, the

coupled domain has been discretized through about 3 · 105 elements, adopting a discretization of the arterial wall consisting in about 4.5 · 104 elements, and a fluid mesh with an average element size less than 2 mm. As regards the time discretization, the overall time-path [0, T] has been subdivided in B = 10 intervals, corresponding to t = 0.1 s (see Section 2). In turn, in order to satisfy the Courant stability requirement, the fluid problem associated to each time interval t has been solved by considering

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE 10

ARTICLE IN PRESS

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

Fig. 6. Hemodynamic risk assessment. (a) Blood streamline patterns at diastole (t = 0, left) and systole (t = t ∗ , right). Frontal and lateral views of the spatial distribution at the fluid-vessel interface for: (b) OSI (Eq. (27)), (c) TBDN (Eq. (31a)) and (d) TBDT (Eq. (31b)).

Fig. 5. Spatial distribution of: (a) average along-the-thickness amplitude H˜ of collagen fibers normalized with respect to Ho at diastole (t = 0, left) and systole (t = t ∗ , right); (b) tangent circumferential stiffness Cϕ at diastole (t = 0, left) and systole (t = t ∗ , right); (c) fluid pressure pf at systole t = t ∗ ; (d) displacement norm us  at systole with respect to the diastolic configuration. The systolic time t = t ∗ is defined in Fig. 3.

adaptative time-steps varying in the range [5−6 , 10−3 ] s. Moreover, the incremental structural problem has been treated by considering Ks = 10, corresponding to dτs = 0.01 s (see Section 2.2). Numerical simulations have been performed by running on a multiprocessor Intel Xeon II HPZ800 workstation with 48 GB of RAM. With reference to the addressed case study, a complete analysis of the fluid-structure interaction problem required a computational cost of about 20 GB of RAM for about 9 hours (with about 60% of the computational cost allocated for solving the incremental structural problem and about 5% for the FSI coupling procedure). 3.2. Numerical results With reference to the exemplary case study previously introduced, some results obtained by solving the FSI problem associated to a cardiac cycle (namely, t ∈ [0, 1] s) via the proposed computational framework, are presented and discussed. In detail, Fig. 3 shows the proximal pressure profiles computed through the adopted 0D Windkessel-type circuits at the outlet RCIA and LCIA boundaries. As expected, different time-dependent pressure patterns are experienced, with a higher pressure peak observed for RCIA than LCIA. Fig. 4(a) depicts the time evolution of the computed WSS signal at the four control points P1 to P4 of the fluid-vessel inter-

face. As expected, moderate values of the shear stress (peak values less than 1 Pa) are obtained along the descending aorta (at P1), while lower, retrogade and oscillating WSS patterns are observed within the aneurysmatic region (namely, at P2 and at P3). At the iliac bifurcation (at P4) a net increase of the WSS peak value is observed, though a physiological trend still occurs (with a small retrogade phase). These features are further confirmed by the spatially-averaged WSS values within the three distinct regions R1 to R3, reported in Fig. 4(b). In particular, a physiological behavior is computed both in R1 and R3, with the highest shear stress peaks experienced at the iliac bifurcation (R3), whereas a totally distorted signal is obtained within the region R2. Fig. 4(c) and (d) provide the quantitative evidence of the average (through the thickness) differential circumferential strain ε ϕ (see Eq. (24)) at the control locations P1 to P4, and of the corresponding spatially-averaged values within the three regions R1 to R3. It is worth noting that, the proposed results highlight the model capability to reproduce the high inhomogeneity of the strain field as strongly dependent on the blood dynamics and on the tissue non-linear response. In detail, the highest variation of the circumferential strain with respect to the diastole (i.e., t = 0) are observed in the healthy segment (R1), confirming the higher compliance of physiological vessel walls with respect to diseased ones (R2). As a matter of fact, since the inhomogeneity of the strain field (also occurring at diastole as a result of the preliminary filling stage, see Section 2.3), collagen fibers in the diseased region (R2) exhibit a reduced crimp in comparison with those in the healthy region (R1), leading to a stiffer response in R2 then in R1. This is fully confirmed also by analysing the results proposed in Fig. 5, where at diastole (t = 0) and near the systolic peak (namely, at t = t ∗ , see Fig. 3) the strict correlation between fiber amplitude and tissue stiffness is clearly shown. In particular, the high spatial heterogeneity of the average along-the-thickness fiber amplitude H˜ (Fig. 5(a), see Eq. (23)), induced by the non-uniform pressure distribution inside the vessel (Fig. 5(c)), leads to similar inhomogeneous patterns of the tangent circumferential tissue stiffness Cϕ (Fig. 5(b), see Eq. (25)), reflecting in turn on the vascular macroscale mechanics in terms of the displacement distri-

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE

ARTICLE IN PRESS D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

bution (Fig. 5(d)). Accordingly, the adopted multiscale constitutive rationale allows to obtain indications on the local evolution of microstructural features within the tissue and on the corresponding evolution of local tissue stiffness. By choosing the healthy zone (R1) as reference and considering results at the systolic peak (t = t ∗ ), regions within the aneurysm sac (R2) exhibit a decrease of about 20% for the average fiber amplitude (corresponding to a higher tangent circumferential stiffness up to about 20%), while significant localization mechanisms associated to a decrease in H˜ of about 40% are observed within the iliac bifurcation (R3, corresponding to higher values of Cϕ up to about 30%). It is worth remarking that, although tissue structural properties are assumed to be uniformly distributed on the CT-based reconstructed tissue domain, the initial deformation state in the arterial wall (namely, the pre-strain field in s (0)) induced by the preliminary filling stage is not uniform in space. Thereby a heterogeneous distribution of the tissue stiffness is also attained at t = 0, since the adopted multiscale approach. This occurrence is clearly confirmed by analysing results depicted in Fig. 5(a) and (b) at diastole, and it is fully in agreement with the evidence reported by Vorp [56]. The proposed computational framework is also able to furnish useful indications on the local blood dynamics inside the vessel. In detail, Fig. 6(a) shows the simulated streamline patterns at diastolic (left) and systolic (right) stages. As expected, large vortex structures occur during diastole within the aneurysm bulge, contributing to a low and retrograde WSS signal. In addition, the pathologic enlargement induces both strong flow deviations and small vortices also during systole, thus contributing to reduce WSS peak values and to introduce shear stress multiple oscillations. These effects further result in an altered hemodynamics at the iliac bifurcation that may experience abnormal mechanical stimuli. Finally, Fig. 6(b–d) provide frontal and lateral views of the spatial distributions at the fluid-vessel interface for the three risk indexes herein considered, i.e. OSI (Eq. (27)), TBDN (Eq. (31a)) and TBDT (Eq. (31b)). As expected, OSI provides a certain risky condition within the aneurysm sac, although fault positive regions are highlighted as risky along the healthy tract R1 (see Fig. 6(b)-right). A more clear risk indication is provided by the combined analysis of the two TBD-based indexes. In particular, the analysis of results associated to TBDN highlights high risk levels within the aneurysm region R2, while TBDT provides a unique significant risky condition between the end of the bulge and the beginning of the iliac bifurcation R3. Such an indication sounds as more robust then the one provided by the OSI-based analysis, since it does not reveal any fault positive area along the healthy tract. 4. Conclusions This work presents a novel computational framework for the FSI analysis of vascular biomechanics based on a novel flow-tissue multiscale strategy. In particular, a multiscale blood flow description, based on a 3D treatment of a vascular segment integrated with a 0D representation of the downstream vasculature, is coupled with a multiscale constitutive model of vessel tissues, furnishing an original and novel contribution to the actual state of the art in this context. The adopted FSI strategy is based on an explicit time-marching approach, solved via a finite-element scheme and consisting in a discretized physical coupling. In detail, a discrete mesh updating is performed to restore fluid-structure compatibility, without specific treatments for advective effects induced by the wall moving within the fluid domain at each time step, that is without considering a rigorous Arbitrary-Lagrangian–Eulerian formulation. Within the limitations of the case study herein numerically addressed, and associated to a patient-specific aortic abdominal aneurysmatic geom-

[m5G;July 6, 2017;16:0] 11

etry, present approach revealed computationally efficient, as well as effective and accurate. Soundness of the proposed approach in predicting realistic vascular macroscale behavior is proven by: the higher outflow pressure peaks obtained for RCIA (right common iliac artery) than LCIA (left common iliac artery), fully in agreement with the evidence reported by Xiao et al. [41] and emphasizing the role of the adopted 3D-0D multiscale modeling strategy for the fluid problem; the evidence that the computed strain-based variation of configuration of the vessel segment reproduces available in vivo observations on arterial wall dynamics experienced during the cardiac cycle [57,58]. The latter outcome is obtained without the need of a case-specific tuning of material parameters. On the contrary, since the non-phenomenological character of the adopted constitutive multiscale description, model parameters are straight set by referring to values of histological, biochemical and biophysical features within well-established experimental ranges. Moreover, present results highlight that the proposed approach allows to accurately assess the mechanical environment at tissue microstructure during cardiac cycles, opening the door to clarify cell-cell signaling pathways which drive tissue remodeling. This occurrence could be exploited for contributing to a better understanding of aneurysm etiology. The methodology is also amenable for the correct evaluation of the shear stress gradient, shown to contribute to pathological endothelial states [2,61]. In this perspective, kinetic-based models [62] of the blood flow can be compared with continuum approaches [63] for the optimal evaluation of endothelial shear stress in cardiovascular pathologies. With respect to existing contributions that adopt such a multiscale constitutive rationale [16,35], vessel mechanics is here predicted on the basis of realistic loading conditions, depending on the bloodstream, coupled with a patient-specific geometrical model. Moreover, recent contributions that consider a multiscale description of blood flow are herein generalized by including a non-linear anisotropic material modeling, with the aim of providing a better assessment of hemodynamic risk factors [20,21]. In particular, presented results show the effectiveness of the proposed strategy in furnishing clinical insights associated to both microstructural features and wall-shear-stress-based risk indexes. Therefore, the proposed computational framework moves towards a refined in silico assessment of vessel mechanics, in terms of reliable quantities to be employed in therapeutic decisions and clinical/surgical planning. Future works will account for microstructural and/or biochemical defects that are associated with the pathological remodeling of aneurysmatic tissues, as well as with possible aneurysm rupture. The employed multiscale constitutive description will allow to straightforwardly address the heterogeneity of the material properties, which plays a key role in the patient-specific variability of aneurismatic features [8]. Moreover, remodeling laws will be defined by considering a chemo-mechano-biological framework, that is by coupling present approach with the description of both molecular diffusion and transport phenomena affecting cell-cell signaling pathways [47]. Furthermore, present approach will be generalized in order to account for dominant damage mechanisms of tissue constituents, in agreement with well-established multiscale strategies [59,60]. Conflict of interest D. Bianchi, E. Monaldo, M. Marino, A. Gizzi, M. Marino, G. Vairo disclose any actual or potential conflict of interest including any financial, personal or other relationships with other people or organizations within three years of beginning the submitted work that could inappropriately influence, or be perceived to influence, their corresponding work.

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE 12

ARTICLE IN PRESS

[m5G;July 6, 2017;16:0]

D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

Acknowledgments D. Bianchi and G. Vairo acknowledge the Italian Minister of University and Research, MIUR (Research program: “Consolidate the Foundations 2015”; Project: BIOART; Grant number (CUP): E82F60 0 0850 0 05). A. Gizzi and S. Filippi wish to thank the partial support of the Italian National Group of Mathematical Physics (GNFM-INdAM). M. Marino acknowledges that this work has been carried out within the framework of the SMART BIOTECS alliance between the Technical University of Braunschweig and the Leibniz University of Hannover. This initiative is financially supported by the Ministry of Economy and Culture (MWK) of Lower Saxony, Germany. References [1] Salsac AV, Sparks SR, Chomaz JM, Lasheras JC. Evolution of the wall shear stress during the progressive enlargement of symmetric abdominal aortic aneurysm. J Fluid Mech 2006;560:19. [2] Chatzizisis YS, Coskun AU, Jonas M, Edelman ER, Feldman CL, Stone PH. Role of endothelial shear stress in the natural history of coronary atherosclerosis and vascular remodeling: molecular, cellular, and vascular behavior. J Am Coll Cardiol 2007;49:2379. [3] Barbour JR, Spinale FG, Ikonomidis JS. Proteinase systems and thoracic aortic aneurysm progression. J Surg Res 2007;139:292–307. [4] Davis FM, Rateri DL, Daugherty A. Mechanisms of aortic aneurysm formation: translating preclinical studies into clinical therapies. Heart 2014;100:1498–505. [5] Nakahashi TK, Hoshina K, Tsao PS, Sho E, Sho M, Karwowski JK, et al. Flow loading induces macrophage antioxidative gene expression in experimental aneurysms. Arterioscler Thromb Vasc Biol 2002;22(12):2017–22. [6] Sho E, Sho M, Hoshina K, Kimura H, Nakahashi TK, Dalman RL. Hemodynamic forces regulate mural macrophage infiltration in experimental aortic aneurysms. Exp Mol Pathol 2004;76(2):108–16. [7] Vardulaki KA, Walker NM, Day NE, Duffy SW, Ashton HA, Scott RAP. Quantifying the risks of hypertension, age, sex and smoking in patients with abdominal aortic aneurysm. Br J Surg 20 0 0;87(2):195–200. [8] Ruddy JM, Jones JA, Spinale FG, Ikonomidis JS. Regional heterogeneity within the aorta: relevance to aneurysm disease. J Thorac Cardiovasc Surg 2008;136. 1123–30 [9] Carmo M, Colombo L, Bruno A, Corsi F, Roncoroni L, Cuttin MS, et al. Alteration of elastin, collagen and their cross-links in abdominal aortic aneurysms. Eur J Vasc Endovasc Surg 2002;23. 543–49 [10] McGloughlin T. Biomechanics and mechanobiology of aneurysms. Berlin Heidelberg: Springer-Verlag; 2011. [11] Robertson AM, Duan X, Aziz KM, Hill MR, Watkins SC, Cebral JR. Diversity in the strength and structure of unruptured cerebral aneurysms. Ann Biomed Eng 2015;43. 1502–15 [12] Fedosov DA, Caswell B, Karniadakis GE. A multiscale red blood cell model with accurate mechanics, rheology, and dynamics. Biophys J 2010;98(10):2215–25. [13] Groen D, Borgdorff J, Bona-Casas C, Hetherington J, Nash RW, Zasada SJ, et al. Flexible composition and execution of high performance, high fidelity multiscale biomedical simulations. Interface Focus 2013;3(2):20120087. [14] Arrieta C, Uribe S, Sing-Long C, Hurtado D, Andia M, Irarrazaval P, et al. Simultaneous left and right ventricle segmentation using topology preserving level sets. Biomed Signal Process Control 2017;33:88–95. [15] Cito S, Mazzeo MD, Badimon L. A review of macroscopic thrombus modeling methods. Thromb Res 2013;131:116. [16] Bianchi D, Marino M, Vairo G. An integrated computational approach for aortic mechanics including geometric, histological and chemico-physical data. J Biomech 2016;49(12):2331–40. [17] Sotelo J, Urbina J, Valverde I, Tejos C, Irarrázaval P, Andia ME, et al. 3d quantification of wall shear stress and oscillatory shear index using a finite-element method in 3d CINE PC-MRI data of the thoracic aorta. IEEE Trans Med Imaging 2016;35(6):1475–87. [18] Bernabeu MO, Nash RW, Groen D, Carver HB, Hetherington J, Krüger T, et al. Impact of blood rheology on wall shear stress in a model of the middle cerebral artery. Interface Focus 2013;3(2):20120094. [19] Gizzi A, Bernaschi M, Bini D, Cherubini C, Filippi S, Melchionna S, et al. Three-band decomposition analysis of wall shear stress in pulsatile flows. Phys Rev E Stat Nonlin Soft Matter Phys 2011;83(3):031902. [20] Nestola MGC, Gizzi A, Cherubini C, Filippi S, Succi S. Novel risk predictor for thrombus deposition in abdominal aortic aneurysms. EPL (Europhys Lett) 2015;112(2):28001. [21] Nestola MG, Gizzi A, Cherubini C, Filippi S. Three-band decomposition analysis in multiscale FSI models of abdominal aortic aneurysms. Int J Mod Phys C 2016;27(02):1650017. [22] Rybicki FJ, Melchionna S, Mitsouras D, Coskun AU, Whitmore AG, Steigner M, et al. Prediction of coronary artery plaque progression and potential rupture from 320-detector row prospectively ECG-gated single heart beat CT angiography: Lattice boltzmann evaluation of endothelial shear stress. Int J Cardiovasc Imaging 2009;25:289–99.

[23] Khan MO, Chnafa C, Gallo D, Molinari F, Morbiducci U, Steinman DA, et al. On the quantification and visualization of transient periodic instabilities in pulsatile flows. J Biomech 2016;52:179–82. [24] Di Martino E, Fumero A, Guadagni G, Ballerini G, Spirito R, Biglioli P, et al. Fluid-structure interaction within realistic three-dimensional models of the aneurysmatic aorta as a guidance to assess the risk of rupture of the aneurysm. Med Eng Phys 2001;23:647–55. [25] Wolters B, Rutten M, Schurink G, Kose U, de Hart JF, van de Vosse F. A patient-specific computational model of fluid-structure interaction in abdominal aortic aneurysms. Med Eng Phys 2005;27:871–83. [26] Auricchio F, Conti M, Morganti S, A R. Patient-specific finite element analysis of carotid artery stenting: a focus on vessel modeling. Int J Numer Method Biomed Eng 2013;29:645–64. [27] Vande Geest J, Schmidt D, Sacks M, Vorp D. The effects of anisotropy on the stress analyses of patient-specific abdominal aortic aneurysms. Ann Biomed Eng 2008;36:921–32. [28] Auricchio F, Ferrara A, Morganti S. Comparison and critical analysis of invariant-based models with respect to their ability in fitting human aortic valve data. Ann Solid Struct Mech 2012;4:1–14. [29] Holzapfel G, Gasser TC, Ogden R. A new constitutive framework for arterial wall mechanics and a comparative study of material models. J Elast 20 0 0;61:1–48. [30] Gasser TD, Ogden RW, Holzapfel GA. Hyperelastic modelling of arterial layers with distributed collagen fibre orientations. J R Soc Interface 2006;3:15–35. [31] Ferrara A, Pandolfi A. Numerical modelling of fracture in human arteries. Comput Methods Appl Mech Eng 2008;11(5):553–67. [32] Vasta M, Gizzi A, Pandolfi A. On three- and two-dimensional fiber distributed models of biological tissues. Prob Eng Mech 2014;37:170–9. [33] Gizzi A, Vasta M, Pandolfi A. Modeling collagen recruitment in hyperelastic bio-material models with statistical distribution of the fiber orientation. Int J Eng Sci 2014;78:48–60. [34] Maceri F, Marino M, Vairo G. A unified multiscale mechanical model for soft collagenous tissues with regular fiber arrangement. J Biomech 2010;43:355–63. [35] Maceri F, Marino M, Vairo G. Age-dependent arterial mechanics via a multiscale approach. Int J Comp Met Eng Sci Mech 2013;14:141–51. [36] Marino M, Vairo G. Computational modeling of soft tissues and ligaments. In: Jin Z, editor. Computational modeling of biomechanics and biotribology in the musculoskeletal system: biomaterials and tissues. woodhead publishing series in biomaterials, vol. 81. Cambridge, UK: Woodhead Publishing/Elsevier; 2014. p. 141–72. [37] Marino M, Wriggers P. Finite strain response of crimped fibers under uniaxial traction: an analytical approach applied to collagen. J Mech Phys Solids 2017;98:429–53. [38] Jansen IGH, Schneiders JJ, Potters WV, van Ooij P, van den Berg R, van Bavel E, et al. Generalized versus patient-specific inflow boundary conditions in computational fluid dynamics simulations of cerebral aneurysmal hemodynamics. AJNR Am J Neuroradiol 2014;35(8):1543–8. [39] Pennati G, Corsini C, Cosentino D, Hsia TY, Luisi VS, Dubini G, et al. Boundary conditions of patient-specific fluid dynamics modelling of cavopulmonary connections: possible adaptation of pulmonary resistances results in a critical issue for a virtual surgical planning. Interface Focus 2011;1(3):297–307. [40] Formaggia L, Quarteroni A, Veneziani A. Cardiovascular mathematics: modeling and simulation of the circulatory system, vol. 1. Springer Science and Business Media; 2010. [41] Xiao N, Humphrey JD, Figueroa CA. Multi-scale computational model of three-dimensional hemodynamics within a deformable full-body arterial network. J Comput Phys 2013;244:22–40. [42] Lagana K, Dubini G, Migliavacca F, Pietrabissa R, Pennati G, Veneziani A, et al. Multiscale modelling as a tool to prescribe realistic boundary conditions for the study of surgical procedures. Biorheology 2002;39(3,4):359–64. [43] Veneziani A, Vergara C. An approximate method for solving incompressible Navier–Stokes problems with flow rate conditions. Comput Methods Appl Mech Eng 2007;196(9):1685–700. [44] Oshima M, Torii R, Tokuda S, Yamada S, Koizumi A. Patient-specific modeling and multi-scale blood simulation for computational hemodynamic study on the human cerebrovascular system. Curr Pharm Biotechnol 2012;13(11):2153–65. [45] Sankaran S, Moghadam ME, Kahn AM, Tseng EE, Guccione JM, Marsden AL. Patient-specific multiscale modeling of blood flow for coronary artery bypass graft surgery. Ann Biomed Eng 2012;40(10):2228–42. [46] Guan D, Liang F, Gremaud PA. Comparison of the windkessel model and structured-tree model applied to prescribe outflow boundary conditions for a one-dimensional arterial tree model. J Biomech 2016;49(9):1583–92. [47] Marino M, Pontrelli G, Vairo G, Wriggers P. Coupling microscale transport and tissue mechanics: modelling strategies for arterial multiphysics. In: Becker S, editor. Modeling of microscale transport in biological processes. Elsevier; 2017. [48] Maceri F, Marino M, Vairo G. Elasto-damage modeling of biopolymer molecules response. Comput Model Eng Sci 2012;87:461–81. [49] O’Connell M, Murthy S, Phan S, Xu C, Buchanan J, Spilker R, et al. The three-dimensional micro- and nanostructure of the aortic medial lamellar unit measured using 3d confocal and electron microscopy imaging. Matrix Biol 2008;27:171–81. [50] Morbiducci U, Gallo D, Cristofanelli S, Ponzini R, Deriu MA, Rizzo G, et al. A rational approach to defining principal axes of multidirectional wall shear stress in realistic vascular geometries, with application to the study of the influence of helical flow on wall shear stress directionality in aorta. J Biomech 2015;48(6):899–906.

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028

JID: JJBE

ARTICLE IN PRESS D. Bianchi et al. / Medical Engineering and Physics 000 (2017) 1–13

[51] Himburg HA, Grzybowski DM, Hazel AL, LaMack JA, Li XM, Friedman MH. Spatial comparison between wall shear stress measures and porcine arterial endothelial permeability. Am J Physiol Heart Circ Physiol 2004;286(5):H1916–22. [52] Corbett SC, Ajdari A, Coskun AU, Nayeb-Hashemi H. Effect of pulsatile blood flow on thrombosis potential with a step wall transition. ASAIO J 2010;56(4):290–5. [53] COMSOL multiphysics®. comsol multiphysics user’s guide. Version 5.2. Stockholm, Sweden: COMSOLAB; 2015. [54] Caselles V, Catte F, Coll T, Dibos F. A geometric model for active contours. Numer Math 1993;66:1–31. [55] Hallock P, Benson IC. Studies on the elastic properties of human isolated aorta. J Clin Invest 1937;16:595–602. [56] Vorp DA. Biomechanics of abdominal aortic aneurysm. J Biomech 2007;40(9):1887–902. [57] Hansen F, Mangell P, Sonesson B, Länne T. Diameter and compliance in the human common carotid artery variations with age and sex. Ultrasound Med Biol 1992;21(1):1–9.

[m5G;July 6, 2017;16:0] 13

[58] Sonesson B, Länne T, Vernersson E, Hansen F. Sex difference in the mechanical properties of the abdominal aorta in human beings. J Vasc Surg 1994;20(6):959–69. [59] Marino M, Vairo G. Influence of inter-molecular interactions on the elasto– damage mechanics of collagen fibrils: a bottom-up approach towards macroscopic tissue modeling. J Mech Phys Solids 2014;73:38–54. [60] Marino M. Molecular and intermolecular effects in collagen fibril mechanics: a multi-scale analytical model compared with atomistic and experimental studies. Biomech Mod Mechanobiol 2016;15:133–54. [61] Dolan JM, Meng H, Singh S, Paluch R, Kolega J. High fluid shear stress and spatial shear stress gradients affect endothelial proliferation, survival, and alignment. Ann Biomed Eng 2011;39(6):1620–31. [62] Pontrelli G, Halliday I, Melchionna S, Spencer TJ, Succi S. The lattice Boltzmann method and multiscale hemodynamics: recent advances and perspectives. IFAC Proc Vol 2012;45(2):30–9. [63] Cherubini C, Filippi S, Gizzi A, Nestola MGC. On the wall shear stress gradient in fluid dynamics. Commun Comput Phys 2015;17(03):808–21.

Please cite this article as: D. Bianchi et al., A FSI computational framework for vascular physiopathology: A novel flow-tissue multiscale strategy, Medical Engineering and Physics (2017), http://dx.doi.org/10.1016/j.medengphy.2017.06.028