Image-based flow simulation of platelet aggregates under different shear rates

Hemodynamics is crucial for the activation and aggregation of platelets in response to flow-induced shear. In this paper, a novel image-based computational model simulating blood flow through and around platelet aggregates is presented. The microstructure of aggregates was captured by two different modalities of microscopy images of in vitro whole blood perfusion experiments in microfluidic chambers coated with collagen. One set of images captured the geometry of the aggregate outline, while the other employed platelet labelling to infer the internal density. The platelet aggregates were modelled as a porous medium, the permeability of which was calculated with the Kozeny-Carman equation. The computational model was subsequently applied to study hemodynamics inside and around the platelet aggregates. The blood flow velocity, shear stress and kinetic force exerted on the aggregates were investigated and compared under 800 s−1, 1600 s−1 and 4000 s−1 wall shear rates. The advection-diffusion balance of agonist transport inside the platelet aggregates was also evaluated by local Péclet number. The findings show that the transport of agonists is not only affected by the shear rate but also significantly influenced by the microstructure of the aggregates. Moreover, large kinetic forces were found at the transition zone from shell to core of the aggregates, which could contribute to identifying the boundary between the shell and the core. The shear rate and the rate of elongation flow were investigated as well. The results imply that the emerging shapes of aggregates are highly correlated to the shear rate and the rate of elongation. The framework provides a way to incorporate the internal microstructure of the aggregates into the computational model and yields a better understanding of the hemodynamics and physiology of platelet aggregates, hence laying the foundation for predicting aggregation and deformation under different flow conditions.


Introduction
Platelet aggregation is the initial step of thrombus formation, that begins with platelet adhesion to the sites of vascular injury [1][2][3]. It is a multistep dynamic process involving distinct receptors and adhesive ligands such as GPIbα and α IIb β 3 receptors, von Willebrand factor (VWF), fibrinogen and fibronectin [4]. Those receptor-ligand interactions regulate tethering, plateletplatelet cohesion and stabilization of the formed aggregates. Blood flow also plays a critical role in the process as distinct shear ranges lead to unique ways to form aggregates [4,5]. Under low shear conditions, fibrinogen and integrin α IIb β 3 are known to be the predominant factors of platelet aggregation while VWF and fibronectin have progressively increasing contribution to the aggregation with the increase of the flow shear rate [6].
To better understand the impact of the hemodynamics on the formation and mechanics of blood clots, several computational models were developed [7][8][9][10][11][12]. Since the porosity of the forming aggregate influences its formation [13-16], a growing number of computational studies considered the thrombus as a porous medium. Tomaiuolo et al.
[17] established a twodimensional computational model, in which the thrombus consists of two homogeneous porous parts. One with highly activated, densely packed platelets was named core and the other with lower levels of activation was named shell. Different permeabilities were assigned for the core and shell, respectively. They have shown that the solute transport behaviour changes with changes of the platelet packing density. Xu et al. [18] developed a similar novel two-dimensional multi-phase computational model to characterise the interplay between the main components of the clot. They showed that the shear force acting on the clot surface is affected by the internal structure of clots. Also, their results indicated that the porosity of the clots decreased as the flow shear rate increased.
The growing interest in the porous structure of thrombus or platelet aggregate has prompted measurements of their permeability [19,20]. Experimental studies have demonstrated several ways to measure the permeability of clots [20][21][22][23][24]. However, owing to the significant differences in the setup, the composition of the flow and the heterogeneous nature of the formation, the experimental results of permeability usually correspond to unique conditions and are difficult to reuse in different scenarios. At present, the Kozeny-Carman formula [25,26] is one of the most widely applied methods to estimate permeability. It approximates the permeability based on the volume fraction and geometric properties of the solid components that form the porous substance [14,27].
Instead of considering the clot geometry as a simplified object such as half ellipse in the model, the geometry of the clot can be obtained from in vivo or in vitro imaging data, contributing to further studies on the microstructure of platelet aggregates [28]. For example, in [29], a two-dimensional model was developed to calculate the shear stress of the near-thrombus region using the fluorescent labeling images captured in mice. They quantitatively measured the thrombus dynamics in the early stages of hemostasis based on the images, but did not consider the internal structure of the thrombus. Taylor et al. [30] and Pinar et al. [31] examined the dynamics of the blood flow around the thrombus based on the in vitro experiments, but both of them considered thrombi as impermeable solids. In [32], the porous geometry of the clot was reconstructed from the experimental images. This setting is closer to reality compared to the impermeable solid and provides the prediction of the transport of inert solutes in the aggregates.
In this paper, an image-based method to capture the internal microstructure is proposed. Based on this method, a three-dimensional computational model is developed, that can simulate the blood flows through the formed platelet aggregates in a microchannel accurately. In this model, the platelet aggregates are considered as porous media and the microstructure of platelet aggregates is based on in vitro experiments of whole blood perfusion in microfluidic chambers coated with collagen. Two sets of images were used to capture the geometries and internal structures of the aggregates. The images without platelet labeling captured the geometry of the aggregate outline, while the images with platelet labeling are employed to infer the internal density of the aggregates and estimate the permeability. The computational model was subsequently applied to study the hemodynamics inside and around the platelet aggregates. The blood flow velocity, shear stress, kinetic forces of blood flow exerted on aggregates, flow elongational rate and advection-diffusion balance of agonist transport inside the platelet aggregates were investigated and compared under 800 s −1 , 1600 s −1 and 4000 s −1 wall shear rates (WSRs). The proposed model incorporating the internal microstructure of the aggregates provides an accurate estimate of the hemodynamics of platelet aggregates, and thus lays the foundation for predicting further aggregation and deformation under different flow conditions.

