Microﬂuidic Simulation and Optimization of Blood Coagulation Factors and Anticoagulants in Polymethyl Methacrylate Microchannels

: Blood coagulation is a critical and complex reaction that involves various chemical substances, such as prothrombin, ﬁbrinogen, and ﬁbrin. The process can be divided into three main steps, namely the formation of the prothrombin activator, conversion of prothrombin to thrombin, and conversion of ﬁbrinogen to ﬁbrin. In this study, an ANSYS simulation is carried out to determine the prothrombin time (PT) of blood, the chemical changes that occur during coagulation and the anticoagulation factor. The addition of deionized water to the microchannels before the addition of blood and reagents results in a two-phase ﬂow. The evaluation of this two-phase ﬂow is necessary, and dynamic simulations are required to determine the PT. The chemical rate constant and order of the chemical reaction are derived from the actual prothrombin time. Moreover, the genetic algorithms in PYTHON and ANSYS are used to estimate chemical reaction parameters for a 20 s PT. The blood and anticoagulant exhibit increased dynamic behavior in the microchannel. In addition, particles are added to the microchannel and the dynamic mesh method is used to simulate the ﬂow behaviors of the red and white blood cells in the microchannel. The PTs for different volumes of blood are also reported.


Introduction
Recent statistics tell us that the aging of the population has already become a major social issue in several countries, making it necessary to monitor the health of older people. An aging population brings with it an increased risk of several illnesses common to older individuals, such as various types of cancer, cardiovascular, bleeding and metabolic diseases, and severe infections [1][2][3][4]. In addition, the severity of such diseases is higher in elders than in younger people. Therefore, the detection of health risk in the early stages is necessary to preserve the health of this portion of the population. Severe bleeding that causes death are bleeding inside the skull (intracranial), subarachnoid, and rupture. Blood coagulation and platelet defects are the causes of bleeding that lead to death [4]. In order to prevent this, understanding the chemical reaction that occurs during blood coagulation and finding the chemical reaction rate involved in the coagulation process are necessary. For this purpose, a new microfluidic device has been devised to detect the blood clotting time (coagulation time) to assist in the early detection of cardiovascular diseases and bleeding disorders.
Several studies have investigated the blood flow in microchannels for drug delivery [5][6][7]. Some studies have examined the microfluidic flow field structures and mixing performance in the microchannels of such devices [8][9][10][11][12]. Jafari et al. reports the optimized parameter of heat sinks fluid flow microchannel using a genetic algorithm, which states the volume flow rate ranges from 0.005 to 0.009 m 3 /s and the entropy generation rate is minimized. Therefore, genetic algorithm can be optimized for unknown parameters to obtain better results [13]. The inertial force of fluid is small due to the diameter of the microfluid channel, which is on the order of 10 −6 m, resulting in a Reynolds number between 0.01 and 100 [8]. Under this Reynolds number condition, laminar flow of the fluid occurs in the microchannel. Fluids can be mixed into this flow only through diffusion. In biomedical testing, increasing the efficiency of the mixing of fluid with the sample and the reagent pretreatment is necessary. In this work, we use species transport and chemical reaction, which involves an unknown reaction rate constant and reaction order, and use a genetic algorithm to calculate these parameters.
Substances involved in the coagulation process are collectively referred to as coagulation factors [14]. Blood coagulation is a complex process which involves 14 different steps that can be roughly divided into three chemical reactions, namely, the formation of a prothrombin activator, prothrombin to thrombin conversion, and fibrinogen to fibrin conversion. Under normal circumstances, the coagulation reaction is rapid, using an internal cascade mechanism to amplify the coagulation signal step by step, and positive feedback links to enhance the coagulation process. In addition, the coagulation mechanisms during hemostasis combine (especially the platelet thrombosis mechanisms) to reinforce and promote each other. To prevent accidental thrombosis under normal conditions and to control the thrombus formed at the bleeding site, the body triggers anticoagulant and fibrinolysis mechanisms against the coagulation mechanism [15,16].
Since generic algorithms have been demonstrated to be an effective technique for optimization and finding better solutions, they are used in this study to calculate the coagulation reaction coefficient. Here, we simulate and calculate the PT chemical reaction rate of the blood, as well as report the point of care (POC) PT time.

