## 2 Materials and Methods  
In our simulations, the blood plasma is represented as a continuous incompressible Newtonian fluid that is computed using LBM with a D3Q19 grid and the Bhatnagar-Gross-Krook (BGK) collision operator which provides second order accuracy in the bulk of the flow. The red blood cells are represented as discrete element method (DEM) membranes which are coupled to the fluid using a validated in-house immersed boundary method (IBM) [25; 38]. This IBM implementation builds on the original formalism by [PERSON] et al. [39]. The mechanical responses of the immersed cells in this work are based on the model description that is explained and validated in [25]. The RBC mesh consists of \(N_{p}=642\) vertices and \(N_{T}=1280\) triangles, whereas a platelet mesh consists of with \(N_{p}=66\) vertices and \(N_{T}=128\) triangles. These meshes are depicted in Figure 1.  
To simulate a viscosity contrast between the cytoplasmic solution of the RBC and the blood plasma, such an approach is chosen that is similar to the method described in [34]. We define the ratio between the internal and external viscosities to be \(\Lambda=\eta_{int}/\eta_{ext}\), with \(\eta_{int}\) the kinematic viscosity of the internal fluid and \(\eta_{ext}\) the kinematic viscosity of the external suspending fluid. To introduce this viscosity contrast to the fluid, one needs to keep track of the fluid nodes residing inside the membranes on the RBCs. These membranes are represented as closed surface polygonal meshes which enclose lattice nodes of the fluid grid. For these grid nodes, the relaxation parameter, \(\tau\), can be altered to express the increase of local viscosity. If a fluid node is located inside the RBC, the relaxation parameter of that node is adjusted based on the viscosity ratio between the internal and external fluid so that:  
\[\tau_{\text{int}}=\Lambda(\tau_{\text{ext}}-1/2)+1/2, \tag{1}\]  
where \(\tau_{\text{int}}\) is the relaxation time for an internal lattice node, \(\tau_{\text{ext}}\) is the relaxation time for the external lattice nodes, and \(\Lambda\) is the viscosity ratio between the external and internal fluids. Experimental work predicts the inner viscosity to be in the range of 3-10 \(\times\) 10\({}^{-3}\) Pas [2; 40]. In this work, we chose \(\Lambda=5\) for all the performed simulations with a viscosity contrast. Determining whether a lattice point is inside an RBC is done using the _even-odd rule_[41]. This rule builds on simple ray-casting where we determine if a point \(\vec{\tau}\) is inside or outside a polygon by testing how many times a ray, starting from the point \(\vec{\tau}\) to any other point outside the domain, intersects with the polygon. If the number of intersections is an even number, the point is outside of the polygon, whilst a point is inside the polygon if the number of intersections is odd. Testing whether an intersection takes place between a ray and a polygon is done by performing the Moller-Trumbore algorithm [42]. It is applied to each of the triangles which make up the RBC membrane. The chosen parameters of the simulations result in the RBC mesh consisting of 1280 triangles, enclosing approximately 720 lattice nodes. To reduce the amount of computations that have to be performed, the ray is always cast towards a fixed direction and the positions of the triangle vertices are structured using an octree data structure [43]. This drastically reduces the amount of ray-triangle intersections that have to be calculated for each lattice node. The algorithm above does not have to be performed at every iteration. The typical time step of our simulations is in the order of 10\({}^{-8}\) s, during which the membrane always moves much less than the size of a lattice node. Furthermore, the fluid nodes close to a given membrane point are already known due to the application of the immersed boundary method. The normal vectors of the surface at each of the vertices that makes up a cell are also known since the constitutive model makes use of the local curvature.