TO: Professor Andreas Linninger/ Chih-Yang Hsu
Date: December 2, 2013
Subject: Project 1– The Circle of Willis
Simulation of blood flows in arteries of the Circle of Willis: normal and diseased condition
Ruchita Shah[1], Andreas Linninger[2],Chih-Yang Hsu[3]
*University of Illinois at Chicago-Bioengineering, 12/01/2013
Abstract:
The Circle of Willis is a part of brain which is also known as the cerebral arterial circle. The Circle of Willis provides blood to the different parts of the brain. : Basilar artery (BA), left interior carotid artery (LICA), and right interior carotid artery (RICA), are three main arteries which provide blood to the brain. These arteries supplies total ~750 mL/min blood to the brain[5]. The blockage in one of these arteries may lead to many serious injuries to the brain. The blockage of an artery in brain mainly causes stroke. It is also known as ischemia stroke. In this paper, pressures, resistance and flow rate of blood of the Circle of Willis are determined.In addition, pressure and blood flow for the diseased brain is also determined by minor changes in modeling the Circle of Willis. The diameters of arteries arereduced to analyze the diseased brain. The diameter is reduced by 25%, 50%, 75%, and 100% of normal vessel. More specifically, the diameters of BA, RICA and LICA are changed to study the diseased brain. In this paper, the diameter of one artery is changed at a time. Furthermore, the result shows that decrease in diameter changes the blood flow and pressure of the different arteries in the Circle of Willis. It also shows that more than 75% of accumulation in an artery drastically changes the pressure flow of the different arteries. To understand change in flow with diameter of an artery, this paper mainly focused on LICA, BA and RICA. The effect of change in diameter of one of three artery on other two arteries can be studied.Furthermore, MATLAB simulation is used to predict and understand the brain condition. It does not give the accurate result.
Keywords: BA, LICA, RICA, diameter, blood flow, pressure,stroke
Introduction
The brain is one of the vital organs of human body. About 15-20% of the cardiac output goes into the brain to maintain and regulate normal body functions [1]. The blood is delivered in brain through ring- like vascular structure called Circle of Willis, shown in figure 1.The main arteries in Circle of Willis which supply blood to brain are shown in table 1, with their diameters and the mean flow rates obtained before previous studies.Moreover, BA, LICA, and RICA provides blood to different part of brain from the cardiac output to nourish the brain, whereas LMCA, LACA, LPCA, RMCA, RACA, and RPCA provide blood to the different part of brain. In this paper, BA, RICA and LICA are taken as inflow and LMCA, LACA, LPCA, RMCA, RACA, and RPCA are taken as outflow.
However, blockage of blood supply in the brain causes stroke in brain. This kind of stoke is known as ischemic stroke. There are many diseases and reasons that can cause a stroke, but one of the main reasons is the deposition of cholesterol, also called atherosclerosis, and other reason is hemorrhage. In hemorrhage, due to rupture of vessels the supply of blood to the vessels is restricted in the brain. High blood pressure and aneurysm causes thehemorrhage. There are
two different types of strokes in the brain: acute and chronic. The deposition of cholesterol blocks70-80 % of arteries causes health issues, and sometimes causes a mini-stroke. If there is a 100% blockage of artery it causes permanent damage to the brain, and may lead to death.
Furthermore, the main focus of this paper is on ischemic stroke. Ischemic stroke is most common and 80 % of the strokes are ischemic stroke. There are three different conditionsform the blockage in arteries of brain: the formation of a clot within a blood vessel of the brain, called thrombosis;the movement of a clot from another part of the body such as the heart to the brain, called embolism; or a severe narrowing of an artery in or leading to the brain, called stenosis.
Figure1. The Circle of Willis. Figure shows the major arteries of Circle of Willis [2].
Main Arteries of Circle of Willis with Abbreviations / dmm / P
mmHg / FR
ml/min
Inflow / basilar artery / BA / 3.1 / 86.7 / 131 ± 40
Left interior carotid artery / LICA / 4.8 / 98 / 264 ± 52
right interior carotid artery / RICA / 4.9 / 96 / 252 ± 52
Outflow / left anterior cerebral artery / LACA / 2.3 / 92 / 85 ± 26
left middle cerebral artery / LMCA / 2.8 / 89 / 150 ± 31
left posterior cerebral artery / LPCA / 2.5 / 81 / 66 ± 14
right anterior cerebral artery / RACA / 2.3 / 93 / 80 ± 28
right middle cerebral artery / RMCA / 2.9 / 88 / 145 ± 27
right posterior cerebral artery / RPCA / 2.3 / 82 / 63 ± 14
Table 1. Average diameters (d), pressures (P) at the inlets and outlets, mean flowrates (FR) in the major arteries of Circle of Willis [1].The values of flow rates are expressed as mean±standard deviation.
The main goal of this project is to characterize the relationship between the pressure, resistance and flow of the blood through the Circle of Willis. The normal brain and the diseased brain can be compared and analyzed by modeling, and modeling can be done by MATLAB simulation. The pressure and flow can be analyzed in blocked artery and also other arteries.
Methods
The unknown values of pressures and blood flow of Circle of Willis can be simulated using MATLAB. The Circle of Willis showed as linear system in figure 2
This paper will only focus on the followingarteries BA, LICA, RICA, LMCA, LPCA, LACA, RACA, RMCA,and RPCA. In given MATLAB files and codes, faces, points, diameter and connectivity of faces are known. The number of inflows and outflows can be determined by given files.
Computation of Blood Flow and Pressure
Using the properties that were provided by the network, the unknown variables were determined. In order to calculate the flow rates, pressure changes and the resistances, the inlet and outlet pressures (boundary conditions) and the viscosity had to be determined, and some assumptions had to be made.
The following are the assumption that were made during the experiment: (1) The blood used is in the form of incompressible Newtonianfluid;(2)There is no accumulation in arteries;(3)Only Focusing on main nine arteries as defined earlier for less complications while ignoring all the other arteries;(4) Blood flow is in steady state;(5)Temperature is held constant at 370C;(6) Arteriesare non-deformable throughout the experiment.
Inlet pressure(mmHg) / Outlet pressure
(mmHg)
LICA=91 / LPCA =87
BA =89 / LMCA=87
RICA=92 / LACA=87
RPCA=87
RACA=87
RMCA=87
The network for the Circle of Willis consists of three inlets and six outlets. The inlet and outlet pressures are predicted based on flows given in literatures.The known variables are diameter and connectivityprovided by the network. The values found in literature forthe inlet and outlet points of the circuit are shown in table 2, and The viscosity of the blood is constant which is 4.375 * 10-7 mmHg*min[3].
In order to determine the blood flow through the network, two prime concepts were utilized: (1) Conservation balance equation; (2) Constitutive equations.
According to the conservation of volumetric blood flow, at a node the volumetric flow in must be equal to the volumetric blood flow out:
The above equation can be re-written in recital form as follows:
Before, solving pressures and flow rate, the resistance of each vessel shall be computed. This resistance (α) value can be calculated using Hagen-Poiseuille equation as defined below
Where µ is the viscosity, set as 4.375 * 10-7 mmHg* min, L is the length and D is the diameter of the blood vessels in the network. The resistance is inversely proposal to the diameter of vessel. Which means that the resistance of vessel increases with decreasing in diameter of vessel. Eventually, the change in pressure is computed using constitutive equations:
Where α is the resistance of flow, F is the flow of blood through the faces and P is the change of pressure at the node which connects two or more faces. Accordingly, the change in pressure at any junction is equal to the product of the resistance and the flow of blood through the vessels connected by the junction. Above equation can also be written as
All the conservation balances, constitutive equations and the boundary conditions were then placed into a matrix form and listed symbolically as:
Ax = b (6)
WhereA is a matrix and x and b are vectors. Matrix A is filled by the coefficients of the conservation balances, constitutive equations and the boundary conditions and vector b consists of the resulting values of this equations. Thus, vector x is calculated using the following equation:
x = A\b (7)
Computation of diameter for diseased brain:
For diseased brain, the diameter of vessels are decreased. The diameter of vessels can be decreased by below equation,
(8)
Where diameter of diseased vessel (Dd), actual diameter of vessel (DN), and fraction of decreased diameter (fd). The fraction of decreased diameter can be calculated by,
(9)
Where fa is percent accumulation of cholesterol in vessels.
For diseased condition, the diameters of arteries are narrowed down in the Circle of Willis. The diameter of BA, LICA, and RICA are decreased by 25%, 50%, 75% and 100% by using above equations. The flow and pressure of diseased brain can be evaluated by MATLAB simulation, keeping the boundary condition constant. During simulation, the diameters of BA, LICA, and RICA are decreased one at a time.
Results
The inlet and outlet pressures are kept constant for both condition in brain: normal and diseased. The unknown pressures and flows can be found using MATLAB simulation. The flows and pressures of BA,RICA, and LICA are shown in table 2, for normal condition.In the given network, the degree of freedom is 209, the total unknown values, which is summation of 109 flows and 100 pressure.
The MATLAB simulation provides the unknown pressure and volumetric flow of the normal human brain. The circuit can be varied by changing the boundary condition or diameter of the vessel similar as diseased brain, and can compute the pressure and flows. In this simulation inlet and outlet pressures are kept constant.
Making small change in the code, different simulations of brain can be done. From those data, normal brain and diseasedbrain can be compared and evaluated. In addition, the physiological values can be compared with the simulated value. The table 4.a,4.b and 4.c shows the evaluated values for diseased brain. The Figure 3.a, 3.b and 3.c are physical representation of the table 4.a,4.b and 4.c respectively. Each bar graph shows the physical representation and comparison of the change in flows of all three arteries with change in only one artery.
Figure.3.c. The change of flow LICA,BA and RICA with decrease percentage of diameter in LICA is shown. Each bar represent the flow rate with change in diameter.
Figure.3.b. The change of flow LICA,BA and RICA with decrease percentage of diameter in BA is shown. Each bar represent the flow rate with change in diameter.
Discussion
Using MATLAB simulation, the variation in blood flow with the change in diameter of arteries can be understood. Occurrence of stenosis in an artery causes the decrease in diameter of an artery increases the blood flow, which cause stroke in the brain. Table 4.a, 4.b and 4.c shows that, the decrease in diameters, varies the inflow. In this simulation, the pressures are constant. These table shows the total inflow of the blood.
In table 4.a, the diameter of LICA is changed, and blood flow in the LICA, RICA, and BA are measured. The value measured in table 4.a, shows that with the decrease in diameter of LICA decreasing the bloodflow in LICA, while the blood flow in BA and RICA increases. Furthermore, as shown in table 4.b, decrease in diameter of BA, increases the blood flow in LICA and RICA. When the accumulation is 100% in the BA, the blood flow is distributed evenly in LICA and RICA. Furthermore, when diameter decrease in RICA, the flow of BA and RICAdecreases, and the blood flow in LICA increases. The decrease in diameter may causes the stroke. Depending on the severity of stroke, temporary or permanent injury occurred in the brain.
Moreover, comparing total flow rate of diseased condition given in table 4.a 4.b, 4.c to normal condition given in table 1, and the total flow rate is also decreased with decrease in diameter.The total blood flow in the Circle of Willis remains similar (725±25 mL/min) for 0-50% of the accumulation. Between 50-75% may affect the flow of blood in arteries but may cause temporary damage to the brain. However, for 75% and more accumulation of plaque decreases 20-40% of the total blood flowin the Circle of Willis, which may lead to permanent brain injury.
The results shownin table 4.a,4.b, and 4.c are the possible outcome for different condition in brain. The results are predicted using simulation, but they are not accurate, because of many reasons. One of the main reason, the boundary conditions are manipulated, and kept constant during simulation. Other reason is the chosen boundary conditions are different than the physiological condition, comparing table 1 with table 2 shows the values are very different.Furthermore, the network only focused on main arteries in steady state condition.
Validation
The validation is done by using conservational balance equation. Using viewer, at two or three different node the conservational balance is applied. The computed values of flow were used after simulation. The inflow and outflows are equal at each node which shows that the simulation is true. This validate the code for both normal and diseased brain.
Conclusion
Brain is very complex organ of human body, and it is very hard to measure the pressure and blood flow in it. Using MATLAB, the pressure and flow of the arteries of the Circle of Willis can be predicted. The goal of this paper is to understand the effect of change in diameter on the arteries of the Circle of Willis in diseased brain. The result of diseased brain is compared with the normal brain, and disclosed that the flows are changed in different arteries when the diameter of an artery changes. The result shows that, change of diameter in different arteries, changes the path of the blood flow. The direction and blood flow mainly depend on the occurrence of stenosis in artery. Stenosis on BA, shows that the blood flow of BA equally divides into RICA and LICA, and make the alternative pathway. In addition, the total amount of blood goes into brain is almost same as normal brain. Whereas occurrence of stenosis in LICA or RICA distribution of blood flow is different than BA, and the total amount of blood goes in the brain is lesser than the normal brain.
References
1.Devault, K, and et.al. “Blood Flow in the Circle of Willis: Modeling and Calibration.”
2.Marieb, Elaine N., and Katja Hoehn. Human Anatomy & Physiology. 9th. 2013.eBook.
3.Ischemic Stroke, American Heart Association and American Stroke Association .
4.Alnæs,M . and et.al., “Computation of Hemodynamics in the Circle of Willis.”Stroke .2007(38): 2500-2505
5.Brain Basics : Preventing stroke, National Institute of Neurological Disorder and Stroke.
6. Zhao. and et.al., “Regional Cerebral Blood Flow Quantitative Using MR Angiograph
Appendix
Appendix. I. The below code is used for normal simulation of brain.
%% Normal Simulation
clear all
close all
clc
PIn = [ 91 89 92];
POut = [87 87 87 87 87 87];
[diameterVector, ptCoordMx, faceMx, groupMx, pointMx] = CASEReaderCY('CoW.cs31');
[MU, nFaces, nPoints] = initializeConstants(faceMx,pointMx);
alpha = computeLengthAndResistances(diameterVector,ptCoordMx,faceMx,nFaces,MU);
[A,b] = fillingAandb(pointMx, faceMx, alpha, nFaces, PIn, POut);
[flows, pressures] = solveForFlowsandPressures(A, b, nFaces); %% Baseline
Appendix. II. The below code is used for diseased brain.
%% Disease Simulation
percentage = 0.25;
listofFaces = [1];
diameterVector = changeDiameters(diameterVector, percentage, listofFaces);
alpha = computeLengthAndResistances(diameterVector,ptCoordMx,faceMx,nFaces,MU);
[A,b] = fillingAandb(pointMx, faceMx, alpha, nFaces, PIn, POut);
[flowsD, pressuresD] = solveForFlowsandPressures(A, b, nFaces); %% diseased
function [MU, nFaces, nPoints] = initializeConstants(faceMx, pointMx)
MU=4.375*10^-7; % converting MU value from pascal*s to mmHg*min
nFaces = size(faceMx,1);
nPoints = size(pointMx,1);
function alpha = computeLengthAndResistances(diameterVector,ptCoordMx,faceMx,nFaces,MU)
% solving for the resistance alpha
faceLength = zeros(nFaces,1);
alpha = zeros(nFaces,1);
for i = 1:nFaces % i is start from 1 and ends at length of faceMx
A1=ptCoordMx(faceMx(i,2),1); % assigninng point coordiantes in the
%face to find length
B1=ptCoordMx(faceMx(i,2),2);
C1=ptCoordMx(faceMx(i,2),3);
A2=ptCoordMx(faceMx(i,3),1);
B2=ptCoordMx(faceMx(i,3),2);
C2=ptCoordMx(faceMx(i,3),3);
faceLength(i)=sqrt(((A2-A1)^2)±((B2-B1)^2)±((C2-C1)^2)); % equation of length
alpha(i)=(128*MU*faceLength(i))/((diameterVector(i)^4)*pi*(10^6));
end% end of loop
function [A,b] = fillingAandb(pointMx, faceMx, alpha, nFaces, PIn, POut)
% general form to fill A matrix with zeros
A_Dim = size(pointMx,1)±size(faceMx,1);
A = zeros(A_Dim,A_Dim);
b = zeros(A_Dim,1);
%for counting the input and output
%Filling Matrix A
nEq = 1;
inlet = [];
outlet = [];
%for conservation eq
for i = 1:size(pointMx,1) %Looping rows
if pointMx(i,2) == 0
if pointMx(i,1) < 0
inlet = [inlet i];
else
outlet = [outlet i];
end
else
for u = 1:size(pointMx,2) % looping cols
aFace = pointMx(i,u);
if aFace == 0
break
end
if aFace < 0
A(nEq,abs(aFace)) = -1;
else
A(nEq,abs(aFace)) = 1;
end
end
nEq = nEq±1;
end
end
%for constitutive eq
for k = 1:nFaces
ln1 = faceMx(k,2);
ln2 = faceMx(k,3);
A(nEq,k) = alpha(k);
A(nEq,nFaces±ln1) = -1;
A(nEq,nFaces±ln2) = 1;
nEq = nEq ± 1;
end
% % %filling inlet and outlet
for i = 1:length(inlet)
A(nEq,nFaces±inlet(i)) = 1;
b(nEq,1) = PIn(i);
nEq = nEq±1;
end
for j=1:length(outlet)
A(nEq,nFaces±outlet(j)) = 1;
b(nEq,1)=POut(j);
nEq=nEq±1;
end
function [flows, pressures] = solveForFlowsandPressures(A,b,nFaces)
%Sol=A\b;
flows = Sol(1:nFaces);
pressures = Sol(nFaces±1:end);
function diameterVector = changeDiameters(diameterVector,percentage,listofFaces)
diameterVector(listofFaces) = diameterVector(listofFaces)*percentage;