Methodology
A microfluidic chip (length 100 mm × width 100 mm) including a microchannel of 1 mm diameter was modeled and analyzed in this work, which was adopted from Manda et al. [17]. The microchannel exhibits two inlet and an outlet, and its substrate is made of polymethyl methacrylate (PMMA). Reagent and blood were used as a reactant for the calculation of PT time and the chemical reaction rate constants in the microchannel. The model of microfluidic chip and its dimension are presented in Figure 1. performance in the microchannels of such devices [8][9][10][11][12]. Jafari et al. reports the optimized parameter of heat sinks fluid flow microchannel using a genetic algorithm, which states the volume flow rate ranges from 0.005 to 0.009 m 3 /s and the entropy generation rate is minimized. Therefore, genetic algorithm can be optimized for unknown parameters to obtain better results [13]. The inertial force of fluid is small due to the diameter of the microfluid channel, which is on the order of 10 -6 m, resulting in a Reynolds number between 0.01 and 100 [8]. Under this Reynolds number condition, laminar flow of the fluid occurs in the microchannel. Fluids can be mixed into this flow only through diffusion. In biomedical testing, increasing the efficiency of the mixing of fluid with the sample and the reagent pretreatment is necessary. In this work, we use species transport and chemical reaction, which involves an unknown reaction rate constant and reaction order, and use a genetic algorithm to calculate these parameters. Substances involved in the coagulation process are collectively referred to as coagulation factors [14]. Blood coagulation is a complex process which involves 14 different steps that can be roughly divided into three chemical reactions, namely, the formation of a prothrombin activator, prothrombin to thrombin conversion, and fibrinogen to fibrin conversion. Under normal circumstances, the coagulation reaction is rapid, using an internal cascade mechanism to amplify the coagulation signal step by step, and positive feedback links to enhance the coagulation process. In addition, the coagulation mechanisms during hemostasis combine (especially the platelet thrombosis mechanisms) to reinforce and promote each other. To prevent accidental thrombosis under normal conditions and to control the thrombus formed at the bleeding site, the body triggers anticoagulant and fibrinolysis mechanisms against the coagulation mechanism [15,16].
Since generic algorithms have been demonstrated to be an effective technique for optimization and finding better solutions, they are used in this study to calculate the coagulation reaction coefficient. Here, we simulate and calculate the PT chemical reaction rate of the blood, as well as report the point of care (POC) PT time.

Methodology
A microfluidic chip (length 100 mm × width 100 mm) including a microchannel of 1 mm diameter was modeled and analyzed in this work, which was adopted from Manda et al. [17]. The microchannel exhibits two inlet and an outlet, and its substrate is made of polymethyl methacrylate (PMMA). Reagent and blood were used as a reactant for the calculation of PT time and the chemical reaction rate constants in the microchannel. The model of microfluidic chip and its dimension are presented in Figure 1.  Computational fluid dynamics [18] are used to simulate and observe the mixing of blood and reagents in the mixed-flow field in a microchannel. The convection diffusion Coatings 2021, 11, 1394 3 of 14 equation and a series of chemical reaction equations for blood coagulation are introduced as the source/sink term for the equation. Moreover, appropriate initial and boundary conditions are applied. A numerical simulation is carried out to analyze the reaction process between the reagent and the blood sample in the microchannel for analysis of the chemical reaction and evaluation of the mixing process.
The genetic algorithm method involves the following steps [19]: (1) generate a set of starting values for a parallel search to avoid the local minimum-or maximum-value loop; (2) use probabilistic selection rules rather than deterministic rules; (3) calculation is mainly performed by the coding of the chromosome rather than the parameters themselves; (4) obtain the fitness score by using the objective function rather than other derivative or auxiliary information. Some studies have suggested using genetic algorithm methods for optimization of the transport behavior of microchannel heat transfer and mass transfer [20,21]. However, since ANSYS Fluent does not provide this feature, we used PYTHON to develop the genetic algorithm. The interactive calculation with ANSYS Fluent is depicted in Figure 2. Computational fluid dynamics [18] are used to simulate and observe the mixing of blood and reagents in the mixed-flow field in a microchannel. The convection diffusion equation and a series of chemical reaction equations for blood coagulation are introduced as the source/sink term for the equation. Moreover, appropriate initial and boundary conditions are applied. A numerical simulation is carried out to analyze the reaction process between the reagent and the blood sample in the microchannel for analysis of the chemical reaction and evaluation of the mixing process.
The genetic algorithm method involves the following steps [19]: (1) generate a set of starting values for a parallel search to avoid the local minimum-or maximum-value loop; (2) use probabilistic selection rules rather than deterministic rules; (3) calculation is mainly performed by the coding of the chromosome rather than the parameters themselves; (4) obtain the fitness score by using the objective function rather than other derivative or auxiliary information. Some studies have suggested using genetic algorithm methods for optimization of the transport behavior of microchannel heat transfer and mass transfer [20,21]. However, since ANSYS Fluent does not provide this feature, we used PYTHON to develop the genetic algorithm. The interactive calculation with ANSYS Fluent is depicted in Figure 2. An efficient dynamic mesh technique [18], which can be used to simulate the problem of flow field variation with respect to time and boundary motion, can be devised using the ANSYS Fluent dynamic mesh model. The motion of the boundary can be rigid or deformed. It can be a predefined motion or a motion that is not defined in advance. During simulation, the mesh update is automatically performed using ANSYS Fluent according to the changes in the boundaries of each iteration step. Dynamic meshes have also been used in several studies to simulate the flow of blood in blood vessels [22][23][24]. An efficient dynamic mesh technique [18], which can be used to simulate the problem of flow field variation with respect to time and boundary motion, can be devised using the ANSYS Fluent dynamic mesh model. The motion of the boundary can be rigid or deformed. It can be a predefined motion or a motion that is not defined in advance. During simulation, the mesh update is automatically performed using ANSYS Fluent according to the changes in the boundaries of each iteration step. Dynamic meshes have also been used in several studies to simulate the flow of blood in blood vessels [22][23][24].