In vitro experiments
All the experiments were performed at the Laboratory for Thrombosis Research at KU Leuven. A schematic diagram of the in vitro experiments of platelet aggregates is shown in Fig 1. Fresh human blood was drawn into a citrated tube (3.8% sodium citrate, citrate:blood = 1:9) and mixed with MitoTracker Deep Red (Invitrogen, USA). Citrate chelates extracellular calcium, and without calcium coagulation cannot occur, hence fibrin cannot be formed. A temperature control (Thermo Fisher Scientific, US) was used to ensure the blood was at a human-body temperature. The chip μ−Slide VI 0.1 (Ibidi, Germany) was positioned under the microscope and connected to the test tube with blood and the pump via tubing. Before the experiments, the complete surface of the perfusion chamber was coated using 100 μg/ml Horm-collagen (Takeda, Austria) and stood overnight, resulting in a uniform distribution of collagen fibrils over the complete perfused area. The pump (Harvard Apparatus, US), which was on the other side of the chip, pulled the blood through the chamber with a preset WSR. After five minutes of blood perfusion, for every one micrometer in the z-direction, the images of cross-sections of the channel were captured under an Axio Observer Z1 inverted fluorescence microscope (Zeiss, Germany) using differential interference contrast (DIC) and fluorescence microscopy at 100x magnification. The perfusion time was chosen to ensure after which a stable clot was formed in vitro. The three-dimensional structure of platelet aggregates was reconstructed by stacking all the images over the z-direction later. Two types of images modality were captured for the platelet aggregates. One set of images captured the shape of the aggregate outline with non-labelled platelets (DIC images), while the other employed platelet labelling to infer the internal density (fluorescent images).

