Abstract
This study proposed an ordinary differential equation model to simulate platelet aggregation by considering the Stokes viscosity resistance and interaction force between each platelet. This model is similar to the car-following traffic model, and its stability condition can be derived via long-wave expansion. The phase diagram can be divided into stable, metastable, and unstable regions. The theoretical results show that the flow viscosity and platelet diameter play essential roles in thrombus formation. The mKdV equation near the critical point is derived to describe the platelet density waves’ evolution properties by applying the reductive perturbation method. Furthermore, through the simulation of the space-time evolution of the distance distributions for platelets, it is shown that a smaller platelet diameter and lower viscosity will lead to thrombus formation, and the analytical results are consistent with the simulation results. The model indicates that thrombosis is similar to traffic congestion.
1. Introduction
Platelets promote wound closure and prevent excessive blood loss when the skin or tissue is injured, but they are also the main factor affecting thrombus formation. If the receptor proteins on a platelet’s surface are activated by platelet-activating factors such as thrombin or ADP, the platelet will firmly adhere to a blood vessel and become the base of a thrombus. Subsequently, a large number of platelets will attract each other and construct a filamentous protein network through platelet-platelet interaction via GPIba-filament protein. At the same time, the platelet membrane will deform and exhibit filopodia, which raises the probability of interconnection due to increased surface area. Abnormal platelet function may lead to bleeding disorders or thrombosis, which can be life threatening. This condition affects more than 2% of the population in the United States [1].
Many studies have carried out numerical simulations of platelets to supplement experimental and clinical results. These simulations typically include processes in the macroscopic scale, coagulation cascade scale, cell scale, and microscopic scale.
The macroscopic scale is equal to the hydrodynamics represented by the commonly recognized blood flow. The characteristic length ranges from 10 to 500 m, including both capillaries and the aorta, where fluid equations are applied. For example, Poziridis [2] used the boundary element method to study the behavior of platelets at the boundary of the blood vessel wall. Wu et al. [3] investigated the flipping of platelets close to the blood vessel via the lattice Boltzmann method. Moreover, capillaries are a significant part of the blood-brain barrier and connect the smallest arteries to the smallest veins. Hence, nanofluid flow through capillary attracts much attention [4–11].
Platelet aggregation corresponds to the coagulation cascade scale. The characteristic length of this scale is near 10 m, so this scale is capable of showing how platelets interact with filamentous proteins in the blood to construct a network structure while enveloping red blood cells and white blood cells and promoting more blood factors to participate in the reaction, finally producing a thrombus. Most studies simplify each platelet as one mass point. Yazdani et al. [12] used the mass point model to simulate thrombus formation and the interactions of platelets with red blood cells under the Morse potential. The phase-field method, percolation theory, and reaction-diffusion equations have also been used to analyze thrombus formation.
The characteristic length of the cell scale is near 1 m, which is close to the size of platelets. This scale focuses on the interaction between platelets and the deformation effect. In one study, the cell scale was used to explore the interactions between platelets, red blood cells, and white blood cells, and it was found that most platelets were distributed near the blood vessel wall [13].
Furthermore, the microscopic scale can be used to investigate receptor and membrane reactions. Many receptor proteins are distributed on the platelet membrane at the molecular level, and biochemical reactions involving these proteins are the fundamental causes of platelet aggregation. Regarding this scale, due to the limitation of the calculation speed of all-atom molecular dynamics, the coarse-grained molecular dynamics (CGMD) method is normally used for research. For example, Zhang [14] conducted a comparative study on platelet surface elasticity and found that the rigid-body model’s surface pressure was higher than the actual results. Prachi [15] proposed the modified Morse potential, which is calibrated with in vitro experiments, to simulate platelet aggregation. Recently, Zhang et al. [16] proposed a new model based on the results of scanning electron microscopy and established a corresponding model to simulate filopodia. Figure 1 and Table 1 show some typical methods and references that are applied to platelet aggregation, and it is easy to observe that the thrombus simulation involves biomedicine, physics, and numerical methods.
According to the above analysis, most numerical simulations for platelets with molecular dynamics require a large amount of calculation and complex integration of multiscale information. Hence, the objective of this study was to develop a simplified model to analyze the thrombosis phenomenon based on the traffic flow model, since traffic congestion and thrombosis have many similar features. For example, a cluster of low-speed cars will lead to traffic congestion; similarly, a thrombus will block blood flow.
This paper is organized as follows. Section 2 introduces Bando’s car-following model and our model considering Stokes’ law and different interaction potential functions. Section 3 presents the stability analysis for this model. Section 4 describes the reductive perturbation method used to derive the nonlinear equation. Section 5 details the numerical simulation performed to validate the theoretical analysis. Section 6 offers the conclusions of this study.
2. Model
In 1998, Bando et al. proposed the optimal velocity model based on ordinary differential equations, which assumed that each driver would set their optimal speed by referring to their distance from the vehicle ahead. The equation of this model is as follows:where is the position of the th vehicle at time and is the distance between the th vehicle and the th vehicle. In addition, the sensitivity coefficient reflects the driver’s reaction time, and the optimal velocity is the driver’s expected velocity. Generally, the optimal velocity is related to the distance between the preceding vehicle and can be written as
As the distance between adjacent vehicles increases, the expected speed increases monotonously with space, but it remains almost straight when the velocity approaches the maximum speed . Therefore, many scholars [28–49] have extended the model to cover scenarios such as intelligent transportation, multiple lanes, and the honk effect.
Here, the optimal velocity function shows the interaction behavior between two cars, which inspires us to convert the optimal velocity function into the interaction force function to analyze the movements of platelets by using Newton’s second law. Then, thrombosis can be treated as traffic congestion because the way in which a thrombus blocks a blood vessel is similar to the way in which a traffic jam blocks a road.
In fluid mechanics, the force formula of an object in fluid (Stokes’ formula of viscous resistance [50]) iswhere is the viscosity coefficient and is the relative velocity. The shape of an unactivated platelet is an ellipsoid, whereas that of an activated platelet is a sphere with filopodium. represents Stokes’ radius, which means the particle with radius has the same dynamics parameters as an irregularly shaped platelet. Hence, platelets can be regarded as fine solids, and Stokes’ law can describe the movement of fine solids in fluid. However, the movement of platelets is different from general free particles in that there is a specific interaction force between platelets. Therefore, after considering the interaction force , we propose the one-dimensional platelet motion model equation as follows:where is the interaction force between platelets and is the mass of one platelet. Yazdani et al. [12] used the Morse potential energy function to study the aggregation behavior of platelets. Therefore, we will also use the Morse force as follows:where is the well depth, is the width of energy, and is the equilibrium distance. For convenience, we set for (5), and the blue curve in Figure 2 represents the Morse force function. By expanding the Morse potential energy function near , one can obtain the standard elastic potential energy, shown by the red dashed line. Meanwhile, we use the hyperbolic tangent format, shown by the black dotted line, corresponding to (2) by approaching the Morse force. In Figure 2, the upper area of each line indicates the repulsive force, and the lower region indicates the attractive force.
3. Linear Stability Analysis
In the ordinary differential equation dynamical system, linear stability analysis has been widely applied to investigate the perturbation response. In other words, stability refers to free flow, which means each platelet moves individually without the cluster. Instability means there exists a thrombus blocking blood flow. In this study, we use the same linear stability analysis as in the classical car-following model to investigate the impact of the evolution of thrombus on blood flow. To do so, we consider blood flow during the stable state of platelet movement. In uniform flow, all the platelets are moving through a blood vessel with the same velocity. The steady-state solution for this equation is as follows:where is the length of the blood vessel and is the number of platelets. In the steady state, the interactive force between each platelet is zero, and all platelets move with the same velocity.
By introducing a small deviation into the steady-state solution, the uniform solution can be written as
Then, we can obtain the linearized equation by substituting (6) and (7) into (4):where and . For convenience, some coefficients in (8) can be replaced with and . By setting for the Fourier models, we can obtainand then we can obtain the first-order and second-order terms of by substituting into (9) as follows:
When , the evolution of uniform steady-state flow still remains stable for long-wavelength modes. In the case of , the flow will become unstable. Finally, we can obtain the neutral stability condition for this simplified dynamical model of platelets as follows:and
To maintain that there exists no thrombosis in a blood vessel, the parameters should satisfy
Here, it can be seen clearly that decreasing viscosity leads to thrombosis. Increasing density with smaller also promotes platelet adhesion, which means thrombosis may occur. These findings are consistent with the findings in [51].
Figure 3 shows the stable region of for different interaction force functions. The lower part below the solid line corresponds to the unstable region, and the upper part over the solid line is the stable region. The region below the red dot-dash line is the infeasible part because has no physical meaning. The space between the coexistence curve and the neutral stability curve is the metastable region mixed with free particles and thrombi.
4. Nonlinear Analysis
The nonlinear behavior near the critical point of the system is analyzed, especially for the function . By introducing a small scaling parameter, the gap between each platelet is set as follows:where is the mean distance between each platelet. For convenience, with and , (4) can also be rewritten as
By expanding each term in equation (16) to the fifth order of , we obtainwhere
By substituting (17)–(21) into (15) and expanding to the fifth order, the partial differential equation can be obtained as follows:
At the critical point , let and take . Then, the second-order and third-order terms of in (22) can be eliminated and written aswhere
To obtain the standard mKdV equation, we introduce the following transformation:and then we obtain the standard mKdV equation aswhere
To determine the wave velocity in (27), we use the following solvability condition:
By integrating (28), we obtain the wave velocity :
Therefore, the standard kink-antikink wave solution of the mKdV equation is
Finally, the kink-antikink wave solution of the platelet spacing can be obtained aswhere the amplitude of the solution can be written as
The kink-antikink solution of the mKdV equation means that there are coexisting phases in the system; in the context of our study, it means a small-scale thrombus and free platelets appear simultaneously. Figure 3 shows the coexistence curve when .
5. Numerical Simulation
To verify the theoretical analysis results and simulate thrombus formation in a one-dimensional system for length and platelet number , we use the fourth-order Runge–Kutta method to solve the equation. Here, we use the periodic boundary condition and set , and we set the initial conditions as follows:
Figures 4–6 show the time-space evolution diagrams corresponding to different functions. Although the types of forces are different, the system becomes increasingly stable with the increase of . According to stability analysis, the phase transition occurs at when . Therefore, the system becomes unstable with the decrease of , so some platelets aggregate, which indicates the formation of a thrombus. For , Figure 6 exhibits the kink-antikink wave structure, a typical solution of the mKdV equation. For the force of , although it shows characteristics of noise, no thrombus exists. Hence, thrombus formation is not related to the viscosity coefficient, and the shape of the platelet is also affected by aspects of the external environment, such as average distribution density.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
Figures 7–9 show the numerical results for the phase space , which verify the stable analysis result. In these figures, we use color to represent when . Here, the step size of each parameter of is 0.05. The dashed lines in Figures 7–9 are consistent with the results of the theoretical analysis, as shown in Figure 3. However, we also find that there exists an unstable area when in Figure 4 because the corresponding stability curve is negative at this region. This result shows that even if the platelet density is low, thrombus formation may still be caused by nonlinear effects, which explains why thrombosis still occurs under the low platelet of vaccine for COVID-19 theoretically, which was reported in [52].
6. Conclusions
In this paper, we introduced additional interaction forces into the Stokes formula of viscous force and analyzed the dynamic mechanism of thrombosis. Theoretical analysis for the car-following model is applied because thrombosis is similar to traffic congestion. In addition, the theoretical analysis results are consistent with the experimental results, indicating that platelet density impacts thrombus formation. Nonlinear analysis shows that in the case of force function similar to optimal velocity function, the model can be reduced to the mKdV equation, and it exhibits the characteristics of kink-antikink wave solutions. Numerical results show that thrombus formation is related to fluid and platelet parameters, and phase-space simulation indicates thrombus formation can still occur when platelet density is low. Some of the conclusions are consistent with experimental data and phenomena. In summary, the one-dimensional model we proposed can effectively simulate thrombosis without requiring excessive calculation.
Data Availability
No data were used to support this study.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This study was supported by the Special Project on High-Performance Computing of the National Key R&D Program (no. 2016YFB0200604) and Nantong Science & Technology Research Plan (no. JC2021133).