Governing Equations and Boundary Conditions
To prevent the microchannel from being contaminated by impurities and dust, pure water (DI water) is first introduced into the microchannel before the addition of the blood and reagent. This creates a two phase flow. The blood and reagent are mixed in one phase, and the pure water flow is the other phase. Therefore, the Eulerian-Eulerian method is used to simulate the two phases as continuous fluids. The conservation equation is mainly applied to the blood sample and reagent phase [18] (ANSYS 2019).

Equation of Continuity
where the density is mixed by volume weight as follows:

Momentum Equation
where the stress tensor τ i is expressed as follows: The viscosity can be obtained by mass weight mixing as follows:

Species Transport Equations
where the diffusion flux J j and the chemical reaction rate R j are expressed as follows: The general formula for the chemical reaction is as follows: For the aforementioned reaction formula, the reaction rate or rate of formation of the substance k in the reaction r can be derived as follows: According to Seegers et al., blood coagulation can be mainly divided into three main stages [16].
As mentioned earlier, the formation reaction of the prothrombin activator is a complex reaction. The prothrombin activator is assumed to have a constant value in the present study. In addition, the prothrombin activator reacts with a reagent (anticoagulant). Therefore, the chemical reaction can be simplified as follows: Prothrombin Activator The chemical reaction rate can be expressed as follows: where R is a reagent (anticoagulant); PA is the prothrombin activator; T is thrombin; and FG is fibrinogen. The objective function is used to estimate the chemical reaction parameters. Therefore, in this study, the objective function is set as the minimum value of the average mass fraction of the prothrombin activator and fibrin at the end of the microchannel, as presented in Equation (18) [18].
where f l is an objective function coefficient corresponding to the i th variable; dA is the average mass fraction of the prothrombin activator and fibrin; and wl is the weighting mass fraction of the prothrombin activator and fibrin. The following parameters are selected for the genetic algorithm: population size of 100; chromosome length of 10; crossover probability of 0.75; and mutation probability of 0.05. In addition, to simulate the flow of red blood cells (RBC) and white blood cells (WBC) inside the microchannel, particles are added to the flow field. Moreover, the transport behavior of the particles in the microchannel is observed using the dynamic mesh technique and Newton's second and third laws of motion. In this study, to reduce the convergence time and complexity of the simulation, the particles are divided into three groups. Each group has approximately 3-10 RBCs, and the groups do not affect each other.