Image data processing
A schematic workflow of the methodology employed in this work to transfer the experimental image data to the computational model is provided in Fig 2. The non-labelled images were segmented with 3D Slicer [33] manually and stacked together over the z direction. As a final step to reconstruct the surface of the platelet aggregate, interpolation was applied to fill between the image slices. Smoothing was applied to further improve the surface quality of the platelet aggregate meshes. However, the surface mesh may still contain low quality mesh elements, therefore it was remeshed with MeshLab [34] using the uniform resampling filter. After proceeding with the surface mesh, a corresponding volume mesh was created in Gmsh [35] for use in the numerical simulations.
The intensity of the fluorescence in the labelled image data is proportional to the local density of the platelet aggregate, and can be further used to infer the porosity of the aggregates. Therefore density measurements have been made to approximate the local density of aggregates. This was done by manually counting and the details of it are described in S1 File. Fig 3 shows the results of this measurement which provide us with a mapping between the fluorescence intensity and density of the entire aggregates. The permeability of the aggregates was subsequently inferred using the Kozeny-Carman equation. Further details on this are presented in the next subsection.
To simulate the blood flow around the platelet aggregates, a rectangular domain was created and uniformly meshed to represent a part of the channel of the chip in the experiments, as shown in Fig 4. The platelet aggregates formed on the collagen coated bottom plate, approximately in the middle of the channel. The intensity data of platelet aggregates were interpolated from the platelet aggregate volume mesh to the volume mesh of the flow domain. This interpolation between two meshes allows us to assign the corresponding fluorescence intensity value of the aggregate to a specific location in the simulated domain. Finally, this volume mesh of the flow domain was then used as the spatial discretization of the finite element method to simulate the blood flow.

Computational blood flow model
Blood flow was modelled as an incompressible Newtonian fluid governed by the Navier-Stokes equations, and the influence of the porous clot on the flow was introduced using a Darcy term:  Here, u is the velocity of the fluid, p is the pressure, ρ is the fluid density, and μ is the dynamic viscosity. The last term on the right-hand side accounts for the momentum exchange between the fluid and the solid phases of the porous medium. It is defined by Darcy's law [36,37], which describes flow through a porous medium. In this term, k is the permeability of the platelet aggregates. Permeability is a function of the volume fraction of pore size and fiber or cell/ cell aggregate size [21] and was estimated using the Kozeny-Carman equation Here F s , � and D p are the sphericity of the platelets, the porosity of the aggregates and the platelet diameter, respectively. The complete list of parameter values used in the model and their literature sources are shown in Table 1. According to the inlet flow settings in experiments, the Reynolds numbers for 800 s −1 , 1600 s −1 and 4000 s −1 WSRs is 0.67, 1.37 and 3.33, respectively. We expect that the flow velocity in the platelet aggregates are significantly smaller than the velocity in the channel due to the porous structure, which fits the Stokes flow assumption of the Kozeny-Carman equation. The time-step size was chosen to satisfy the Courant-Friedrichs-Lewy condition. A mesh convergence study was carried out to select the appropriate mesh resolution for the simulation. Considering the computational cost and the convergence of the simulation, a mesh with 3.0 million tetrahedral elements was chosen. The average element size of the mesh is 0.837 μm 3 , which is significantly smaller than the average size of a platelet (7.2 − 11.7 fL) [38]. As shown in Fig 4, the blood flowed through the domain from left to right (along the positive y direction). Since the blood flow domain was considered as a small part of the channel in the x  direction, constant and parabolic velocity profiles were prescribed in the x and z directions respectively at the inlet (for the coordinate system see Fig 4). A constant-pressure condition was prescribed at the outlet. No-slip boundary conditions were imposed at the top and bottom walls. The flow field was solved using the finite element method (FEM) [43] and implemented in FreeFEM [44]. An iterative solver based on flexible generalized minimal residual (FGMRES) method [45] was used and implemented via PETSc [46]. This work was carried out on the Dutch national supercomputer Snelllius (SURF, Netherlands). Each simulation was performed on AMD Rome 7H12 CPU × 2 and parallelized on 128 cores. The corresponding code of the simulation is available at https://github.com/UvaCsl/AggregateFlowSimu/releases/tag/v1.0.0.

