Numerical and experimental investigations of Fractional Flow Reserve (FFR) in a stenosed coronary artery.

. The coronary artery is numerically investigated based on CFD techniques for measuring the severity of stenosis. In mild cases, medication is often preferred whereas for severe cases surgical intervention is `required but most of the cases fall in between. Thus, it poses a problem for clinicians in choosing an appropriate action. The Fractional Flow Reserve (FFR) is a number which helps to predict the functional significance of stenosis in this scenario. In this study, various cases of occlusion percent ranging between 40 to 70 are considered using different models for predicting FFR in the stenosed coronary artery. The lesion length is also varied between 1 cm and 5 cm based on patient-­specific data. The experimental investigation of FFR in the coronary stenosis using silicon model is also carried out in this study and compared with numerical results. The effect of occlusion percentage and lesion length on the FFR value are quantified .


1.INTRODUCTION
Coronary heart disease (CHD) is the one of disease state which is responsible for approximately 20 percent deaths around the world. The number of CHD related deaths in India has increased from 17% (26% in adults) in [2001][2002][2003] to 23% (32% in adults) in 2010-2013 [1]. Thus, any fast and accurate method of CHD diagnosis will improve its management, including treatment in India and elsewhere.
In the case of severe coronary occlusions, making a decision to carryout angioplasty or more robust repair process is often straightforward. In extremely mild cases of CHD, medication is often preferred. However, in a large number of cases, CHD falls between mild and severe. This poses a problem for clinicians in choosing a right course of action. Fractional flow reserve (FFR) is the ratio of pressure across a coronary stenosis that is used as a clinical measure to determine the severity of CHD [2]. = � is the distal pressure and � is the proximal pressure of the stenotic lesion.The technique followed to diagnose patients for finding out the severity of stenosis is invasive. The procedure involves inserting a catheter into a patient's artery and taking it to the stenosed location for investigation. Eventually the values of pressures are measured across the stenosis quantify the severity of stenosis. But due to its invasive nature, FFR can be associated with medical complications, such as coronary dissections, and is not entirely suitable for the follow-up of medically treated stenoses, in addition to cost implications [3]. Now-days, a Computed Tomography Angiogram (CTA) is performed on patients, which obtains 3D pictures of the moving heart and bigger arterial vessels. This 3D image is reconstructed, and numerical simulations are performed on the generated 3D model of arteries to predict the FFR value. This procedure does not involve any risk as this is non-invasive. Though the 3D modeling predicts FFR value close to that of the value obtained using the catheterization procedure, the time taken for simulations is enormous because it accounts for the actual geometry of the stenosis(non-symmetric) and it also accounts for the coupling of lumped model terminal vessels into the 3D computational domain, which adds to the overall computational cost. Hence several studies focusing on 1D simulations of the heart and arteries were done to reduce the time taken for simulation. The numerical simulation is performed to carry out full arterial tree with diseased state in one of the artery and 1D simulation can predict the severity with low computational cost with fairly good accuracy. The governing equations of blood flow through artery was given by Sherwin et al. [4].A complete model of the 1D arterial circulation by considering interactions of the systemic circulation with coronary circulation provided by Mynard and Nithiarasu [5].They incorporated a heart model in the arterial circulation and carried out simulations. Most of the studies done earlier used 3D model for the calculation of FFR. but the 3D simulations require a lot of computational time, whereas, 1D simulations require comparatively lesser time and computational cost.
The impact of lesion length on functional severity for intermediate coronary lesions are investigated clinically by Iguchi et al [6] .They found that lesion length has physiological significance in intermediate-grade coronary. Several authors have clinically investigated the impact of length length on severity of stenosis but there is no significant numerical study reported. There is a lack of understanding of the intermediate grade stenosis due to impact of the lesion length. There is need of detailed analysis of hemodynamics parameter due to the impact of occlusion percent and lesion length coupled simultaneously.
A 1D model of the arterial network and simulations were done. In this study ,the coupled effect of occlusion percentage and length of lesion on haemodynamics parameters has been studied.The coupled geometric influence of occlusion on FFR is explicitly investigated in this study. The comparison of 1D model is also done with patient specific 3D model.The experimental FFR using patient specific silicon model is quantified and compared with 3D numerical results.

Governing equations
The artery is considered as a cylindrical vessel and its wall considered to be elastic. The 1D equations as given in Sherwin et al. [4] for the continuity and momentum conservation are where is the cross-sectional area, is the mean velocity is the pressure inside the artery, is the density of the blood, is the friction force per unit length and is modelled based upon the assumption of a steady, laminar and Poiseuille flow Three unknown variables are present in the two Equations & (2) for calculating the pressure, the system is closed by providing an equation that relates the pressure to the crosssectional area which is based on elasticity of the vessel by where is the pressure exerted by the surrounding tissues, is the area when transmural pressure is zero (i.e. ) and the material properties of the vessel are accounted by .
The unknown variables ( ) and ( ) are obtained from the values of forward ( ) and backward ( ) characteristics as The Locally Conservative Galerkin (LCG) method is used for the numerical simulation. The LCG numerical technique was introduced by [7], in which each element is treated as a sub-domain with its own boundaries.

Arterial tree modeling
The entire arterial network is modeled with systemic as well as coronary circulation and is shown in Fig.1 .The arteries and their properties for simulation are obtained from Mynard et al [5]. The left epicardial artery is modelled as stenosed left anterior descending along with all parts of tree to be normal.  The Equation representing the stenosed region for a straight duct is given by [8]. This has been modified for a duct which is tapering and the modified expression is

Geometry of stenosed arterial segment
, ��� is the length of the arterial segment, � is the radius at stenosed region, � is the distance from the beginning of the arterial segment to the place where is a maximum constriction of the cross-sectional area, � is the length of the stenosed region in artery, � is the inlet diameter , � is the outlet diameter and severity is the percentage of stenosis that has occurred in the artery.

Validation study
The artery considered for validation is right carotid artery. The trend of flow rate and pressure are observed to be almost identical for normal conditions and functioning of the heart as shown in Fig. 3 and Fig. 4.

Grid independence study
Three grid sizes of 50 nodes, 100 nodes and 150 nodes are taken and simulations are performed. The axial velocity value variation along the length of the artery (Left Endocardial artery) is compared for the three grid sizes in Fig. 5. It is found that both the 100 and 150 grid sizes produce similar results. Hence, 100 nodes per artery is considered for the simulation.

Governing equations
Blood flow is considered as Newtonian fluid and its flow through the blood vessel is assumed to be laminar and incompressible. Here the density of blood is constant so the reduced equation is used in simulation are A finite volume-based solver named PimpleFOAM in OpenFOAM is used to solve the set of equations with second order time, gradient and divergence schemes. The CT scan data obtained from hospital is segmented using snake algorithm and converted into CAD model. It is further meshed in ICEM CFD. One of the CAD model is shown in Fig.6.  The solver is validated using sonu et al. [8] and nondimensional variables are used for validation. The results are found to be showing good agreement with literature results as shown in Fig.7 and Fig.8. The variables used in the figures defined as U=(u/ ��� ) and Y=(y/D). Where U is the non-dimensional velocity, u is the dimensional velocity, ��� is the average velocity given at inlet. D is the diameter at inlet, y and Y are the dimensional and nondimensional cross-stream length respectively.

Boundary Condition
The Womersley velocity profile is obtained from volume flow rate waveform of left anterior descending artery from 1-D simulation. The velocity profiles are calculated from the result of FFT performed on volume flow rate waveform. Womersley velocity profile is given as inlet and for outlet, constant pressure(p=0) boundary condition.

Grid Independence study
The grid independence study is performed for all cases but here it is shown for only one cases. The error is computed with respect to finest mesh that is considered in the simulation. The simulation is performed using 205433 elements in the computational domain.  The experimental setup of flow through coronary artery has been developed in the present study and it is shown in Fig.9. The components used in the experiment are a pulsatile pump, silicon tubes, patient-specific stenosed segment and pressure sensors Water-glycerol at specific ratio is used as blood mimicking fluid [10]. The pressure sensor is calibrated using a two-point calibration technique. The pulsatile pump creates physiological flow relevant for mammalian blood circulation. The patient-specific coronary tree is generated form CT scan data and printed as a 3-D model. The Y-connector is used for insertion of wire catheter and it is taken to the stenosed region for measurement of the pressure in proximal and distal location. Bridge amplifier is used for amplification of the pressure signals. Data acquisition system (DAQ) is used for data visualization of the measured data. Repeatability of the flow measurement was ascertained from experiment carried over a threemonth period.

RESULTS AND DISCUSSION
The flow rate waveform observed in the stenosed coronary artery for various percentage of occlusion is showed in Fig.11. A flow rise is observed in the systole for diseased state condition.The flow rise during systole is reminiscent of the observations of Rammos et al. [11]. This phenomenon is referred in the literature as systolic flow rise.The wave reflections increase with the increase in occlusion of stenosis. The velocity value(u ∝ ( � + � )) is proportional to forward and backward travelling characteristics. As the occlusion increases, the crosssection becomes narrower at the mid-length of the stenosed section causing an increase in the backward traveling wave (w2), resulting in reduction of velocity value.Due to this, the backward travelling characteristics, � will increase with increasing occlusion .Hence the velocity magnitude dependent on the characteristics will become lesser as the percent blockage increases. Thus, the flow rate decreases with increase in occlusion. For the same occlusion , it is seen that with increasing lesion length, the same phenomenon is observed as shown in Fig.12.This observation indicates that not only occlusion percent is important parameter but also lesion length plays a vital role in the intermediate grade artery.    The FFR value is decreasing with increasing occlusion percent and lesion length .For same occlusion percent , the coronary lesion with higher lesion length is more severe and showing a similar trend with clinical investigation [6].The comparison of FFR for full scale 3D numerical simulation with reduced order i.e 1D model is done and is shown in Fig.15. Various parameters like coronary tortuosity, bifurcation and wall deformation may be responsible for deviation between two models.
A correlation is suggested based on numerical simulation results by regression analysis. The results obtained using numerical simulation and correlation are in good agreement and it can be used instead of carrying out 3-D simulation for computing FFR.
= + * + * � + * � + * � + * + * � + ℎ * � + * � + * � where, * is blockage percentage, * is lesion length in mm, a=76.1380, b=-0.9624, c=0.69393, d=-0.1825, e=-6.8064, f=0.24259, g=-0.0043, h=-0.00426, i=0.0000369 and j=-0.0000001262. It is observed that � value for FFR is .9843 for correlation model are adequately significant at 99 confidence. Fig. 9 shows the goodness of the regression model obtained for FFR, where the confrontation between the values predicted by the correlation model and numerical simulation is shown. It is observed that the errors between the numerical values for simulation and the ones given by correlation have been calculated and the deviation from the simulated values is within ± 2.

CONCLUSION
 The flow inside the coronary artery is influenced by backward and forward travelling characteristics.
 The flow rise in systolic phase of cardiac cycle is observed in stenosed artery as compared to healthy state condition.  The pressure drop is increasing with occlusion percent and decreasing with lesion length .  The FFR values continuously decrease as the occlusion percent and lesion lesion length increases.
 The lesion length has a physiological influence on the functional severity of stenosis.  A correlation of FFR is suggested based on statistical analysis are in reasonably good agreement and the results from proposed correlation may be used instead of numerical simulation.  A better prediction of the FFR values can be obtained by considering terminal vessels as windkessel elements as opposed to the resistance model in 1D models considered in the present study.