Boundary Equations
Since the microfluidic channel is filled with DI water at the beginning of the process, if the inlet flow rate is less than a certain value, blood or reagent cannot flow in the microfluidic channel. Therefore, the inlet speed of the blood and reagent is set to 0.005 m/s. The condition between the microchannel and the wall is assumed to be no-slip. In addition, the prothrombin activator and fibrinogen contents in the blood are assumed to be 0.01 and 0.04 wt%, and the other inlet is a reagent, respectively.

Results and Discussion
The reference value for the prothrombin time (PT) varies depending on the analytical method used but is generally approximately 12-13 s. For subjects who have not received anticoagulation therapy, the international normalized ratio (INR) reference value is 0.8-1.2. If anticoagulant drugs are used, the INR value can be 2-3 [25]. Therefore, in this study, a dynamic simulation is used to investigate changes in the blood and anticoagulant in the microchannel. Figure 2 shows a flow chart of the genetic algorithm. We first discuss the distribution of the mass fraction for each item in the flow channel with various k 1 -k 3 and n 1 -n 5 values after 20 s, because different k 1 -k 3 and n 1 -n 5 values affect the prothrombin activator, thrombin, and fibrin contents in the microchannel. The small area in the panel on the left-hand side of Figure 3 indicates that the blue zone contains the blood inlet side and mixing area. The brown zone is the reagent inlet side, which indicates that the mass fraction of the prothrombin activator decreases because the k 1 and n 1 values are small. The rate of decrease of n 1 is higher than that of k 1 . Figure 3 shows that initially the prothrombin activator in blood side and mass fraction is 100% before the reaction; then, prothrombin gradually descends in the microfluidic channel. At x = 0, it immediately reacts with the reagent, so the reaction is complete before x < 5. Since the prothrombin reaction is related to R1 (Equation (18)), the k 1 and n 1 will affect the time of prothrombin activator in the microchannel. From Figure 3, it can be found that when n 1 is equal to 2.0, the reagent will flow to the blood side and react with prothrombin rapidly, hence it will complete reaction when x < 5. As can be seen in Figure 4, the tendency of the thrombin content to decrease after an increase in response does not occur in the mixed zone. Therefore, it is speculated that the prothrombin activator has completely reacted with thrombin before it enters the mixed zone. The larger the k 1 value and the smaller the n 3 value, the slower the rate of thrombin formation in the flow channel. As indicated on the left of Figure 4, the increase in magnitude when k 1 = 0.1 and 0.001 is almost the same as for 1.0 and 0.1, indicating that there should be a limiting value for k 1 . As depicted on the right of Figure 4, the maximum value is obtained at x = −10 mm when n 3 = 2, which indicates that the values of k 1 and n 3 have a considerable effect on the distribution of thrombin in the microchannel. According to Equation (17), the thrombin content is proportional to the fibrin content. Thrombin is the intermediate product of prothrombin and fibrinogen. Prothrombin activator starts to produce thrombin in the blood side of the microchannel at the beginning. When the reagent and the prothrombin activator react completely, the thrombin will react only with fibrinogen, hence thrombin will decrease in the microchannel. When n 3 = 2, which indicates that thrombin will react with fibrinogen to produce fibrin, hence, the thrombin will decrease at x > 0. As depicted on the left in Figure 5, when k 2 = 0.001, the rate of thrombin generation is low, which results in the fibrin mass fraction being almost equal to 0. Conversely, part of the thrombin content continues to react to form fibrin in the mixed zone when k 2 = 0.1. When x < 0, the fibrin mass fraction increases along with the microchannel, and when it is not mixed with reagents and chemical reaction (x > 0), it starts to decrease along the microchannel.     In this study, the parameters k1-k3 obtained using the genetic algorithm were 0.987, 0.43, and 0.877, respectively, and those obtained for n1-n5 were 0.548, 1.631, 1.23, 0.459, and 1.341, respectively. These results are indicative of the behavior of the prothrombin activator and fibrin in the microchannel. As can be seen in Figure 6, the chemical reaction In this study, the parameters k 1 -k 3 obtained using the genetic algorithm were 0.987, 0.43, and 0.877, respectively, and those obtained for n 1 -n 5 were 0.548, 1.631, 1.23, 0.459, and 1.341, respectively. These results are indicative of the behavior of the prothrombin activator and fibrin in the microchannel. As can be seen in Figure 6, the chemical reaction with the anticoagulant does not occur at 1 s, because the prothrombin activator has yet to enter the mixed zone at this time. The prothrombin activator reacts with the anticoagulant within 3-5 s, and the reaction rate is high. Although a high prothrombin activator content is present in the blood at the side channel at 3 s, the activity of the prothrombin activator decreases after 5 s. Finally, the mass fraction of the prothrombin activator is nearly stable. In this study, the parameters k1-k3 obtained using the genetic algorithm were 0.987, 0.43, and 0.877, respectively, and those obtained for n1-n5 were 0.548, 1.631, 1.23, 0.459, and 1.341, respectively. These results are indicative of the behavior of the prothrombin activator and fibrin in the microchannel. As can be seen in Figure 6, the chemical reaction with the anticoagulant does not occur at 1 s, because the prothrombin activator has yet to enter the mixed zone at this time. The prothrombin activator reacts with the anticoagulant within 3-5 s, and the reaction rate is high. Although a high prothrombin activator content is present in the blood at the side channel at 3 s, the activity of the prothrombin activator decreases after 5 s. Finally, the mass fraction of the prothrombin activator is nearly stable.  (16) and (17), thrombin and fibrin in the blood undergo spontaneous reactions and a proportional relationship exists between them. The fibrin reaction has already begun before entering the mixed zone (Figure 7). A slowdown in the fibrin flow speed occurs at 3 s in the microfluidic channel because of the high viscosity of fibrin [21], which causes thrombin and fibrin to flow out of the microchannel after 5 s. In addition, the reaction extends to the reagent channel side, except in the mixed region.  (16) and (17), thrombin and fibrin in the blood undergo spontaneous reactions and a proportional relationship exists between them. The fibrin reaction has already begun before entering the mixed zone (Figure 7). A slowdown in the fibrin flow speed occurs at 3 s in the microfluidic channel because of the high viscosity of fibrin [21], which causes thrombin and fibrin to flow out of the microchannel after 5 s. In addition, the reaction extends to the reagent channel side, except in the mixed region. Figure 8 depicts the average mass fraction of each component in the microchannel as a function of time. The prothrombin activator enters the mixing zone and reacts rapidly with the anticoagulant to produce a chemical reaction at 1.2 s, and its content then quickly decreases with time. Almost no chemical reaction occurs in the microchannel. In addition, the fibrin reaction is rapid and tends to stabilize at approximately 4 s. (16) and (17), thrombin and fibrin in the blood undergo spontaneous reactions and a proportional relationship exists between them. The fibrin reaction has already begun before entering the mixed zone (Figure 7). A slowdown in the fibrin flow speed occurs at 3 s in the microfluidic channel because of the high viscosity of fibrin [21], which causes thrombin and fibrin to flow out of the microchannel after 5 s. In addition, the reaction extends to the reagent channel side, except in the mixed region.  The prothrombin activator enters the mixing zone and reacts rapidly with the anticoagulant to produce a chemical reaction at 1.2 s, and its content then quickly decreases with time. Almost no chemical reaction occurs in the microchannel. In addition, the fibrin reaction is rapid and tends to stabilize at approximately 4 s.