Results
The distribution of the fluorescence intensity obtained from the experimental images under three WSRs is presented in Fig 5(a)-5(c). With the increase of the shear rate, the ratio of fluorescence intensity clustering at the tail of the distribution increases, which means more parts in the aggregates have a higher density. Note that the values of this fluorescence intensity are in arbitrary units and not comparable under different scenarios due to the fact that they are influenced by other factors such as initial labeling concentration and light intensity. By computing the corresponding aggregate density from the fluorescence intensity via regression, the average density and the distribution of this density inside the aggregates are demonstrated in To further study how fluorescence intensity, porosity and permeability distribute inside the platelet aggregate, the corresponding results on the cross-sections of the aggregate formed under 1600 s −1 WSR are shown in Fig 6. In the core of the aggregate, the intensity is higher than that in the shell. Correspondingly, the platelets in the core are positioned denser, which leads to much lower permeability in this part. This non-homogeneity results in complex blood flow behavior inside the platelet aggregates.

Flow within and around the platelet aggregate
The steady flow pattern on a cross-section under 1600 s −1 WSR is shown in Fig 7(a). Although the platelet aggregate is considered as a porous medium, the permeability of the aggregates is extremely small. Consequently, most of the blood flows over the aggregate and only a small amount of the blood permeates the platelet aggregate. Fig 7(b) shows the corresponding blood flow velocity inside the platelet aggregates on cross-sections. Although in the shell with relatively high permeability, the blood flow velocity is significantly lower than the extra-thrombus flow, it is still considerably higher than the velocity in the core with lower permeability. However, in the transition from the shell to the core, the magnitude of the velocity is greater than that of some regions of the shell. After blood enters the core, its velocity becomes extremely  low. The minimum velocity in the interior core reaches 10 −11 mm/s, which means barely any blood flow.

Stress analysis of the blood flow
Since the platelet aggregate is considered as a porous medium, the momentum of the flow exerts a kinetic force inside the aggregate, which can lead to partial or complete aggregate embolism [32]. These forces indicate the interaction between blood flow and platelet aggregates under the assumption that the deformation of the aggregate is insignificant, and induce stresses in the aggregate structure. However, these are not discussed in the current work. This kinetic force was calculated by Darcy's law, f ¼ À m k u, and the result is shown in Fig 8(a). The highest forces appear on the top outer layer of the platelet aggregate, due to the relatively high fluid velocity over these parts. Also, similar to the blood velocity in the transition from the shell to the core, there exists an area with higher kinetic force compared to the force in the shell and core.
Furthermore, the local shear stress has a significant influence on platelet activation, aggregation and adhesion [47]. It has been shown in [48] and [49] that high levels of shear stress could initiate platelet aggregation. Fig 8(b) demonstrates the fluid shear stress on the surface of such platelet aggregate under 1600 s −1 WSR condition. These fluid-induced shear stress magnitude σ is estimated by: where _ g is the shear rate magnitude obtained by: Here, D represents the strain rate tensor which is defined as D ¼ 1 2 ðru þ ru T Þ. As expected, shear stress at the top part of the aggregate is several-fold higher than at the bottom. This is because the aggregate protrudes into the center of the channel and is exposed to a higher blood flow velocity. This result agrees with the shear stress distribution reported in [28] and [32].

Advection versus diffusion
Biochemical processes involving different agonists such as adenosine-5 0 -diphosphate (ADP), require a certain amount of time to take place. In the advection-dominated region, if such a process takes more time than the characteristic time of agonist advection (i.e. the average time of the agonist to be convected over a distance of platelet), this process is hindered since the flow velocity is high enough to wash away the agonists before any reaction takes place. Therefore, identifying the advection-dominated region is critical for understanding the dynamics of the aggregation process. The advection-and diffusion-dominated regions of inner chemicals were evaluated by Péclet number, a ratio of advection and diffusion time: where u denotes the magnitude of local flow velocity, D is the mass diffusion coefficient and L denotes the characteristic length. In our simulation, the characteristic length is defined as the platelet diameter, therefore the advection time in the denominator of Pe L indicates the time required to move the distance of a platelet diameter. When the Péclet number equals one, advection and diffusion contribute to mass transport equally. If Pe L > 1, the diffusion time of the chemical is longer than the advection time in such area, which means the motion of the chemical is advection-dominated. Otherwise, it is diffusion-dominated. Three agonists, Calcium (Ca 2+ ), ADP and Factor X, were evaluated in this work (see Table 2). Ca 2+ contributes to multiple stages of cellular activation in platelets [50]. It has the lowest molecular weight of these three chemicals but the highest diffusion coefficient. The second one is ADP, an agonist that plays an important role in platelet activation. Finally, Factor X, an enzyme in the coagulation cascade that can increase the thrombosis propensity [51]. It has the highest molecular weight and the lowest diffusion coefficient among the three agonists. The results demonstrate that the advection-dominated volume in the platelet aggregate increases with the decrease of the diffusion velocities of agonists. Furthermore, Fig 9  visualises the distribution of the advection-and diffusion-dominated areas under the same flow condition. It can be observed that for the highly diffusive agonists only minuscule parts of the platelet aggregate are advection-dominated, even in the shell region. However for Factor X, advection takes over the transport at almost the entire outer layer of the aggregate. In conclusion, with the decrease of the diffusion coefficient, advection gradually dominates the transport of chemicals in the shell region of the platelet aggregate. This aligns well with previous findings predicting such behaviour within the core-shell structure of platelet aggregates [13, 17].

Effects of increasing WSR
Under all three shear rate conditions, the average intrathrombus velocity is at least two orders of magnitude lower than the average extra-thrombus velocity, as shown in Table 3. This is  54]. Furthermore, the increase in inlet blood velocity leads to higher average intrathrombus velocity and stronger kinetic forces on platelet aggregates. In order to explore how the inlet WSRs affect the solute transport, the average Péclet number and advection-dominated volume inside the platelet aggregates are also computed. The result indicates that the average Péclet number inside the aggregates has a positive correlation with WSRs, while the relation between advection-dominated volume and WSRs is not clear.
To further understand the influence of the blood flow velocity on the platelet aggregate formation, the shapes of the platelet aggregates formed under three shear rates are shown in Fig  10. It can be observed that the platelet aggregate formed under a high shear rate condition (4000 s −1 ) has a distinctly different shape. It is the tallest out of the three, while the shape of the aggregate formed under low shear rates is the flattest.
Platelet adhesion has been demonstrated to be strongly influenced by the mechanosensitivity of VWF [55][56][57][58]. This protein uncoils and activates only under certain flow conditions [59]. Previous studies have shown that VWF will unfold when the shear rate or the rate of elongation exceeds a threshold [59,60]. Therefore, high shear rate and elongation rate facilitate the binding of platelets and can lead to platelet activation [61]. The appearing shear rate and elongational rate in the three cases are investigated in Fig 10. The shear rate magnitude _ g is  computed based on the corresponding flow rates: where u, v and w are the flow velocity in the x-, y-and z-directions, respectively. The rate of the uni-axial elongation is defined as the magnitude of the diagonal elements of the rate of the flow strain tensor [62,63]: ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi @u @x According to [59,60,64], the threshold of elongation rate for VWF to unfold is much smaller than that of the shear rate: when the shear rates reach 5000s −1 or elongation flow strain rates are higher than 1000 s −1 , the compact conformation of VWF will change [60,65]. Hence, to observe the region where VWF unfolds in our cases, the color bar maxima for results of the shear rate and the elongational flow are set to coincide with the thresholds reported by these studies, namely 5000 s −1 and 1000 s −1 , respectively. The figure indicates that under a high shear rate condition (4000 s −1 ), the region where the shear rate or the elongation rate exceeds the threshold, is considerably larger than the region in the case of low shear rate flow conditions. This leads to a higher potential for platelets to aggregate and adhere under high shear rates.