As indicated by Equations
The flow condition in which the RBCs and WBCs pass through the microchannel is illustrated in Figure 9. The RBCs and WBCs are ready to enter the mixing zone at approximately 4.3 s, as illustrated in Figure 8. The chemical reaction rapidly generates a large amount of fibrin, thereby decreasing the flow velocity of particles in the microchannel. In addition, the time it takes for microparticles to flow out of the micro passage increases, because the passing particles collide with the curved surface of the wall and rebound. As can be seen in the figure, the particles flow completely out of the micro flow channel after 11.7 s.  The flow condition in which the RBCs and WBCs pass through the microchannel is illustrated in Figure 9. The RBCs and WBCs are ready to enter the mixing zone at approximately 4.3 s, as illustrated in Figure 8. The chemical reaction rapidly generates a large amount of fibrin, thereby decreasing the flow velocity of particles in the microchannel. In addition, the time it takes for microparticles to flow out of the micro passage increases, because the passing particles collide with the curved surface of the wall and rebound. As can be seen in the figure, the particles flow completely out of the micro flow channel after 11.7 s. We have calculated the prothrombin time and compared the values with those found in past studies [26][27][28][29]. Table 1 shows the comparison between the point of care We have calculated the prothrombin time and compared the values with those found in past studies [26][27][28][29]. Table 1 shows the comparison between the point of care (POC) and lab PT results from the literature [30] and the simulated PT results and different blood volumes (40-100%). The simulated PT and POC PT are similar, and when the volume of the blood is low, both the POC PT and simulated PT are high. Figure 10 shows a graph of the POC PT and simulated PT results for comparison of the difference between the mean POC PT and simulated PT results. When the blood volume is more than 55%, the PT time is less than 16.2 s on average. When the volume of the blood increases the PT time is shorter. When the volume is 100%, the simulated PT time is 12.3 s and the lab POC PT result is 12 s. (POC) and lab PT results from the literature [30] and the simulated PT results and different blood volumes (40-100%). The simulated PT and POC PT are similar, and when the volume of the blood is low, both the POC PT and simulated PT are high. Figure 10 shows a graph of the POC PT and simulated PT results for comparison of the difference between the mean POC PT and simulated PT results. When the blood volume is more than 55%, the PT time is less than 16.2 s on average. When the volume of the blood increases the PT time is shorter. When the volume is 100%, the simulated PT time is 12.3 s and the lab POC PT result is 12 s.

Conclusions
In this study, a genetic algorithm method was used to estimate the various parameters in the chemical reaction to obtain the actual behaviors of blood and anticoagulant in the microchannel. The following conclusions can be drawn: 1. Through simulation, it can be found that the chemical reaction of prothrombin in blood is a factor. When the reaction rate is k1, and reaction order n2, n3 increase, prothrombin will be quickly depleted, causing the drug to flow to the blood end, but too slow the time of PT is extended. The reaction rate and reaction order can be obtained more accurately by genetic algorithm. The chemical reaction rate constants

Conclusions
In this study, a genetic algorithm method was used to estimate the various parameters in the chemical reaction to obtain the actual behaviors of blood and anticoagulant in the microchannel. The following conclusions can be drawn:

1.
Through simulation, it can be found that the chemical reaction of prothrombin in blood is a factor. When the reaction rate is k 1 , and reaction order n 2 , n 3 increase, prothrombin will be quickly depleted, causing the drug to flow to the blood end, but too slow the time of PT is extended. The reaction rate and reaction order can be obtained more accurately by genetic algorithm. The chemical reaction rate constants indicate that the reaction rate of prothrombin activator is faster than thrombin, fibrin activator is faster than thrombin, and prothrombin activator is faster than fibrin.

2.
Predicted by simulation and experimental results, high blood concentration (>65%) is more accurate in predicting PT. 3.
In order to observe whether RBCs and WBCs are obstructed in the microchannel, the PT time cannot be accurately and effectively evaluated. Therefore, the dynamic mesh method can clearly determine the time for RBCs and WBCs to pass through the microchannel and can predict the dynamic behavior of blood and coagulants. The white blood cell flow time is 11.7 s and there are no obstructions. Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.
Data Availability Statement: Data available on request due to restrictions eg privacy or ethical. The data presented in this study are available on request from the corresponding author. The data are not publicly available due to their containing information that could compromise the privacy of research participants.

Conflicts of Interest:
The authors declare no conflict of interest.

Nomenclature
Nomenclature V velocity . m ij mass from phase i to phase j X, Y mass fraction M w molecular weight P static pressure J diffusion flux R net rate of production by chemical reaction D j,m diffusion coefficient for species j in the mixture N number of chemical species in the system N r number of chemical species in reaction r v j,r stoichiometric coefficient for reactant j in reaction r v j,r stoichiometric coefficient for product j in reaction r M j symbol denoting species j k f,r forward rate constant for reaction r k b,r backward rate constant for reaction r C k,r molar concentration of each reactant and product species k in reaction r η k,r forward rate exponent for each reactant and product species k in reaction r η k,r backward rate exponent for each reactant and product species k in reaction r t time g gravity acceleration k 1 -k 3 chemical reaction rate n 1 -n 5 chemical reaction order w l weighting mass fraction f l objective function coefficient dA average mass fraction Greek letters α volume fraction β interphase momentum exchange coefficient µ viscosity τ stress tensor ρ density Subscripts i phase, species j, k, l species r reaction