Comparison of different ranges of porosity
In this work, the porosity ranges are inferred from the local density of the platelets. The measurement requires one to count the number of platelets manually in multiple regions of the aggregate, which causes uncertainty (for further details see S1 File). Furthermore, aggregates formed under different conditions might display different porosities. Therefore, a series of simulations based on different porosity ranges were performed to quantify the possible influence of potential porosity ranges on the quantities of interest. Different porosity values and different ranges of porosity were chosen and considered to have a linear relation with fluorescence intensity. Fig 11 visualises the changes in average blood flow velocity and advection-dominated volume inside the aggregates in terms of a fixed range of porosity under three WSRs. With the same porosity settings, both the average velocity and the advection-dominated volume are higher under a higher WSR, which is expected. Simultaneously, the increase of the porosity, which means more permeable aggregates, leads to an increase of not only the average velocity but also more advection-dominated areas inside the aggregates. Moreover, the blood flow behaviour under various ranges of porosity was also investigated and the results are shown in Fig 12. With a wider porosity setting, both the average velocity within the aggregates and the effect of advection on chemical transport increase. Overall, these changes in the porosity lead to quantitative changes in the flow field, however, the observed trends remain unaffected, implying the wider generality of the results.

Limitations
The proposed image-based computational model is based on several assumptions and has some limitations. First, the porosity of the aggregates was inferred from the fluorescence intensity of the images. However, the exposure time and exposure level in the experiments can have an impact on the fluorescence intensity. How the fluorescence intensity changes with the exposure time and level should also be investigated in more detail to provide a more precise estimation of the platelet density. Second, the porosity of the platelet aggregates is approximated by the ratio of the area occupied by platelets on the images. This process requires manually counting the number of platelets, which involves uncertainty and contributes to the uncertainty in the model outputs. However, the evaluation of interobserver variability was not performed. Therefore, a comparison study based on different porosity settings was performed, and demonstrated how the potential uncertainty in porosity ranges would affect quantities of interest, such as average velocities and advection-dominated volume of agonists under three WSRs. No qualitative changes have been observed due to the variety of porosity ranges. Furthermore, the volume fraction result shown in Fig 6 follows our expectation on platelet aggregates formed under different shear rates, i.e. the aggregates tend to be denser when formed under higher shear rates. Moreover, in the computational model of blood flow, the viscosity in the lumen was assumed to be constant. However, in actual flow, red blood cells in the microchannel migrate toward the center of the lumen resulting in a red blood cell-free layer close to the vessel wall [66,67]. Such a phenomenon contributes to heterogeneity of viscosity in the lumen and therefore influences the flow resistance and biological transport [68]. Also, the platelets and red blood cells cannot penetrate platelet aggregates, which means the fluid inside the aggregates is pure plasma. Its viscosity is smaller than whole blood [69,70]. A simulation with plasma viscosity has been carried out for comparison. The result (S3 Fig) shows that there is no significant influence to the flow dynamics inside the aggregates. Finally, the platelet aggregate pore wall presents barriers to the agonists, resulting in a smaller agonist diffusion coefficient, which contributes to a larger Péclet number. Since the knowledge of how large the effect of the porous on the diffusion of agonists is unknown, we used the diffusion coefficients that were measured in free flow.

Discussion
In this work, two different microscopy image modalities were used to capture the morphology and microstructure of platelet aggregates formed under in vitro blood perfusion experiments. The information extracted from the images yields the geometry and local porosity of the aggregates, that are subsequently applied in the computational model. The image-based computational model enables an accurate prediction of intrathrombus flow under various flow conditions. A separation of the core and shell of the aggregates in terms of the platelet density can be observed in Fig 6. The core is composed of densely packed platelets, where the transport of agonists is diffusion-dominated, while the outer layer is a loose association of platelets.
Furthermore, the blood flow velocity inside platelet aggregates observed has a qualitative agreement with the result reported in [32]. However, a strip area with relatively high flow velocity was found at the boundary between the shell and the core. To understand what happens here, the magnitude of the blood flow velocity under three WSRs was investigated. As Fig  13 shows, this situation occurs under all three WSRs. This might be because the platelets in the core were excessively dense, causing the circumferential flow around the core and contributing to the increase of the blood flow velocity at the transition zone from shell to core. Correspondingly, the large kinetic forces in the area were generated by the resistance when blood flowed from the shell to the core, a much denser medium. This phenomenon may help to identify the boundary between the shell and the core.
In Table 3, no clear relation is found between advection-dominated volume and WSRs. This is mainly due to the density difference of the aggregates formed under three WSRs. As mentioned before, the platelet aggregate is denser under 4000 s −1 WSR, implying that there are more regions in the aggregate that large proteins, e.g. Factor X, cannot advect by the fluid. Therefore, this leads to a smaller advection-dominated volume for Factor X under 4000 s −1 WSR than that under 800 s −1 WSR. In contrast, for the lighter agonists such as Ca 2+ , the transport is not significantly hampered even if the aggregate becomes relatively denser in the core. As a result, advection dominates more area for the transport of lighter agonists inside the platelet aggregates with the increase of the blood shear rate.
As mentioned before, if the shear rate or the rate of elongation exceeds a threshold, the VWF will unfold and expose the A1 domain, hence supporting platelet attachment. As shown in Fig 10, under high shear rate flow conditions, the shear rate and the rate of elongation of most of the area around the aggregate exceed the threshold, which means the aggregates are exposed to a large amount of unfolded VWF. This might contribute to the easy attachment of VWF to the surface of aggregates, especially on the top of the aggregates. It further explains the phenomenon that the aggregates are prone to grow in parallel to the flow under a low shear rate while the ones formed under a high shear rate show an aggregation tendency perpendicular to the flow direction. The results of the shear rate and the rate of elongation might be further used to predict the aggregation process over time. However, this observation on clot geometries is based on a limited number of experiments. More experiments are required to show if it has statistical significance.
In this work, we considered the aggregates as porous media and mainly focused on the blood flow behaviour and its influence on the aggregates. However the aggregate might deform, contract, or even break away during this dynamic process. A model could be further developed to study the mechanical behaviour of the platelet aggregates. The kinetic force exerted on the platelet aggregates by the blood flow calculated from our model can be applied to this mechanical model. Combined with the experiments, the break-away dynamics of the aggregates and the mechanical stress forces under different environments can be further simulated and analyzed.

Conclusion
In this work, we present an image-based computational blood flow model to simulate the blood flowing through the platelet aggregates which are considered as porous media. The shapes and the microstructures of the platelet aggregates are extracted from microscopy image data of platelet aggregates formed in the blood perfusion experiments. The fluorescence intensity of the platelet labelling is used to infer the porosity of the aggregates and the corresponding permeability is calculated via the Kozeny-Carman equation. The imaged-base model allows us to study the details of hemodynamic transport, such as the velocity, shear stress, and advection-diffusion of agonist transport inside and around the platelet aggregates. Relatively high blood flow velocity and forces were found at the transition from the shell to the core due to the compact distribution of platelets in the interior core. The large forces in the area could contribute to the activation of the platelets in the shell. Furthermore, the results demonstrate that both the flow shear rate and aggregate microstructure have a substantial impact on the transport of agonists. Finally, the shear rate and the rate of elongation flow have also been investigated and the findings imply that there is a strong correlation between these values and the shapes of growing platelet aggregates. Overall, the proposed computational model incorporates the internal microstructure of the aggregates and enables a more precise prediction of the hemodynamics in the platelet aggregates. Our work lays the foundation for predicting shear-dependent aggregation and deformation under different flow conditions, which could subsequently contribute to developing shear-selective anti-thrombotic or anti-platelet drugs.