Next Article in Journal
Recent Developments and Future Directions in Flow Visualization: Experiments and Techniques
Next Article in Special Issue
Spurious Aeroacoustic Emissions in Lattice Boltzmann Simulations on Non-Uniform Grids
Previous Article in Journal
Pressure Behavior in a Linear Porous Media for Partially Miscible Displacement of Oil by Gas
Previous Article in Special Issue
Numerical Simulation of First-Order Surface Reaction in Open Cavity Using Lattice Boltzmann Method
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

The Lattice Boltzmann Method with Deformable Boundary for Colonic Flow Due to Segmental Circular Contractions

INRAE, Centre IDF Jouy-En-Josas/Antony, Université Paris-Saclay, 78352 Jouy-en-Josas, France
Fluids 2025, 10(2), 22; https://doi.org/10.3390/fluids10020022
Submission received: 2 December 2024 / Revised: 29 December 2024 / Accepted: 10 January 2025 / Published: 21 January 2025
(This article belongs to the Special Issue Lattice Boltzmann Methods: Fundamentals and Applications)

Abstract

:
We extend the 3D Lattice Boltzmann method with a deformable boundary (LBM-DB) for the computations of the full-volume colonic flow of the Newtonian fluid driven by the peristaltic segmented circular contractions which obey the three-step “intestinal law”: (i) deflation, (ii) inflation, and (iii) elastic relaxation. The key point is that the LBM-DB accurately prescribes a curved deforming surface on the regular computational grid through precise and compact Dirichlet velocity schemes, without the need to recover for an adaptive boundary mesh or surface remesh, and without constraint of fluid volume conservation. The population “refill” of “fresh” fluid nodes, including sharp corners, is reformulated with the improved reconstruction algorithms by combining bulk and advanced boundary LBM steps with a local sub-iterative collision update. The efficient parallel LBM-DB simulations in silico then extend the physical experiments performed in vitro on the Dynamic Colon Model (DCM, 2020) to highly occlusive contractile waves. The motility scenarios are modeled both in a cylindrical tube and in a new geometry of “parabolic” transverse shape, which mimics the dynamics of realistic triangular lumen aperture. We examine the role of cross-sectional shape, motility pattern, occlusion scenario, peristaltic wave speed, elasticity effect, kinematic viscosity, inlet/outlet conditions and numerical compressibility on the temporal localization of pressure and velocity oscillations, and especially the ratio of retrograde vs antegrade velocity amplitudes, in relation to the major contractile events. The developed numerical approach could contribute to a better understanding of the intestinal physiology and pathology due to a possibility of its straightforward extension to the non-Newtonian chyme rheology and anatomical geometry.

1. Introduction

“The human colon is one of the least understood organs of the body, largely because of difficulty with access, both in vivo and in terms of availability of human resection specimens for experimental studies”, according to a relatively recent consensus [1]. However, it has been realized that the hydrodynamics of the mass transfer through the digestive tract plays a crucial role in an understanding of its functioning and disease [2,3]. On the one hand, measurements and numerical predictions of complex multi-directional velocity fields and shear rates make it possible to quantify water absorption [4], to estimate the dispersion of nutrients and microbiota [5,6], to predict the bolus mixing in the gut [7,8], and to design more efficient drug release [9,10,11]; on the other hand, they are also relevant for explaining colonic dysmotility [12,13] as well as constipation and defecation anomalies [14,15,16]. At the same time, the in vivo measurements [17,18,19,20,21], the in vitro experiments [3,9,22], the predictions [23] and the modeling [7,24,25,26,27,28,29] of pressure waves inform us about muscular contractions due to the wall motility governing digestive mechanics [29,30]. In this regard, computational fluid dynamics (CFD) could provide a number of insights relating to all these physical fields due to its ability to vary characteristic parameters, to prescribe wall motility scenarios, and to visualize flow patterns, which can hardly be realized experimentally [31].
The undulatory propulsive movement of smooth muscles is commonly called peristalsis, which is considered a primary motility that spreads food through the digestive system, starting in the throat and then continuing through esophagus, stomach, small intestine, and large intestine (colon). Realistic computational models of the intestinal flow and mixing are mostly developed in the small intestine employing its geometric prototypes for both animals [24,25,32,33] and humans [3,5,26,27], and also recently using an anatomical model [29].

1.1. The Motility in Small Intestine

The small intestine is composed of the three parts called duodenum, jejunum, and ileum, and plays a principal role in water and nutrient adsorption from chyme. The peristalsis combines the high-amplitude travelling waves due to circular (inner, internal) muscle contractions with the relaxation of the longitudinal (outer, external) muscle on the wave front, e.g., [29]. Two additional motility mechanisms of the intestinal walls are (i) segmentation [34] and (ii) pendular movement [35]. The segmentation, isolated or rhythmic, is the contraction by pinching of the circular smooth muscle which breaks down chyme and brings it to the epithelial walls, making it more likely to encounter the absorptive villi. The villi, for their part, are regular structures of the mucosa folds which amplify absorption processes due to an increase in the intestinal surface area; on the hand, they also promote mixing because of their pendular movements, involving contraction and relaxation of the longitudinal muscle layer. An important dependence has been established in a human ileum: “each contraction within the clusters consisted of a larger-amplitude contraction, which was a continuation of a lower-amplitude ripple” [18].
In vitro, the “Small Intestinal Simulator” [22] imitates segmentation contractions, where pressure and velocity fields, along with mixing effects, are measured as a function of fluid viscosity. In silico, the phase–amplitude coupling mechanism [23] successfully modulates the segmentation motor pattern by waxing and waning of propulsive wave amplitude in the small intestine of a mouse. Based on these achievements, a recent one-dimensional numerical approach [36] of the nutrient absorption and bacterial growth models the radius of intestine via a superposition of two sine waves, with a high frequency for the segmentation and a low frequency for the peristalsis.
On the multi-dimensional numerical level, Fullard et al. [24,25] describe digestive flow by the Newtonian fluid, in comparison [24] with the non-Newtonian Bird–Carreau model mostly used for food and beverages. In these works, the Navier–Stokes flow equation is solved numerically by the finite-element Ansys POLYFLOW software (https://www.ansys.com/products/fluids/ansys-polyflow accessed on 29 November 2024) in a cylindrical tube subjected to the longitudinal displacement of an elastic boundary, while its pendular non-propagating longitudinal contraction is carried out according to the recorded experimental markers in the rabbit ileum [37]. Further, the longitudinal contraction is coupled [25] with the segmented circumferential contraction, which first slightly relaxes and then sharply reduces tube radius up to 75%, followed by a slower backward relaxation. In particular, it is found [25] that the segmented circular contraction alone enhances mixing in a highly viscous flow above the level of longitudinal contraction and even above that of their combination.
R. A. Agbesi and N. Chevalier [28] examine single, colinear, and counter-propagating colliding waves due to circular contractions in rat duodenum. In their numerical experiment, an axi-symmetric creeping Stokes flow is modeled with the help of the COMSOL Multiphysics software (https://www.comsol.com/, accessed on 9 January 2025); in that, the contractile waves of lumen occlusion ratio between 25 % and 88 % are induced when a linearly-elastic thin wall freely moves and deforms in a radial direction due to prescribed Gaussian-shaped displacement force density. It is found that both single and colliding waves significantly increase solute dispersion and absorption due to the formation of re-circulating vortex upstream of the contraction waves only provided that these waves are highly occluding. The pressure between two colliding waves exceeds the single wave by a factor of 10 to 100, while the jet velocity reaches 140 times the speed of the contractile wave.
N. Palmada et al. [29] extend their experimental and computational studies [3], first performed with the C-tube representation of the duodenum to the peristaltic model in an anatomical description of a human duodenum. The three-dimensional lumen flow is described by the OpenFoam software (https://www.openfoam.com/, accessed on 9 January 2025), where an elastic boundary displacement is performed according to the electrical model of the rhythmic slow waves produced by the pacemaker cells, called interstitial cells of Cajal (ICCs) [30,38]. On the one hand, their model mimics an interaction between the circular and longitudinal muscle layers; on the other hand, it takes into account their coordination by the pacemaker cells; these three types of muscle type layers are spatially described by radial subdivisions. It is found that the retrograde waves, which propagate in the opposite direction to the contractile wave, do not promote mixing; in contrast, the colliding waves (which make up about 15 % of all) produce local flow gradients sufficient for mixing enhancement, in agreement with [28]. Finally, it is estimated [29] that the multi-scale computational effort can become unrealistic to further describe the 3D pendular micro villi flow at its own scale of ≈500 μm, which would, however, make it possible to more precisely predict the flow and dispersion near the walls.
The two groups, of Brasseur [7,26,27] and of Loubens [5,32,39], apply an alternative numerical approach for solving the digestive flow and mixing, represented by the Lattice Boltzmann Method (LBM) [40,41,42,43,44,45,46,47,48]. The pioneering work [7] develops the 2D model to describe segmental contractions in an antral geometry of stomach through the boundary velocity prescribed on the moving surface, which is reproduced from the results of Magnetic Resonance Imaging (MRI). This work records the spatio-temporal structure of the contraction-induced retrograde and circular flows vs. gastric fluid pressure variation and concludes that the most occluding contractions best improve the mixing, in agreement with the conclusion [28]. However, the simulated pressure waves have much lower amplitudes than the manometric data, which is suggested [7] to happen due to the close contact of the manometer with a contracting wall.
Y. Wang et al. [26] apply the LBM to model the 2D pendular coaxial micro-flow, coupled with the macro-scale lumen eddy flow driven by a lid, in regular arrays of villi. More recently, Brasseur and Wang [27] extend villi movement to the 3D lateral oscillations. These two works examine the role of villi arrays and conclude that (i) the pendular motion is more significant for adsorption enhancement than an increase in an epithelial surface area due to villi’s networks, and (ii) the multi-scale flow interaction further increases adsorption due to creation of the 3D micro-mixing.
In parallel works, Loubens et al. [32] employ the LBM to simulate the 2D microscopic creeping fluid flow through villi’s arrays in the proximal duodenum of rats and guinea pigs. They model the oscillations of villi’s surface according to either a cyclic variation of the longitudinal strain-rate or a recorded spatio-temporal evolution of the longitudinal contractions; it is concluded that the pendular motion induces a sufficiently high local shear to promote mixing in a sufficiently viscous fluid. A further extension [5] to the 3D Newtonian and non-Newtonian (pseudo-plastic) fluid flow, which is also coupled with the advection–diffusion equation (ADE) for nutrient mixing and absorption, observes enhancement in the transport of solid particles due to the dynamics of the villi.
Y. Zhang et al. [4] simulate the 2D multi-scale channel flow on the moving mesh; their model combines ideas [27,32] and employs the COMSOL MultiPhysics software; the results indicate a 72 % increase in the mass transfer and adsorption in the rat duodenum due to pendular motion of the villi. Zha et al. [8] develop the 2D axi-symmetric multi-scale parametric study with the help of the Ansys Fluent Fluid Simulation software (https://simtec-europe.com/?gad_source=1&gclid=Cj0KCQiAqL28BhCrARIsACYJvkfSo_d0r96A3Hg0XkNAY8mguYks472mS7nM9Hvu_zkUAav4gRhBrZEaAtj9EALw_wcB, accessed on 9 January 2025). Their work concerns the role of amplitude, frequency, and wavelength of a wall contraction on the flow and mixing in a human duodenum, when an inner layer of circular submucosa folds, approximately 8 mm in height, is covered by villi 1 mm in height. The results show that segmentation due to the contraction of circular folders can intensify mixing.

1.2. The Colonic Motility

According to the textbook by Christensen [49], the human large intestine is “surprisingly more complex than that of other omnivores”, and being about 1.5 m long, it consists of five components, “appendix, cecum, colon, rectum, and anal canal”. In turn, the colon is subdivided into the proximal colon, partitioned into ascending and transverse colons, and the distal colon, partitioned into descending and sigmoid colons (see Figure 1). Like the rabbit’s colon, an external (outer) longitudinal smooth muscle is separated into three taeniae coli (taenia) [50,51], which follow the entire colon and then merge to form a single longitudinal muscle layer that continues into the rectum. The three taenia interact, via the interstitial cells of Cajal (ICC) located in intramuscalar layer, with the internal (inner) circular smooth muscle layer [13,52,53]. The rings of circular contractions subdivide the bulging intertaenial walls into a chain of segments approximately 2 cm long, called haustra. During contractions, the haustral membrane occludes colonic cross-section and gives it a triangular shape; the capsule endoscopy images [53] visualize spectacularly the opening of the lumen and even its complete occlusion.
The seminal work [52] identifies the four types of motility in the rabbit colon using the high-definition spatio-temporal mapping of muscle membrane movement. These four types of motility are composed of two groups, first the two components of mass (neural) peristalsis, which are (i) “deep circular contractions” and (ii) regular “longitudinal taenian contractile activity”, and second the (iii) phasic contractions of the circular muscles called “ripples” and (iv) “fast phasic contractions” due to the interaction between the longitudinal and circular muscles. Following the methodology [52], the phase coupling between the motor patterns of different frequency and amplitude was quantified [55] from the relationships that exist between changes in diameter (video image) and corresponding changes in intraluminal pressure (manometry). In turn, phase-coupling modulation [51] extends the methodology [23] to the rabbit colon. This work shows that slowly spreading circular contractions subdivide the colon into the pockets (haustra) whose rhythmic membrane contractions increase and decrease amplitude at the frequency of ripples; the latter being identified as low-amplitude, rapidly propagating contractions, in agreement with [52].
Consensus [1] classifies all reported colon motor functions and delineates the similarities between those observed in animals and humans in vivo, but also in in vitro, then concludes that “the most obvious similarity is that between neural peristalsis that propels the contents in animals when distended and the high-amplitude propagation contractions (HAPCs) in humans”.
Over the years, the HAPC pattern has been defined as pressure waves mostly originating in the proximal colon, with a high amplitude varying between 50 and about 150 mmHg, and extending during 10–30 s at least 10–30 cm along the length of the colon [51]. Over the past decade, high-resolution colonic manometry (HRCM) has made available pressure wave measurements, which are typically recorded by 1 cm distanced sensors along a stripe 30–60 cm long [17,19,21,53,56,57,58]. It was observed that the high-frequency cyclic motor model coexists with the low-frequency one [57], and that a quasi-permanent low-frequency retrograde cyclic motor pattern generates a braking mechanism to transit through the sigmoid colon [16,19,57,58]. In turn, the cyclic myogenic motor pattern (CYMP), due to ripple activity, induces mixed pressure waves, which propagate in both antegrade (oral) and retrograde (aboral) directions, while the colliding waves annihilate each other [53].
In addition, there is a chaotic interhaustral flow showing different velocity patterns in the three intertaenial regions [51]. This vortex flow should explain a surprisingly long residence time of microbiota in the colon. In silico, simplistic models of gastric flow account for this effect via a huge increase in the molecular diffusion coefficient [59], a creation of an artificial vortex [6], or even an artificial viscosity variation inside the mucus boundary layer [60].
However, a relatively small number of numerical works undertakes a direct discretization of the hydrodynamic equations in the colon. So far, an early work [61] employs the COMSOL Multiphysics software to model the Newtonian flow through a lumen represented by a cylindrical tube. In this work, the localized reductions in tube diameter are induced by applying the phasic Gaussian pressure pulse. Remarkably, the measured and simulated pressure profiles [61] show qualitatively similar responses, with pressure increasing rapidly from the onset of wall deformation; in addition, the largest peaks in pressure magnitude then correspond to a higher deformation rate and/or higher fluid viscosity. The authors conclude [61] that these results disapprove of a popular opinion that the pressure fluid manometry receives a weaker signal when the lumen is not obstructed.
Recent work [62] depicts, based on anatomical data, a long colonic slice measuring 60 cm in length through a three-dimensional curvilinear canal, where the mean radius decreases linearly from the cecum to the rectum, while a sine-wave local variation of the radius imitates a peristaltic wave subjected to parametric control over the occlusion range. In this work, the Navier–Stokes equation integrates the non-Newtonian power law rheology of chyme flow, while the model is simulated by the Ansys FEM software with dynamic mesh reconstruction. The results [62] report that the mixed flow of two similar amplitudes reaching a velocity interval of [ 0.2 ,   0.33 ] cm / s occurs mainly in an ascending colon, while in the narrow descending colon, the propulsive antegrade flow prevails, and its amplitude increases with the peristaltic wave speed. However, these numerical results do not seem to detect the retrograde cyclic motor pattern observed experimentally, which is considered responsible for the braking mechanism in the descending colon [19,58].

1.3. The Dynamic Colon Model (DCM)

The “Dynamic Colon Model” (DCM) [9,63,64] in vitro provides functionality of the ascending colon from the Caecum to the Hypatic flexure; the DCM device is composed of 10 haustral units in the horizontal cylindrical tube of 20 cm long (see Figure 2). The DCM is based on clinical data obtained in vivo in the proximal colon from Magnetic Resonance Imaging (MRI); this means that the DCM reproduces an anatomically representative membrane motility and lumen occlusion scenarios due to peristaltic circular contractions, but neglects the longitudinal ones. A synchronized segmented motility pattern of the DCM is composed of four stages, such as (i) deflation, where a neutral value of lumen occlusion decreases to the minimum (relaxed) degree, (ii) inflation, where the wall contracts and occlusion value increases to its maximum degree, (iii) pause, and (iv) backward relaxation to the neutral state with a speed four to five times slower than the contractile speed, imitating an effect of membrane viscoelasticity. The motility pattern obeys an “intestine law”, meaning that when a given segment reaches its most occluded position, its downstream neighbor becomes most relaxed. A high intraluminal pressure peak of about 25 mmHg is registered [9] with a 40 % occlusion increase at peristaltic wave speed of v w = 2 cm / s ( V M A X hereafter), which is close to maximum wave speed value observed in vivo. Subsequent work [64] performs the two types of MRI analysis following [31], where it addresses wall motion, luminal flow pattern, and stream-wise velocity of the contents; these data are recorded in several segment-central tagged locations at different motility stages by increasing occlusion to 60 ± 5 % but reducing v w to 0.4 cm / s ( S L O W ) and 0.8 cm / s ( F A S T ). The degree of mixing is estimated [64] from the predominant direction of the tag displacement (antegrade or retrograde); interestingly, this indicator is practically independent of the low/high volume level of the liquid in lumen (hereafter, L V O L 60 % and H V O L 80 % , respectively), but it decreases from low-viscosity fluid ( L O V I S ) to high-viscosity ( H I V I S ) and increases with wave speed v w .

1.4. The DCM “Digital Twin” (DCMDT)

Unlike the direct discretization methods (DDMs hereafter), the Smooth Particle Hydrodynamics ( S P H ) is a mesh-free numerical method which approximates solutions of the Navier–Stokes equation via a simplified molecular dynamics of the fluid and solid particles. S P H is adapted [65] for full ( 100 % filled) simulations in lumen using the “Mass-Spring Model” for solid particle dynamics, which aims to reproduce the peristaltic contractions of the membrane.
“The Discrete Multiphysics Method” [66] then extends the S P H model to a partially filled colon, also extending the mechanics of membrane deformation through stretching and bending; this method is validated by comparing the tracer tracking with the DCM [63]. Further, the hydrodynamic effects on mixing are compared through simulations [67] of single-phase (full filled), two-phase (gas-liquid), and free-interface (partially filled) fluid. The highest amplitudes in velocity and mixing are then observed in the single-phase fluid, while no significant impact is observed between the two-phase and free-interface models, concluding that a detailed description of the gas phase is not necessary.
Recently, the DCM “digital twin” ( D C M D T ) [11] replicated the segmented motility pattern applying the Hookeen force to solid particles, which move according to the “Lattice-Spring Model” and respect the prescribed visco-elasticity rate. The D C M D T is validated [68] through the Magnetic Resonance Imaging [64] of the DCM results, reporting good agreement for the transient distributions of the mean and peak longitudinal velocity in the central and bounding segments. These reported measurements [68] are called MRI hereafter; they address the S L O W and F A S T wave speeds, the L O V I S and H I V I S kinematic viscosities, all combined with the L V O L and H V O L filling degrees. Unfortunately, no experimental (MRI) or simulated ( D C M D T ) pressure distribution appears to be reported in association with velocity measurements.

1.5. Summary on the Numerical Aspects

The above analysis of small intestine modeling suggests that the flow simulation software which is based on the finite-elements and finite-volumes methods is limited in the multi-scale to an axisymmetric representation of the intestine because of resorting to the specific meshing tools on the moving surface, where the computational effort does not yet allow the coupling of the 3D anatomic discretization of lumen with the mobile layer of micro-structures.
Concerning the colonic modeling discussed above, and according to our definition of the degree of occlusion, the DCM motility patterns appear not to conserve lumen volume, which might become possible due to the presence of an air layer in the device (see Figure 3), and above all a possibility of fluid mass to freely rise upward the “hepatic flexure” (see Figure 2). In turn, the D C M D T [68] allows fluid “to escape” from the modeled colon strip. Furthermore, the S P H method does not seem to strictly interconnect the fluid mass with its volume, unlike the direct discretization methods for incompressible flows ( D D M ), which generally constrain mesh deformation to the volume preservation. This issue is expected to pose difficulties for the D D M to cope with a temporary accumulation of chyme content or a volume variation of the modeled colonic stripe. Another D D M difficulty could be caused by mesh generation in the case of only piece-wise continuous bounding contour due to haustra segmentation.
The LBM might be regarded as an “intermediate” approach between the S P H and D D M ; on the one hand, unlike the S P H , the LBM operates an evolution of the real numbers (populations) on the regular grid according to their prescribed discrete velocities; on the other hand, unlike the D D M , the macroscopic equations represent the second-order accurate approximation of the exact microscopic conservation laws satisfied by the populations in time and space. Indeed, alternatively with the D D M , the LBM is explicit in time, and therefore does not require any linear or non-linear solver; moreover, the LBM’s main bulk (collision) operations are local. Another key point is that this method proceeds both static and moving solid surfaces staying on the Cartesian grid; the LBM is then expected to support the extreme occlusion scenarios more robustly than the D D M due to no need for remeshing. All these factors make the LBM particularly suitable for multi-processor computations, e.g., [69], as well as for an efficient implicit coupling between the macro and micro flows, as recognized for simulations in heterogeneous porous media, e.g., [70].
Furthermore, the idea of LBM solvers operating directly on voxel images of a given porous structure was judiciously extended towards the LBM resolution of different biological mechanisms on the regular grids formed from pixels of medical images. This approach is sequentially developed in works [71,72,73] for a segmented description of thrombosis growth in a cavity of a giant cerebral aneurysm, which is modeled by LBM coupling [74] of anisotropic diffusion and blood-flow solvers. Indeed, several alternative LBM approaches for hemodynamics are under intensive development [75,76,77].
On the other hand, as pulsatility presents one of the key concerns that distinguishes physiological flows from other branches of fluid mechanics, such as, for example, flows in the arteries due to the beating of the heart, or a flow of air in the lung when breathing, the LBM’s emphasis in this work on colonic flow should extend to other pulsatile mechanisms.

1.6. Towards the LBM Colonic Modeling

In this work, we first formulate and validate the 3D LBM algorithms with deformable boundary ( L B M - D B ) for solving the Navier–Stokes equation of the Newtonian fluid. These schemes extend the moving boundary algorithms [78] called there “without LBM in the solid”. We then apply the L B M - D B to model lumen flow driven by the D C M -type segmented colonic motility due to propagating circular contractions. In our approach, an occlusion protocol makes it possible to model the deformation of the membrane surface according to the Dirichlet velocity condition as a function of time and space. In this, the L B M - D B prescribes the boundary closure in the “fluid” nodes and reconstructs the population solution in the “fresh” nodes, emerging from a deformable solid surface, using the most accurate Dirichlet boundary schemes [44,79]; their appropriate combinations are pre-selected through a specially designed reference simulation in the deforming channel flow in the presence of acute angles on the moving front.
The L B M - D B does not sacrifice the compactness of the previously employed moving boundary schemes for gastric flow due to the avoidance of interpolations for the unknown population and macroscopic fields, but also due to keeping away from an adaptive grid refinement in the vicinity of the wall. The latter might become required using the low-accurate bounce-back reflection, e.g., [5], where the discretization error due to the stair-case surface approximation only decreases linearly with the spatial resolution, even on static walls [80,81]. On the back side, the advanced boundary rules do not preserve the mass of outgoing populations [78,82], which is consistent with their analytical (parabolic) exactness in grid-inclined channel flow [80], while an artificial ensuring of mass conservation is known to deteriorate boundary accuracy on static and deforming walls [83]. At the same time, only linear accuracy is expected from an alternative approach of the “Immersed Boundary Method” ( I B M ) [84,85] due to its approximate localization of solid surface and neglect of the flow-curvature and local pressure gradients [86,87]. The I B M -based time-explicit fluid–structure interaction, e.g., [39,88,89], however, rivals in popularity the advanced direct Dirichlet LBM schemes on the moving fronts due to its robustness and solution smoothness [90,91].
The L B M - D B is expected to enable much more freedom in terms of the deformation and volume variation than the D D M because the proposed algorithms do not require conservation of the fluid volume by their construction. Hence, the L B M - D B supports any deformation prescribed for the surrounding walls and, in particular, any motility protocol imposed on the bounding colonic membrane. However, according to the L B M - D B reconstruction, the population mass in the “fresh” nodes does not necessarily scales exactly to the change in volume with a reference value of the fluid density. A “physical” reason for this is that the L B M - D B produces a weakly compressible fluid, where the Mach number is controlled by an adjustable numerical parameter, known as the “speed of sound in square” c s 2 , which linearly relates the local density (population mass) to modeled pressure.
Therefore, the reconstruction algorithm should produce a very smooth transition of velocity and pressure (density) towards the moving wall. For this objective, first, the L B M - D B pre-selects the Dirichlet boundary schemes that respect the non-equilibrium pressure–gradient component in population solution. They were originally represented by the two-node parabolic-accurate “multi-reflexion” ( M R ) class [48,78] and the “magic linear” M L I class [44,81,92], but these two sets of infinite members have recently been extended to a single node, linear but pressure–gradient accurate, L I 4 sub-class [79]. Second, following [78], the L B M - D B algorithms keep the definition of the local density consistent in the “fresh” nodes. To quantify the conservation and compressibility effects, a temporal evolution of the colonic pressure and longitudinal velocity is reported in parallel with the corresponding dynamics of lumen volume and fluid mass.
L B M - D B numerical studies are conducted in the segmented colonic strip whose membrane’s motion imitates the motility patterns of the DCM [9] and MRI [68] experiments. These two motility models are combined with the five scenarios of increasing occlusion degree, and hence increasing numerical difficulty, moving far away from their in vitro counterparts. Each of these protocols is conducted in a full-volume lumen flow with the two viscosity values [68] ( L O V I S and H I V I S ), and three peristaltic speeds: S L O W & F A S T [68] and VMAX [9]. The inlet/outlet configurations correspond to either periodic or closed conditions, which can be seen as an attempt to mimic either an intermediate fragment of the ascending colon or a permanent reflux in the sigmoid colon. Additionally, all these configurations are applied with the two distinct transverse shapes of the colonic tube: the traditional shape is circular (“C”), while an alternative new opening is called “parabolic” (“P”). The idea of the “P” cross-section is that the three parabolic contours meet each other at the fixed vertices of an inscribed triangle, which mimics a triangular-like lumen aperture formed by the three taenia bisections with the contracting circular muscle [53].
In what follows, Section 2.1 analyses the DCM experiments [9,68]; Section 2.2 describes and validates L B M - D B ; Section 3 introduces and performs the L B M - D B colonic experiments; Section 4 summarizes the findings and outlines perspectives. Appendix A, Appendix B and Appendix C complement Section 2.1, Section 2.2 and Section 3, respectively.

2. Materials and Methods

This section discusses the DCM experiments [9,64,68] and introduces the Lattice Boltzmann Method with the Deformable Boundary ( L B M - D B ).

2.1. The DCM, MRI and LBM-DB: Geometry, Luminal Occlusion, Motility and Flow

This section provides our understanding of the DCM experiments [9,64], presents their L B M - D B counterparts, and analyses the MRI velocity profiles [68].

2.1.1. Geometry

The DCM device of the total volume V = 290 cm 3 is segmented into ten haustral units of 2 cm in length, and its 3D transverse shape is presented in Figure S1 [9]. The 2D cross-section of the DCM colonic tube forms a “triangular” opening of the lateral length l Δ   =   3.8 cm , which is delimited by the three symmetric internal parabolic-type contours inside the lumen (see Figure 2 and Figure 3). These three parabolic contours are separated by thin fixed circumferential segments, 0.5 cm thick, which imitate the bisection of the three taenia colis (longitudinal muscles) and circular smooth muscle. The haustral units are separated by the “semilunar folders” of 0.2 cm thickness.
Numerical models D C M D T and L B M - D B subdivide the cylindrical tube into ten haustra segments of length l = 2 cm . At the neutral degree of occlusion, D C M D T represents the contracting contour by the circular ring composed of the 25 solid particles; among them, the three symmetrically located particles mimic the bisection with the longitudinal taenia colis and stay immobile, while the other solid particles relax and constrain according to the applied radial forces.
L B M - D B starts to construct the 2D cross-section from an equidistant triangle of lateral length l Δ   =   3.8 cm ; it is inscribed in the circle of radius r = l Δ × 3 3 2.19 cm representing a delimiting boundary of the cylindrical tube, having for volume V = π r 2 × l × 10 302.43 cm 3 . The L B M - D B then builds an inner “parabolic” membrane contour where the lumen aperture is delimited by the three deforming parabolic curves, sharing the immobile vertices of an inscribed triangle; in this work, the parabolic curves are set non-overlapping and symmetric. Alternatively, an interior circular membrane contour deforms uniformly along its radius without any fixed point. The two cross-sections are dynamically built from the prescribed (current) value of lumen occlusion degree (see Appendix A). An occluded part of the colon is located between its external (outer) wall and the internal deforming surface, delimited by either a “parabolic” (“P”) or a circular (“C”) contour.

2.1.2. Luminal Occlusion

An opening of the colon is dynamically occluded by the deforming membrane whose motion mimics the circular smooth muscle contractions between the smallest and the highest occlusion degrees called, respectively, relax and forward ; the contractions start from an intermediate ( neutral ) occlusion value. As one example, Figure S1b [9] specifies the distance of l r , n = 0.6 cm between the relax and neutral surfaces, and the distance of l r , f = 1 cm between the relax and forward surfaces; in this case, all three contours are convex and their openings have a lower degree of occlusion than that of the inscribed triangle. Additionally, both convex and concave membrane contours are depicted in Figure S4b [9], increasing “the % reduction of the cross-section of the DCM unit” from 11 % to 94 % . Figure 5 [9] presents the modeled relative (with respect to neutral occlusion) contractile interval from 0.4 ( relax , convex) to 0.4 ( forward , concave), whereas it is reported to vary between 0.2 and 0.6 with the MRI and D C M D T experiments [68]; however, it seems to us that these works do not provide any formal definition of “occlusion degree” and do not give an exact value of neutral occlusion applied in their experiments.
In the L B M - D B colonic model, an occlusion degree o ( t ) expresses the ratio of the 2D occluded transverse surface to its whole circular surface value π r 2 ; an open (lumen) fraction of the surface, available for flow, is then defined as s ( t ) = 1 - o ( t ) . Hence, an occlusion scenario is specified by the three occlusion degrees { o r , o n , o f } which define the three principal membrane positions { relax , neutral , forward } .

2.1.3. Motility

The whole membrane contracts according to the motility protocol which prescribes (i) motility pattern for each segment and (ii) how the motility pattern expands over the haustral segments. Starting and ending in the neutral position, the DCM motility pattern can be composed from the five subsequent time steps:
1.
t r : relax from o n to o r ;
2.
t r , p : pause r from o r to o r ;
3.
t f : forward from o r to o f ;
4.
t f , p : pause f from o f to o f ;
5.
t b : back from o f to o n .
In order to mimics the viscoelasticity of the intestinal wall in vivo, the DCM performs the back relaxation with a slower speed than when it contracts. The motility pattern aims to reproduce what is called the “intestine law”. This law constrains a downstream neighbor segment n + 1 of a given segment n to reach the most relaxed position o n + 1 ( t ) = o r when segment n arrives at its most occluded position, o n ( t ) = o f ; the motility pattern moves downstream the colon with the prescribed speed v w of the “propagating pressure wave” ( P P W ). The S L O W and F A S T regimes [68] apply v w = 0.4 cm / s and v w = 0.8 cm / s , respectively, while the DCM [9] operates with v w = 2 cm / s , hereafter referred to as V M A X ; this value reflects the highest realistic travel P P W speed amplitude.
In the L B M - D B simulations, all segments are initially found in the neutral position, o ( t = 0 ) = o n . The first segment n = 1 starts from the forward step following the DCM [9], where it adapts its occlusion speed to reach the highest occlusion degree o n = 1 ( t = t r + t r , p ) = o f , exactly at the moment when the downstream neighbor n = 2 leaves its relax position, with o n = 1 = o r when t [ t r , t r + t r , p ] . From the second segment and downward, all segments perform the five steps from t r to t b obeying the “intestinal law” during their individual motility cycle t [ 0 , t w ] , such that segment n + 1 starts its relaxation step t = t f later than segment n. Hence, when segment n arrives to its forward position at its relative motility time t ( n ) = t r + t r , p + t f , segment n + 1 ends its relaxation at relative time t ( n + 1 ) = t r + t r , p ; the motility cycle over N s haustral units then lasts T N s = t w + ( N s 2 ) t f .
L B M - D B quantifies the viscoelastic effect through ratio e b , f = t b t f . D C M D T then aims to match e b , f = 3 with the help of the radial forces, whereas the DCM [9] directly applies e b , f 4.57 (if we set e b , f equal to the ratio of their forward occlusion speed 1.6 cm / s to their back relaxation speed 0.35 cm/s).
L B M - D B applies and compares the two motility protocols from Table 1, M O T - I and M O T - I I , which advance with the three propagating speeds v w . The characteristic propagation time t f = l / v w is defined with haustra length l = 2 cm . M O T - I applies the three contraction steps without pauses ( t r , p = t f , p = 0 ) with t r = t f , t b = 3 t f following [68]; segment n then starts its back relaxation immediately after reaching the highest occlusion degree, when the downstream segment n + 1 starts its contraction. In Figure 4 left, M O T - I S L O W offers the relative to neutral interval [68] [ o r , n , o f , n ] = [ 1 5 , 3 5 ] .
M O T - I I sets t r = l r , n / v r = 0.375 s with the VMAX, providing (a) distance l r , n = 0.6 cm between neutral and relax positions from Figure S1 [9], and (b) v r = 1.6 cm / s in the V M A X regime [9]; M O T - I I then relaxes faster than M O T - I . According to [9], M O T - I I then sets t f = 1 s and t b = 4.625 , and hence e b , f = 4.625 (to simplify data, this value is set slightly larger than its DCM estimate above); the fastest M O T - I I V M A X cycle then lasts t w = 7 s in agreement with Figure 5 [9]. Since M O T - I I makes break t f , p = t f during occlusion peak o ( t ) = o f , segment n then starts its back relaxation exactly at the moment when its downstream neighbor n + 1 reaches occlusion peak o ( t ) = o f . M O T - I I V M A X is exemplified in Figure 4 right with [ o r , n , o f , n ] = [ 3 5 , 2 5 ] .
These two motility patterns do not preserve an open volume except of the particular combinations of motility times and occlusion degree; the corresponding relative open volume variation δ V ( t ) = V ( t ) V ( o = o n ) 1 with respect to the neutral position is additionally displayed in Figure 4.

2.1.4. Fluid

The DCM( M R I ) and D C M D T experiments [68] are conducted with the low fluid volume of 60% ( L V O L ) and the high fluid volume of 80% ( H V O L ). The DCM fluids [9] are represented by the aqueous NaCMC solutions of two viscosity values, both obeying the power law rheological behavior in response to shear stress. D C M D T simulates the Newtonian fluid where the two viscosity values, H I V I S and L O V I S , are established according to the “linear viscoelastic regime” with a shear rate below 40 [s−1]:
1.
H I V I S : μ = 85 m P a s , ρ = 1020 k g / m 3 , ν H = 8.3 ( 3 ) × 10 5 m 2 / s .
2.
L O V I S : μ = 26 m P a s , ρ = 1017 k g / m 3 , ν L = 2.45 × 10 5 m 2 / s .
L B M - D B simulates two kinematic viscosities ν p h = μ / ρ = { ν H , ν L } with ratio ν H ν L 3.4 .

2.1.5. Entry/Exit

According to the DCM description [66], the fluid rises and moves through the “hepatic flexure”, located at the end of the device, and according to [68], “the hepatic flexure is modeled as a reduction to create a back pressure, guided by the in vivo situation”. In turn, D C M D T operates with a relatively long drive-tank of about 40 cm; it is located at the end of the segmented colon and allows the fluid to escape it.
L B M - D B considers two bounding configurations: closed and periodic. In the first case, the segmented colon is closed with zero (no-slip) velocity on the tube’s entry and exit, as, for example, imitating a permanent mixed wave propagation in a sigmoid column [19,58]. In the second case, the two ends are inter-connected by periodic conditions with the aim of modeling an intermediate colonic fragment, for example, in an ascending or transverse colon, without any interference with the entry and exit conditions. Short neutral static buffers of 0.4–0.6 cm length are inserted before the first and after the last haustral unit.

2.1.6. The MRI Velocity Measurements in the DCM

Table 2 reports the mean velocity interval [ v ¯ m i n , v ¯ m a x ] and the peak velocity interval [ v m i n , v m a x ] which we estimate in Segment 6 ( s . 6 hereafter) over an entire device motility cycle according to the “ M R I ” and “ M R I ( p e a k ) ” graphical output from Figures 7 and 9 [68] ( S L O W ) and Figures 10 and 11 [68] ( F A S T ).
Mean velocity distribution v ¯ ( t ) is obtained via an averaging across the tagged image midway the cross-section of a given segment, while peak distribution v ( t ) is averaged over several central voxels. The velocity profiles are displayed [68] together with the relative occlusion distribution o ( t ) o n o n [ 0.2 , 0.6 ] [recall, we do not know their formal definition o ( t ) and their value o n ].
All reported velocity profiles in s . 6 , s . 2 , and s . 10 demonstrate one distinctive negative minimum value v ¯ m i n = min t ( v ¯ ( t ) ) and v m i n = min t ( v ( t ) ) , occurring at either the most constrained position of the given segment or quickly after it, when the downstream segment only starts its contraction. Table 2 provides the (relative) moment value t m i n [ t r + t f , t w ] of the largest retrograde amplitude v ¯ m i n . Table 3 shows that normalized values t m i n ( t f + t r ) t b [ 0 , 0.4 ] , confirming that the highest retrograde event happens at either maximum contraction or the first half of back relaxation.
Table 2 shows that amplitudes | v ¯ m i n | , and especially | v m i n | , exceed their antegrade (stream-wise, positive) counterparts, v ¯ m a x = max t ( v ¯ ( t ) ) and v m a x = max t ( v ( t ) ) , respectively. Typically, the antegrade maximum amplitudes occur first during the upstream contraction phase, and then during extreme relaxation in the last segment. Table 3 normalizes the mean-and peak velocity intervals by v w ; these data suggest that v ¯ m i n and v m i n scale almost linearly with v w , while their dependence on kinematic viscosity diminishes (i) with wave speed and (ii) fluid volume, in that | v ¯ m i n | increases very slightly from L V O L to H V O L with S L O W , and this trend reverses with F A S T . v ¯ m a x also scales approximately with v w and shows a very weak dependence on viscosity. These estimates quantify observation [68] that “faster wave produces higher velocity”; moreover, they suggest an almost a linear scaling of the highest mean velocity amplitudes with wave speed.
Additionally, Table 3 presents the ratios of the retrograde to antegrade velocity amplitude, R ( v ¯ ) = | v ¯ m i n | v ¯ m a x and R ( v ) = | v m i n | v m a x ; they confirm that the retrograde amplitude dominates the antegrade one, with the largest contrast R ( v ¯ ) = | v ¯ m i n | v ¯ m a x 5 and R ( v ) = | v m i n | v m a x 7.5 in the L V O L H I V I S fluid; H V O L displays smaller disparity, with R ( v ¯ ) = | v ¯ m i n | v ¯ m a x [ 1.675 , 2.5 ] and R ( v ) = | v m i n | v m a x [ 3 , 3.75 ] from L O V I S to H I V I S .
Finally, Reynolds number Re ( | v ¯ m i n | ) = | v ¯ m i n | l ν p h is computed with the highest retrograde mean-velocity amplitude | v ¯ m i n | and haustra length l = 2 cm. An estimated range is Re ( | v ¯ m i n | ) [ 0.327 , 6.53 ] , with the largest value in the case L V O L H I V I S F A S T .

2.2. The LBM with the Deformable Boundary (LBM-DB)

Originally, the LBM operates the Newtonian fluid inside a fixed physical volume bounded by a non-deforming solid surface. Several moving boundary algorithms, applied within the framework of the Parallel Particles code (ParPac [93]), have been described and validated [78]; we extend here one of them, called there “without lbm in solid”, and denote all novel schemes by L B M - D B . The L B M - D B allows solid bodies to move and the fluid–solid interface to deform, while the fluid domain is not strictly limited to the conservation of mass and volume on the solution of weakly compressible flow to the Navier–Stokes equation. L B M - D B applies the “collide and stream” steps only inside a fluid sub-domain, while the Dirichlet boundary conditions and, based on them, the reconstruction step define unknown populations that would enter the fluid domain from the static and moving boundary, respectively; these last two steps are improved in the present algorithms according to more advanced schemes [44,79,81], with emphasis on irregular geometries (corners).

2.2.1. L B M - D B : Main Steps

L B M - D B applies on the d- dimensional equidistant computational grid composed of fluid nodes r V f ( t ) ; the nodes are interconnected by the discrete anti-symmetric velocity set q = Q m / 2 Q m / 2 c q , while an immobile velocity vector is c 0 = 0 ; the Q m = Q 1 non-zero velocities have either zero or ± 1 components along the coordinate axes; convention c q = c q , q = 1 , , Q m / 2 is adopted for every pair of the opposite velocities. The minimal suitable hydrodynamic velocity sets are d 2 Q 9 ( Q m = 8 ) in 2D and d 3 Q 15 ( Q m = 14 ), d 3 Q 19 ( Q m = 18 ) in 3D. We operate with d 3 Q 15 , but L B M - D B is transparent for spatial dimension and discrete velocity sets, in principle at least.
A time step of the L B M - D B iterative update includes the following main steps:
(a)
c - s - s t e p . The local variables of the method are real numbers q = Q m / 2 Q m / 2 f q ( r , t ) called populations; their post-collision values q = Q m / 2 Q m / 2 f ^ q ( r , t ) are computed via a local linear collision transformation, then they are streamed to the directional neighbors: f q ( r + c q , t + 1 ) = f ^ q ( r , t ) .
(b)
bc - s t e p . We let S ( t ) denote a (virtual) fluid–solid interface at time t; S ( t ) can be composed of static and moving, or deforming, components. A fluid node  r b V b ( t ) V f ( t ) is called a boundary node if it has at least one cut link c q such that r b + δ q ( r b , t ) c q S ( t ) with δ q ( r b , t ) [ 0 , 1 ] .
bc - s t e p defines with the individually prescribed boundary rule for a cut link c q population solution f q ( r b , t + 1 ) .
(c)
f - s t e p . When c - s - bc - s t e p is performed on fluid set V f ( t ) , the fluid-solid interface moves towards S ( t + 1 ) , where each solid node transformed into fluid boundary node is called a “fresh” node r n ( t + 1 ) R ( t + 1 ) .
(re)fill  f - s t e p defines population solution q = Q m / 2 Q m / 2 f q ( r n , t + 1 ) on all newborn nodes.
(d)
r - s t e p . f - s t e p is the main component of reconstruction  r - s t e p , which in addition can also include several preparatory operations, such as population and/or macroscopic field interpolations; and r - s t e p may iterate f - s t e p adapting ideas [94,95,96]. L B M - D B tries to avoid population interpolations, having a preference to the use of bc - s t e p as much as possible in the “fresh” nodes.
In summary, L B M - D B iteratively performs the following update:
1.
Collision-Stream Boundary step to define q = Q m / 2 Q m / 2 f q ( r , t + 1 ) for r V f ( t ) ;
2.
Advancement from S ( t ) to S ( t + 1 ) ;
3.
Identification of “fresh” nodes r n ( t + 1 ) R ( t + 1 ) ;
4.
Reconstruction to define q = Q m / 2 Q m / 2 f q ( r n , t + 1 ) for r n ( t + 1 ) R ( t + 1 ) ;
5.
Identification of fluid list V f ( t + 1 ) ;
6.
Identification of boundary list V b ( t + 1 ) ;
7.
Return to Step 1.
Details on code organization and parallelization can be found in Appendix B.1

2.2.2. The c-s-Step

This step computes post-collision solution q = Q m / 2 Q m / 2 f ^ q ( r , t ) in Equation (1a) and shifts it to the directional fluid neighbor in Equation (1b):
collide : f ^ q ( r , t ) = f q ( r , t ) + n ^ q + ( r , t ) + n ^ q ( r , t ) , r V f ( t ) .
stream : f q ( r + c q , t + 1 ) = f ^ q ( r , t ) , i f r + c q V f ( t ) .
Hereafter, n ^ q + = n ^ q + , q = 0 , 1 , , Q m / 2 , denotes the post-collision symmetric component, while n ^ q = n ^ q , q = 1 , , Q m / 2 , is its anti-symmetric counterpart. The simplest consistent linear collision (1a) for modeling an incompressible or weakly compressible Navier–Stokes equation ( N S E ) is the “Two Relaxation Times” ( T R T ) operator [44,48,70] which computes post-collision with two locally prescribed rates τ ± ( r , t ) :
n ^ q ± = 1 τ ± ( f q ± e q ± ) , f q ± = 1 2 ( f q ± f q ) , τ ± > 1 2 .
We operate with the two positive combinations Λ ± ( r , t ) and their product Λ ( r , t ) :
Λ ± ( r , t ) : = ( τ ± ( r , t ) 1 2 ) , Λ ( r , t ) : = Λ + Λ .
The symmetric equilibrium component is denoted e q + ( r , t ) , while its anti-symmetric counterpart is e q ( r , t ) :
e q + ( r , t ) = t q P ( r , t ) + I n s e E q + ( r , t ) , P = c s 2 ρ ,
e 0 ( r , t ) = ρ ( r , t ) 2 q = 1 Q m / 2 e q + , ρ = q = 0 Q 1 f q ( r , t )
E q + ( r , t ) = t q 3 ( j · c q ) 2 | | j | | 2 2 ρ ^ , I n s e = { 0 , 1 } ,
e q ( r , t ) = t q j · c q + E q ( f ) , j ( r , t ) = J + 1 2 F ,
E q ( f ) ( r , t ) = t q Λ F · c q , J ( r , t ) = q = 1 Q 1 f q c q .
The squared speed of sound  c s 2 relates pressure P ( r , t ) to the quantity of local mass ρ ( r , t ) ; the value of c s [0, 1] is adjustable inside its model-dependent stability interval [81,97] but c s 2 = 1 3 improves the anisotropy of the third-order (advection) truncation correction [41,42]; this choice is applied unless indicated. The so-called Navier–Stokes ( N S E ) term E q + ( r , t ) is computed with a given (reference) density value ρ ^ ( r , t ) : = ρ 0 for the modeling of incompressible or weakly compressible flow, while the (linear) Stokes equation is modeled without this term ( I n s e = 0 ). The local momentum quantity is J ( r , t ) , while the macroscopic momentum solution j ( r , t ) sums J ( r , t ) with the half forcing 1 2 F ( r , t ) ; taking into account this correction, external force F ( r , t ) is incorporated into Equation (2) via the term of E q ( f ) ( r , t ) using Λ of Equation (3).
The two weight components t q = { t c , t d } for the “coordinate” (parallel to an axis) and the “diagonal” (others) discrete speeds are fixed by the two hydrodynamic isotropic conditions:
2 q = 1 Q m / 2 t q c q α c q β = δ α β , α , β ; 2 q = 1 Q m / 2 t q c q α 2 c q β 2 = 1 3 , α β .
They provide t q = { 1 3 , 1 12 } in d 2 Q 9 , t q = { 1 3 , 1 24 } in d 3 Q 15 , and t q = { 1 6 , 1 12 } in d 3 Q 19 . The macroscopic equations are derived with the help of Chapman–Enskog non-equilibrium expansion [44,47] from the solvability conditions of the microscopic, mass- and momentum conservation equations satisfied by Equations (2)–(5):
n ^ 0 + ( r , t ) + 2 q = 1 Q m / 2 n ^ q + ( r , t ) = 0 ;
2 q = 1 Q m / 2 n ^ q ( r , t ) c q = F ( r , t ) .
Giving characteristic length L and a small perturbative parameter ε = 1 / L , the second-order accurate T R T approximation of the weakly compressible Navier–Stokes equation is derived [44] in the form
t ρ + · j = O ( ε 3 ) , t j + I n s e · ( j j ρ ^ ) = P + F + · ( ν j ) + ( · ν ξ j ) + O ( ε 3 ) + O ( u 3 ) , P = c s 2 ρ , ν = Λ + 3 , ν ξ = Λ + ( 2 3 c s 2 ) .
An incompressible or weakly compressible N S E becomes expressed in terms of velocity solution u ( r , t ) = j / ρ 0 providing ρ ^ = ρ 0 . Kinematic viscosity ν ( r , t ) is defined via local relaxation value Λ + ( r , t ) in Equation (3), while Λ ( r , t ) , and hence Λ ( r , t ) , is freely adjustable. The single-relaxation rate ( B G K ) operator [41] operates with τ ± = τ ; the B G K is then depleted from relaxation freedom Λ in Equation (3) and cannot improve accuracy, consistency, and stability with its help [44,81]. However, T R T sets the kinematic ( ν ) and bulk ( ν ξ ) viscosities with the same relaxation rate τ + ; this restriction is relaxed with the multiple-relaxation time ( M R T ) collision operators [42,43,78,81,98]. M R T computes collision directly in the Q × Q eigenbasis. It is composed of the (i) Q m / 2 + 1 even-order (symmetric) orthogonal polynomials of the velocity components c q α and (ii) their Q m / 2 odd-order (anti-symmetric) counterparts; corrections n ^ q + and n ^ q in Equation (1a) then sum the post-collision projections on the symmetric and anti-symmetric eigenvectors, respectively. M R T can also include E q ( f ) in e q in Equation (4d) by providing the same value τ for the relaxation rates of all anti-symmetric eigenvectors, including the d conserved ones, presented by the components of the given d dimensional discrete velocity set (see Equation (B.1) [79]).
In this work, the computations are performed with the d 3 Q 15 M R T model in the T R T configuration [78] [respective notations are τ + : = 1 / λ e , τ : = 1 / λ o , Λ : = 3 4 Λ 2 in terms of their parameter Λ 2 : = 4 3 ( 1 λ e + 1 2 ) ( 1 λ o + 1 2 ) ].

2.2.3. The bc-Step

Recall that c q C ( r b , t ) is called cut link in a fluid boundary node r b ( t ) V b ( t ) when r q ( t ) = r b + δ q ( r b , t ) c q S ( t ) with δ q ( r b , t ) [ 0 , 1 ] . The bc - s t e p aims to prescribe the Dirichlet velocity condition at the bisection point r q ( t ) as
u ( r q ( t ) , t ^ ) = u b ( r q ( t ) , t ^ ) , q C ( r b , t ) , t ^ [ t , t + 1 ] .
The task is to determine f q ( r b , t + 1 ) such that Equation (8) is satisfied with a given accuracy in terms of directional Taylor expansion. A popular, local and mass preserving boundary scheme, known as the “bounce-back rule” B B , simply reflects an outgoing population, and thus it does not need to compute δ q :
B B : f q ( r b , t + 1 ) = f ^ q ( r b , t ) 2 t q ρ 0 u b ( r q ( t ) , t ^ ) · c q .
The B B satisfies Equation (8) exactly halfway at r q ( t ) = r b + 1 2 c q only in a linear straight Couette flow, while a parabolic channel (Poiseuille) profile is then only modeled exactly provided that Λ = 3 16 , [44,79,92,99]. In general, the B B produces Λ -dependent velocity solution, which is accurate to second-order, as the best, in a regular geometry, but only first-order accurate in a grid-inclined channel or in irregular porous structure [81]. The linear and parabolic, “multi-reflection” ( M R q ) Dirichlet rules [44,48,78,79] advance the B B in accuracy by computing f q ( r , t + 1 ) via a prescribed linear combination of population components that move downstream and upstream cut link c q (cf. Equation (21) [79]):
f q ( r b , t + 1 ) = M R q ( r b ) + W q ( r q , t ^ ) , M R q ( r b ) = α ^ f ^ q ( r b , t ) + β ^ f ^ q ( r b , t ) + β ( I 1 f q ( r b , t + 1 ) + ( 1 I 1 ) f q ( r b , t ) ) + γ ( I 2 f q ( r b c q , t + 1 ) + ( 1 I 2 ) f q ( r b c q , t + 1 ) ) + γ ^ f ^ q ( r b c q , t ) + F ^ q ( r b , t ) , { I 1 , I 2 } { 0 , 1 } .
Here, solutions f q ( r b , t + 1 ) = f ^ q ( r b c q , t ) ( I 1 = 1 ) and f q ( r b c q , t + 1 ) = f ^ q ( r b 2 c q , t ) ( I 2 = 1 ) are propagated with Equation (1b) before the bc - s t e p from the upstream fluid neighbors r b c q and r b 2 c q , respectively. Coefficients { α ^ , β , β ^ , γ ^ , γ } q ( r b , t ) are assigned per cut link according to value δ q ( r b , t ) and a pre-selected boundary rule; their underscript q is omitted hereafter. M R q is defined [79] as a single-node rule when γ ^ = γ = 0 and I 1 = 1 , and as a local single- node rule when γ ^ = γ = 0 and I 1 = 0 (or β = 0 ); these two families reduce to B B with α ^ = 1 , β = β ^ = 0 . The local post-collision correction F ^ q ( r b , t ) is computed with Equation (1a) prescribing two coefficients K ± = K k ± ( r b , t ) per cut link c q :
F ^ q ( r b , t ) = K + n ^ q + ( r b , t ) + K n ^ q ( r b , t ) .
The term W q ( r q , t ^ ) in Equation (10) is prescribed as
W q ( r q , t ^ ) = α ( p ) e q + ( r b , t ) α ( u ) e q ( ρ ^ u b ( t ^ ) ) .
The first term aims to make vanish the symmetric (pressure) contribution from the underlying expansion of the Dirichlet M R condition, while the second term prescribes the Dirichlet velocity value from Equation (8) providing ρ ^ = ρ 0 for incompressible or weakly compressible flow. The two coefficients α ( p ) = α q ( p ) ( r b , t ) and α ( u ) = α q ( u ) ( r b , t ) are defined per cut link c q as
α ( p ) = α ^ + β + γ + β ^ + γ ^ 1 .
α ( u ) = α ^ + β + γ β ^ γ ^ + 1 .
In practice, α ( p ) is set by a pre-selected class of the boundary schemes, while α ( u ) is freely tunable and it scales the Dirichlet velocity condition prescribed by Equations (10) with (12), which explains an infinite number of formally equivalent Dirichlet boundary rules [44,48,79,94]. They have recently been classified according to accuracy criteria [79], and the main boundary classes are described and extended in Appendix B.2. First, they satisfy the parametrization condition, which means that Reynolds number R e and relaxation freedom Λ fix the steady-state momentum solution on a given grid. Second, they obey such properties as (i) an exactness in an arbitrarily inclined linear Stokes profile (the c - flow criteria satisfied by all considered schemes starting from the L I 1 + class) and (ii) its parabolic counterpart (the p - flow criteria by a single-node class L I 3 + ); (iii) a removal of the first gradient q e q + from the closure relation ( p - pressure criteria by a single-node class L I 4 + ); a combination of the last two properties ( p - p = p - pressure + p - flow by the two-node classes M R and L I 5 ), and an additional removal of the second gradient q 2 e q + ( p - p - p criteria by the most accurate class A V M R ). Each class L I k + contains both the single-node class L I k and the local single-node class E L I k + which extends the schemes [94]. L I 5 is formally introduced in this work; it updates class M L I [44] in the form which complements class L I 4 with directional finite-difference parabolic correction 1 2 α ( u ) δ q 2 q 2 e q ( r b , t ^ ) ; further details can be found in Appendix B.2.
Before performing Equations (10)–(12) at a given iteration, currently cut link c q C ( r b , t ) are sub-divided into the three groups:
1.
the no - nb ( corner )   cut   link has an upstream solid neighbor r b c q ;
2.
the one - nb   cut   link has an upstream fluid neighbor r b c q but the solid second neighbor r b 2 c q ;
3.
the two - nb   cut   link has two upstream fluid neighbors r b c q & & r b 2 c q .
Each group can prescribe its own boundary rule, typically as follows:
1.
a no - nb ( corner )   cut   link : L I k + with I 1 = 0 in Equation (10).
2.
a one - nb   link : L I k + or L I 5 with I 1 = 1 ;
3.
a two - nb   link : L I 5 with I 1 = 1 or A V M R with I 1 = I 2 = 1 .

2.2.4. The Reconstruction r-Step

The “fresh” nodes R ( t + 1 ) = { r n ( t + 1 ) } are new fluid boundary nodes, and their population solution q = Q m / 2 Q m / 2 f q ( r n , t + 1 ) is to be defined. For this purpose, their Q discrete velocities are in turn sub-divided into three groups:
1.
a tangential link c q T ( r n ) has an upstream “fresh” neighbor r n c q R ( t + 1 ) ; in particular, an immobile velocity c 0 T ( r n ) .
2.
a no - nb ( corner )   cut   link c q has an upstream solid neighbor, hence r n + δ ± q c ± q S ( t + 1 ) with δ ± q [ 0 , 1 ] .
3.
the remaining normal links are sub-divided into the stream links , where solution is to be obtained by propagation from an upstream neighbor, and those which are opposed to the cut links, where solution is to be reconstructed with the bc - s t e p .
The “(re)fill” f - s t e p employs the following techniques:
1.
the f - stream sub-step: f q ( r n , t + 1 ) = f ^ q ( r n c q , t ) ;
2.
the f - bc sub-step: f q ( r n , t + 1 ) is defined by Equation (10) using directional approximations from Appendix B.3.1 for some unknown components;
3.
the f - equil sub-step: f q ( r n , t + 1 ) is set equal to an equilibrium value with Equation (A2), where ρ ( r n , t + 1 ) is defined in Equation (A3) and it sums the reconstructed solution, ρ ( r n , t + 1 ) = q = 0 Q m f q ( r n , t + 1 ) ;
4.
the f - extr sub-step: f q ( r n , t + 1 ) is set equal to an averaged extrapolated solution f q ( r n , t ^ ) detailed in Appendix B.3.3;
5.
the f - bc - corner sub-step: the f ± q ( r n , t + 1 ) is defined with the single-node bc -rule using an extrapolated solution { f ± q ( r n , t + 1 ) } and { f ^ ± q ( r n , t ) } at the first step of an iterative update.
The three reconstruction r - s t e p s , (i) the r - equil , (ii) the r - extr , and (iii) the r - bc - corner , combine the first two sub-steps f - stream and f - bc with an optional combination of the last three sub-steps applied in tangential links and no - nb links :
1.
r - equil applies the f - equil sub-step with Equations (A2) and (A3) for all tangential and corner links following [78]. This algorithm requires velocity interpolation to “fresh” nodes as an initial estimate in an iterative update (see a - vel sub-step in Appendix B.3.1).
2.
r - extr assigns an extrapolated solution f q ( r n , t + 1 ) to the moving tangential and corner links with f - extr , then computes an actual velocity value u ( r n , t + 1 ) from the momentum of the reconstructed populations and applies f - equil with Equations (A2) and (A3) only to compute f 0 ( r n , t + 1 ) .
3.
r - bc - corner applies:
(a)
the f - extr sub-step for the moving tangential links ;
(b)
the f - bc - corner sub-step for the no - nb links using f - extr as an initial estimate in an iterative update;
(c)
the f - equil sub-step for f 0 ( r n , t + 1 ) using Equations (A2) and (A3) with an actual momentum value j ( r n , t + 1 ) of all moving links, already reconstructed.
The local iterative update ( r - iter ) is optionally applied in the “fresh” nodes with the ideas of the Local Iterative Refill ( L I R ) algorithms [95] but completely independently in one node compared to another. In our implementation, r - iter iteratively performs the local collision operator followed by the f - s t e p , without any streaming nor to or between new nodes. So, r - iter essentially serves us to update { f ^ q ( r n , t ) } and { n ^ q ± ( r n , t ) } , then to iterate their solution { f q ( r n , t + 1 ) } with the f - s t e p .
The robustness of the f - bc - corner improves when an (unknown) post-collision correction F ^ q in Equation (11) is set to zero in an initial step, then r - iter applies its updated values. Note that the populations obtained with the f - extr sub-step are not updated by the r - iter , because the extrapolations only involve the “old” fluid nodes; at the same time, an equilibrium reconstruction r - equil can benefit, in the presence of the moving tangential links and/or no - nb links , from an iterative update of momentum in Equations (A2) and (A3).
Further details on the reconstruction step can be found in Appendix B.3.

2.2.5. Pre-Selected Algorithms

The benchmark simulations in the deformable channel flow are conducted in Appendix B.4. The three reconstruction algorithms, r - equil , r - extr and r - bc - corner , are applied there with the three sub-steps of the r - iter update, first confirming that r - iter effectively reduces the sharp negative velocity deviations near horizontal walls. The different combinations of the boundary rules compete in these simulations, from the least accurate bounce-back rule ( B B ) to the most accurate p - p - p class A V M R , through c - flow B F L 0 L I 0 + [100], p - pressure B F L 4 L I 4 and p - p B F L 5 L I 5 . B F L k fixes coefficients α ^ γ ^ in Equation (10) according to the fractional B F L rule (see Table III [79] and Equation (A1)), with F ^ q in Equation (11) corresponding to class L I k and F ^ q = 0 in the original (non-parametrized) B F L L I 0 + [100]. Indeed, our preliminary tests indicate that B F L k behaves more accurately than Y L I k [101] and Z L I k [102] (using α ( u ) = 1 ) from Table III [79], at least when Λ = 3 16 where B F L ( δ = 1 2 ) = B B is exact for a mid-grid location of the horizontal walls, as it is prescribed in this benchmark simulation.
The reported results guide us to give preference to the reconstruction algorithm r - bc - corner combined with the 1–3 iterations of the local collision update r - iter , because in these experiments it produces most accurate velocity fields among the three reconstruction algorithms. These simulations also suggest that when developing the method for a given problem, the simplest strategy is to first apply B F L 0 for all cut links, because this linearly accurate c - flow scheme [100] involves no post-collision correction and then applies equally with any collision operator. Once the simulations produce reasonable results, a switch to B F L 5 - B F L 5 - B F L 4 or A V M R - B F L 5 - B F L 4 could be approved, respectively, for two - nb link , one - nb link and no - nb ( corner ) cut link .
Indeed, based on the above results, B F L 4 applied in the corner links improves the accuracy, advising to omit its post-collision correction K 4 n ^ q from Equation (11) during the first reconstruction step in new corner nodes with r - bc - corner . In turn, B F L 5 and A V M R are the most accurate for velocity and pressure in regular boundary nodes. B F L 5 involves a velocity solution from the two opposite directional neighbors, where it may benefit from the prescribed Dirichlet values even in corners. In our simulations below, we switch B F L 5 to B F L 4 in “fresh” corners r n , because our linear velocity interpolation for u ( r n ) , used at first sub-step of r - iter , (optionally) implies the same Dirichlet values as the parabolic correction of B F L 5 ; however, in principle, it should be possible to apply B F L 5 to the new corner links in the following sub-steps.

3. Results

This section describes the L B M - D B simulations of the colonic circular contractions in the DCM-type numerical experiments. Section 3.1 defines geometry and five occlusion scenarios from Table 4; Section 3.2 provides implementation of L B M - D B ; Section 3.3 estimates the parameter space; Section 3.4 describes the numerical experiments; Section 3.5 visualizes and analyzes the results obtained with motility pattern M O T - I ; Section 3.6 compares them to the results obtained with motility pattern M O T - I I ; they are described in detail in Appendix C.1.

3.1. Occlusion Scenarios and Geometry

L B M - D B operates inside an open (lumen) counterpart of the colonic cylindrical tube segmented along the longitudinal z-axis (which is the x-axis in DCM notations) into 10 haustral units of individual length l = 2 cm . The tube is closed by either an impermeable (zero velocity) Dirichlet condition or a periodic condition at the inlet and outlet; the two small neutral occlusion buffers, each of total length l n = 0.4 cm (or l n = 0.6 cm), are placed at the impermeable (or closed, respectively) ends. The transverse shape is dynamically occluded by deformation of the bounding membrane, circular (“C”) or parabolic (“P”), according to one of the two motility protocols in Table 1, M O T - I and M O T - I I . The “C” shape can vary its degree of occlusion o ( t ) between zero (open pipe) and one (closed pipe), while L B M - D B increases the experimental degree of occlusion to o ( t ) 0.938 . However, the “P” cross-section restricts o ( t ) to interval [ o min = 1 s m a x , o max = 1 s m i n ] , where the minimum occlusion degree o min 0.035 ( s m a x = 7 3 4 π 0.965 ) occurs when the three parabolic curves touch an outer circle, while the maximum occlusion degree o max 0.862 ( s m i n = 3 4 π 0.138 ) meets their limit case without overlapping.
The five occlusion scenarios, from the O C L - 1 to the O C L - 5 , are ordered in Table 4 by increasing the forward degree of occlusion o f ; they are shown in Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9 with the “P” and “C” shapes, respectively, and described as follows:
1.
O C L - 1 in Figure 5 sets the relative occlusion interval according to Figure 1 [68], freely prescribing o f = o Δ , which is probably slightly smaller than o f in Figure 5 [9].
2.
O C L - 2 in Figure 6 constructs the “P” contour from the distances provided in Figure S1 [9]. According to Table A1, distance | l ( t ) | to the parabolic boundary along the median linearly expresses through the sector value S ( t ) of the open surface, with S ( t ) = 1 3 ( 1 o ( t ) ) in the considered symmetric case where the solution reads
| l ( t ) | = 1 24 | 3 + 4 3 ( o ( t ) 1 ) π ) r | .
Replacing o ( t ) in Equation (14) by o f , o n and o r , the two distances from Figure S1 [9], (i) | l r | | l f | = 1 [ cm ] , and (ii) | l r | | l n | = 2 5 [ cm ] , give two equations with respect to three unknown fractions, like (i) 5 3 ( o f o n ) π r = 12 and (ii) 3 ( o f o r ) π r = 6 , where we substitute r = l Δ × 3 3 . Two other conditions [9] are (iii) o f , n = 2 5 and (iv) o r , n = 2 5 . These four equations cannot be satisfied simultaneously, O C L - 2 applies the appropriate solution to Equations(i)–(iii), where o n = 30 19 π < o Δ and then [ o r , n , o f , n ] = [ 3 5 , 2 5 ] ; a deviation from [ o r , n , o f , n ] = [ 2 5 , 2 5 ] [9] should occur due to a (slight) difference in the cross-section shape and surface value, but may also apply non-identical definitions of the occlusion degree.
3.
O C L - 3 in Figure 7 sets o f = o max ; this value must exceed the highest degree of occlusion in the DCM and MRI experiments; o r and o f are then assigned providing interval [ o r , n , o f , n ] = [ 1 5 , 3 5 ] [68].
4.
O C L - 4 in Figure 8 applies o f = o max , as O C L - 3 , but it shares interval [9] [ o r , n , o f , n ] = [ 3 5 , 2 5 ] with O C L - 2 . The neutral occlusion o n is highest than in all other cases, and it is slightly larger than the triangular one, with relative difference o n , Δ = 1 / 20 .
5.
O C L - 5 in Figure 9 sets o n = o Δ and relative interval [68] [ o r , n , o f , n ] = [ 1 5 , 3 5 ] in the circular shape, as O C L - 1 and O C L - 3 ; “C” then operates the highest degree of occlusion o f = 0.938 > o max .
In summary, the most challenging scenario in the “P” tube is O C L - 4 due to its largest variation in degree of occlusion from o r 0.246 to o f = o max 0.862 , where the three parabolas form acute angles touching each other. The most difficult scenario in the “C” tube is O C L - 5 with o f 0.938 , which certainly exceeds the highest experimental occlusion degree [9,68].

3.2. Membrane Deformation with the LBM-DB

The collision operator d 3 Q 15 M R T is applied in the two-relaxation time configuration [44,78] to model a single-phase, full-volume, weakly compressible Newtonian fluid through the lumen aperture where the motility protocol prescribes the membrane position at each time step. When the membrane contour is circular (“C”), boundary velocity u b ( r ( t + 1 ) r ( t ) ) n is defined by giving radius values r ( t ) and r ( t + 1 ) before and after membrane displacement, respectively; the unit normal vector n = ( n x , n y , 0 ) is drawn along the radius towards the bisection point ( r b + δ q c q , t + 1 ) .
When the membrane’s contours p ( t ) and p ( t + 1 ) are parabolic (“P”), the boundary velocity is defined approximately, as u b = | P ( t + 1 ) P ( t ) | n by giving the cut link bisection point P ( t + 1 ) = ( r b + δ q c q , t + 1 ) p ( t + 1 ) , drawing the normal from P ( t + 1 ) to the side of an inscribed triangle sharing the two vertices with P ( t + 1 ) , up to the bisection at P ( t ) p ( t ) , and defining the unit normal vector n = ( n x , n y , 0 ) to be directed from P ( t ) to P ( t + 1 ) . Once the Dirichlet velocity u b ( r b + δ q c q , t ^ ) is defined, Equations (10) and (11) are performed routinely with the pre-selected boundary rule; the combination B F L 5 B F L 5 B F L 4 is applied with Equation (A1) in two - nb , one - nb and no - nb cut link, respectively, in both fluid and “fresh” boundary nodes.
In the regular situation, when the cut links bisect the membrane contour at r b + δ q c q with δ q [ 0 , 1 ] , the Dirichlet velocity u b = ( u x , u y , 0 ) of the membrane displacement at this wall point is to be prescribed with Equation (8). However, it may happen in the first and the last sections of a given segment that cut link c q , pointing into the occluded domain in the neighbor segment, does not bisect the membrane contour at r b + δ q c q with δ q [ 0 , 1 ] , because the contractions are uniform per segment, and hence only piece-wise continuous; the membrane displacement in the neighbor segment may then exceed 1 l . u with respect to the location of boundary point r b in a given segment. In such a case, outgoing population f ^ q ( r b , t ) is reflected back with the B B rule in Equation (9); this mass conserving rule approximately traces an impermeable barrier halfway between the exit of the more relaxed membrane and its more occluded neighbor.

3.3. Parameter Space

The size of the cuboid computational grid is set equal to 1 l . u (lattice unit); the length scale L = L l from the physical to the lattice units applies with L = 10 [ l . u ] ( l . u = space step) per haustra length l = 2 cm ; an outer cylinder’s radius r = l Δ × 3 3 ( l = 3.8 cm ) is then discretized with R = r L 10.97 [ l . u ] . The characteristic velocity is set equal to P P W speed v w = l t f [ cm / s ] ; accordingly, L B M - D B prescribes V w = L T f giving haustra length L and T f . The Mach number M a = V w / c s controls the degree of the compressibility level with value c s of the sound speed in Equation (4a).
Unless indicated, the results are reported with T f = 1280 [ l . u ] ( l . u = time step) and c s 2 = 1 3 ; these parameters yield V w = 7.8125 × 10 3 [ l . u ] and M a = 1.35 × 10 2 . Accordingly, the three motility regimes in Table 1, S L O W , F A S T and VMAX, operate with the following scale factors:
  • the time scale from the physical to lattice units:
    T = T f t f = { 256 , 512 , 1280 } [ l . u s ] when t f = { 5 , 2.5 , 1 } [ s ] ;
  • the velocity scale from the lattice to physical units: U 1 = v w V w = { 51.2 , 102.4 , 256 } [ cm / s l . u ] when v w = { 0.4 , 0.8 , 2 } [ cm / s ] .
LBM kinematic viscosity ν = 1 3 ( τ + 1 2 ) [l.u] is prescribed equating Reynolds number R e = L V w ν to its physical counterpart Re p h = l v w ν p h , with ν p h = { ν H , ν L } . The H I V I S lattice values are then ν { 8.14 , 4.07 , 1.63 } × 10 2 [ l . u ] , and the L O V I S lattice values are ν { 2.39 , 1.2 , 0.48 } × 10 2 [ l . u ] , respectively, with S L O W , F A S T and V M A X wave speeds. These simulations then apply τ + [ 0.514 , 0.744 ] from V M A X L O V I S to S L O W H I V I S . τ is computed from condition Λ = 3 16 in Equation (3); this “magic” relationship is best known for providing an analytical solution in straight Poiseuille flow using the bounce-back rule and mid-grid wall location, [78,79,81,92]. The local pressure solution is defined in lattice units as P ( r , t ) = c s 2 ( ρ ( r , t ) ρ m e a n ( t ) ) , where ρ m e a n ( t ) averages the density field ρ ( r , t ) over all fluid nodes. The pressure physical solution p [ mmHg ] = ρ ρ 0 × v w 2 V w 2 P × 7.500615 × 10 3 [l.u] is computed with H I V I S and L O V I S density value ρ [ k g / m 3 ] [68] providing ρ 0 = 1 [ l . u ] .

3.4. Numerical Experiments

The neutral occlusion is initialized with zero fluid velocity and constant pressure P = c s 2 ρ 0 ; the motility cycle over the ten haustral units lasts t = T 10 [ s ] , and it is followed by an equilibration, with no membrane contractions, towards an initial state for T e q = { 20 , 10 , 4 } [ s ] , with S L O W , F A S T and V M A X , respectively (see Table 1). The data are output each 128 iterations, or { 0.5 , 0.25 , 0.1 } s , respectively. To compare, temporal resolution [68] is approximately 2 s (MRI) and 0.25 s ( D C M D T ) with S L O W and F A S T P P W .
Each numerical experiment combines one of two motility patterns from Table 1 with the five occlusion protocols from Table 4; the results of these ten experiments are provided in Table A2, Table A3, Table A4, Table A5, Table A6, Table A7, Table A8, Table A9, Table A10 and Table A11 in the (central) segment s . 6 . Each of these experiments is run six times combining the two viscosity values ν p h , H I V I S and L O V I S , with the three travel velocities v w , S L O W , F A S T and VMAX. Finally, each of these 10 × 6 configurations was simulated in the four geometries as follows: (i) P - C L O and (ii) C - C L O , in the closed pipe, and (iii) P - P E R and (iv) C - P E R , in the periodic pipe.
Our analysis addresses the distributions of stream-wise velocity v ( r , t ) = u z ( r , t ) and pressure p ( r , t ) , where Table A2, Table A3, Table A4, Table A5, Table A6, Table A7, Table A8, Table A9, Table A10 and Table A11 provide the following intervals: (i) [ v ¯ m i n , v ¯ m a x ] v w , (ii) [ v m i n , v m a x ] v w and (iii) [ p m i n , p m a x ] , achieved over the whole motility cycle by (i) segment-averaged velocity v ¯ ; the segment minimum and maximum values of (ii) velocity and (iii) pressure. Additionally, R ( v ¯ ) = | v ¯ m i n | v ¯ m a x and R ( v ) = | v m i n | v m a x are reported; they represent the ratios of the highest retrograde amplitudes  | v ¯ m i n | = | min t ( v ¯ ( t ) ) | and | v m i n | = | min t ( v m i n ( t ) ) | to their antegrade counterparts, v ¯ m a x = max t ( v ¯ ( t ) ) and v m a x = max t ( v m a x ( t ) ) , respectively; Reynolds number R e ( | v ¯ m i n | ) = | v ¯ m i n | l ν p h is then computed with | v ¯ m i n | . The last column reports | δ M | which is the value of the relative mass change at the end of the motility cycle, t = T e n d = T 10 + T e q .
Table 5 summarizes results from Table A2, Table A3, Table A4, Table A5, Table A6, Table A7, Table A8, Table A9, Table A10 and Table A11 separately in four geometries and reports the extreme values over the six runs as follows: (a) min ( R e ( | v ¯ m i n | ) ) ; (b) max ( R ( v ¯ ) ) ; (c) max ( R ( v ) ) ; (d) max ( R e ( | v ¯ m i n | ) ) ; (e) max ( [ p m i n , p m a x ] ) ; (f) max ( | δ M | ) . Table 5 allows us for a quick comparison between (i) the two motility patterns, (ii) the “P” and “C” transverse shapes; (iii) the closed and periodic boundary conditions.
Herein, from top to bottom and right to left, the first and second diagrams display the segment-averaged (“s-mean”) distributions of pressure p ¯ ( t ) and the stream-wise velocity v ¯ ( t ) together in the ten haustral units, respectively. The third diagram displays together the pressure distributions in s . 6 : segment-averaged (“s-mean”, p ¯ ( t ) ), segment-maximum (“s-max”, p m a x ( t ) ), and segment-minimum (“s-min”, p m i n ( t ) ). The fourth diagram reports together the mean, maximum, and minimum, normalized stream-wise velocity distributions in s . 6 : v ¯ ( t ) / v w , v m a x ( t ) / v w and v m i n ( t ) / v w , respectively. The fifth diagram depicts the dynamics of relative occlusion o ( t ) / o n 1 and helps the visual synchronization of the physical and motility events. The sixth diagram displays the relative, exact, and approximate (measured via summation of the fluid grid nodes) evolution of the lumen volume with respect to its neutral value, δ V ( t ) = V ( t ) V ( o = o n ) 1 and δ V n u m = V n u m ( t ) V n u m ( o = o n ) 1 , respectively, also reporting relative mass variation δ M ( t ) = ρ ( t ) ρ ( t = 0 ) 1 .

3.5. Motility Pattern MOT-I

The results of motility pattern M O T - I are reported in Table A2, Table A3, Table A4, Table A5 and Table A6, illustrated in Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15, Figure 16, Figure 17, Figure 18, Figure 19, Figure 20 and Figure 21 and resumed in Table 5. Recall, M O T - I operates without pauses, with relax and forward times equal between them, t r = t f , and back relaxation time t b = 3 t f . According to Table 5, the maximum amplitude of the relative volume variation | min t δ V ( t ) | varies between | 6.55 % | ( O C L - 1 ) and | 16.4 % | ( O C L - 5 ), equally in the “parabolic” and cylindrical pipes; a very slight difference between the closed and periodic pipes is due to the presence of the two “solid” entry/exit units in the closed box.

3.5.1. The MOT-I ∪ OCL-1

Table A2 and Figure 10, Figure 11, Figure 12 and Figure 13 address the combination of M O T - I with the least occluding protocol O C L - 1 . Figure 10 displays the case H I V I S S L O W in the closed geometry with the “parabolic” cross-section ( P - C L O ). The two first diagrams show that the non-peripheral segments, s . 3 s . 8 , display similar segment-averaged distributions for pressure p ¯ ( t ) and velocity v ¯ ( t ) . Segment nine ( s . 9 ) displays negative value v ¯ m i n ( t ) of the largest amplitude, and hence the highest mean retrograde flux, during the forward contraction of its downstream neighbor s . 10 when t [ 45 , 50 ] s ; in turn, s . 10 noticeably increases its mean pressure amplitude p ¯ ( t ) during contraction, followed by the pressure drop of the same amplitude at the start of its back (return) relaxation at t = 50 s . The fourth diagram shows that the three minimum local velocity peaks v m i n ( t ) in s . 6 occur at times of the highest occlusion in the upstream ( s . 5 ), local ( s . 6 , highest amplitude) and downstream ( s . 7 ) segments when t { 25 , 31 , 35 } s , respectively. The third diagram confirms that p m i n ( t ) shows the local peaks just before the high reflux events; at the same time, the sharp peaks of p m a x ( t ) also appear during upstream and local contractions. After downstream contraction, the pressure and velocity distributions display a staircase relaxation composed of decaying amplitude average velocity jumps during subsequent downstream occlusion events. Towards the end of the motility cycle, the local antegrade (positive) velocity reaches its maximum amplitude when the pressure drops drastically at the last segments, at t = 50 s (see first diagram). Finally, an entire tube converges to a uniform pressure value and a zero velocity during equilibration at the neutral occlusion state.
The pressure and the v w normalised velocity distributions have similar shapes for all six ν p h v w configurations in the closed pipes of the two transverse shapes, P - C L O and C - C L O ; however, the faster wave speeds or less viscous flows exhibit dynamic profiles that are less “staircase” but more “nervous” (oscillating). Figure 11 exemplifies the case L O V I S F A S T in the closed circular pipe ( C - C L O ); here, v w increases by a factor of two compared to the previous figure, and hence the motility intervals all become twice as short. Overall, when increasing v w and/or decreasing viscosity value, the pressure distributions resemble a sequence of abrupt oscillations occurring during occlusion events; the retrograde distribution v m i n ( t ) drops abruptly when the upstream neighbor only begins to relax, and v ¯ ( t ) reaches its highest retrograde amplitude shortly after the local contraction, when the downstream segment starts to constrain. These two events are associated with the mean pressure increase and decrease, respectively, in the given segment s . 6 . The highest antegrade amplitude  v ¯ m a x ( t ) reaches its maximum value during the last contraction of s . 10 due to the pressure drop at that location.
Figure 12 and Figure 13 address the periodic pipes, C - P E R H I V I S S L O W and P - P E R L O V I S F A S T , respectively. Remarkably, in the periodic case, s . 9 shows only slightly higher pressure and retrograde velocity amplitudes compared to other segments; v m a x ( t ) now produces the two largest antegrade peaks shortly before the upstream and local occlusion events, and they meet with the p m i n ( t ) decrease in s . 6 . Finally, closed pipes alike, the two transverse shapes show similar solutions.
Table A2 quantifies these effects: it shows that the ratios of the highest retrograde and antegrade amplitudes, R ( v ¯ ) = | v ¯ m i n | v ¯ m a x , and especially R ( v ) = | v m i n | v m a x , vary weakly with ν p h and v w ; except in a few cases, their values smoothly increase from (i) V M A X to S L O W ; (ii) L O V I S to H I V I S ; (iii) “P” to “C”, and (iv) C L O to P E R . According to Table 5, the case H I V I S S L O W achieves in the three geometries (from the four) the largest values max ( R ( v ¯ ) ) [ 1.39 , 2.12 ] and max ( R ( v ) ) [ 1.57 , 1.69 ] . M O T - I O C L - 1 holds max ( R e ( | v ¯ m i n | ) ) [ 4.88 , 7.13 ] , with the minimum and maximum values in P - P E R and the C - C L O , respectively. In turn, the case L O V I S V M A X commonly displays the largest peak pressure interval max ( [ p m i n , p m a x ] ) and the highest Reynolds value Re ( | v ¯ m i n | ) .

3.5.2. MOT-I ∪ OCL-2

Table A3, Figure 14 and Figure 15 address motility pattern M O T - I subjected to occlusion scenario O C L - 2 , where the o n and o f increase relative to O C L - 1 , but o r decreases, and therefore the relaxation and contraction speeds become greater.
Figure 14 displays the case H I V I S S L O W in the periodic “parabolic” pipe P - P E R ; the velocity profiles in its circular counterpart C - P E R are very similar, but the pressure peaks are sharper and larger in the parabolic pipes. The velocity and pressure profiles retain their O C L - 1 shapes, but the two antegrade peaks of v ¯ m a x ( t ) , which occur just before the upstream and local events of the largest contraction, grow relative to the maximum amplitude | min ( v ¯ m i n ( t ) ) | of retrograde velocity, which again takes place at the sequential moments of the highest occlusion (fourth diagram). In turn, pressure peaks p m a x ( t ) increase and reproduce at the times of local and downstream contractions (third diagram).
Figure 15 addresses an “opposite limit” to H I V I S S L O W and displays the case P - P E R L O V I S V M A X . The velocity and pressure profiles become much more oscillatory between main occlusion events; velocity distribution retains the principal features of the slow viscous flow; on the other hand, the two large pressure peaks p m a x ( t ) disappear in the same way as in the case of L O V I S V M A X with O C L - 1 (cf. Figure 10, Figure 11, Figure 12 and Figure 13).
Table A3 confirms that the average velocity continues to increase with v w , whereas the velocity amplitude contrasts, R ( v ¯ ) and R ( v ) , show approximately the same values over the six ν p h v w configurations in a given geometry. In closed pipes, the largest ratio of the average amplitudes max ( R ( v ¯ ) ) [ 1.36 , 1.48 ] remains approximately the same as with O C L - 1 , but the largest ratio of peak amplitudes max ( R ( v ¯ ) ) 1 is reduced by a factor of 1.5 . This reduction becomes even more notable in the two periodic pipes P - P E R and C - P E R , where max ( R ( v ¯ ) ) and max ( R ( v ) ) are reduced by a factor of 1.5 to 3, respectively, and the value of max ( R ( v ) ) 0.65 confirms that the antegrade velocity peaks largely predominate over the retrograde peaks. Hence, we note that although O C L - 2 operates with larger relaxation and occlusion rates, this scenario reduces, or even partly reverses, the contrast between retrogade and antegrade amplitudes, particularly in periodic pipes.
The case L O V I S V M A X reaches the highest retrograde velocity amplitude, then the largest interval max ( R e ( | v ¯ m i n | ) ) [ 7.66 , 13.6 ] from P - P E R to C - C L O , which increases max ( Re ( | v ¯ m i n | ) ) with the previous scenario O C L - 1 by a factor of about two. The case L O V I S V M A X again reaches the largest pressure intervals, and they also increase their values by a factor of about two relative to O C L - 1 .

3.5.3. MOT-I ∪ OCL-3

Table A4, Figure 16 and Figure 17 combine motility pattern M O T - I with occlusion scenario O C L - 3 , where the highest occlusion degree o f = o max reaches its maximum available value in the “P” shape before the parabolic contours overlap; at the same time, o n and o r increase compared to the previous two cases according to occlusion interval [ o r , n , o f , n ] = [ 0.2 , 0.6 ] , which is common with O C L - 1 . Table 5 reports that the highest retrograde velocity amplitude increases and case L O V I S V M A X produces Re ( | v ¯ m i n | ) [ 31.2 , 32 ] and Re ( | v ¯ m i n | ) [ 22 , 16.9 ] in the closed and periodic pipes, respectively. Now, there is a notable difference between periodic and closed pipes, where P - C L O and C - C L O produce much higher ratios max ( R ( v ¯ ) ) and max ( R ( v ) ) due to their larger retrograde averaged and peak values, | v ¯ m i n | and | v m i n | , respectively, especially in the H I V I S flow. On the other hand, closed and periodic “P” pipes reach the largest pressure range, such that p m a x p m i n 1.6 mmHg with the H I V I S VMAX, or about five times greater than the largest pressure intervals reached, as previously in the case of L O V I S VMAX, in circular pipes.
Figure 16 and Figure 17 display the case H I V I S V M A X in P - C L O and P - P E R , respectively. In these two geometries, s . 6 shows a very high positive pressure peak p m a x ( t ) shortly before the highest occlusion o f ; this peak should explain the largest antegrade velocity peak. In turn, the (negative) peaks of p m i n ( t ) reach their largest amplitudes immediately after the local and downstream contractions, and they should explain the local (negative) peaks in retrograde velocity. Herein, the antegrade peaks v m a x ( t ) reach the same amplitude in closed and periodic pipes, but the amplitude of the retrograde flow is higher in the closed pipe; this explains the several-times-larger contrast ratios max ( R ( v ) ) and max ( R ( v ¯ ) ) in P - C L O , which further increases max ( R ( v ¯ ) ) 3.1 in the circular pipe ( C - C L O ) according to Table 5.
Finally, an amplitude of change in relative mass | δ M | increases by an order of magnitude compared to O C L - 2 and remains several times larger in “parabolic” pipes than in circular ones. At the same time, maximum amplitude | min t δ V ( t ) | of the relative (negative) volume variation doubles that of O C L - 2 , equally in all four geometries, and reaches the level of 13 % (cf. last diagrams in Figure 16 and Figure 17, Table 5 and Table A4). Hence, | δ M | grows faster than | min t δ V ( t ) | with the degree of occlusion and especially when the pressure oscillation becomes important due to a non-smooth form of the parabolic opening.
In summary, this example teaches us that the “parabolic” shapes significantly increase their pressure intervals compared to previous occlusion scenarios on the one hand and relative to cylindrical tube on the other hand, when the degree of contraction occlusion o f reaches its maximum value o max with the three non-superimposed parabolas (cf. Figure 7). In these cases, the resulting pressure peaks p m a x ( t ) have similar amplitudes in the closed and periodic “parabolic” pipes, but the retrograde flux at the highest contraction moment is larger in the closed pipe.

3.5.4. MOT-I ∪ OCL-4

Table A4, Figure 18 and Figure 19 report results when M O T - I is combined with O C L - 4 . According to Table 1, O C L - 4 shares o f = o max with O C L - 3 and the relative interval [ o r , n , o f , n ] = [ 0.6 , 0.4 ] with O C L - 2 ; it operates then the highest neutral degree of occlusion o n and the fastest contraction/relaxation speeds compared to the three previous occlusion scenarios. However, Table 5 reports that this combination only very slightly increases the Reynolds range, and hence the largest amplitude of the averaged retrograde velocity compared to O C L - 3 , as the case L O V I S V M A X reaches almost the same interval max ( R e ( | v ¯ m i n | ) ) [ 31 , 34 ] in the closed pipes, and a bit larger in the periodic tubes, as max ( R e ( | v ¯ m i n | ) ) [ 17 , 22 ] ( O C L - 3 ) versus max ( R e ( | v ¯ m i n | ) ) [ 22 , 29 ] ( O C L - 4 ), and Re ( | v ¯ m i n | ) mostly decreases from the circular to the “parabolic” shape. However, both velocity ratios max ( R ( v ¯ ) ) and max ( R ( v ) ) decrease relative to O C L - 3 , similar to what was observed for their less occluded counterparts, O C L - 2 and O C L - 1 , respectively. At the same time, and even stronger than between O C L - 1 and O C L - 3 , both velocity contrasts increase their values with an increase in the degree of neutral occlusion at the same relative intervals [ o r , n , o f , n ] , as here by the factor of approximately 1.5 from O C L - 2 to O C L - 4 . Finally, the case H I V I S V M A X reaches the highest pressure intervals p m a x p m i n 2.2 mmHg over all considered cases; as for O C L - 3 , the largest pressure oscillations take place in “parabolic” pipes with the closed and periodic ends.
Figure 18 and Figure 19 exemplify this last case H I V I S V M A X in the closed, circular and “parabolic” pipes, respectively. As before, the closed circular pipe displays the velocity shape characterized by the two large antegrade peaks v m a x ( t ) that occur shortly before the two local occlusion events, as (i) when the upstream neighbor contracts and propels the flow downstream into (relaxed) s . 6 and (ii) when the s . 6 contracts itself but s . 7 relaxes. At the same time, retrograde peaks v m i n ( t ) remain distinctive over a series of occlusion events, and they display the global minimum value min t ( v ¯ ) shortly after the moment of the highest contraction due to sudden drop in p m i n ( t ) . These occlusion events are all more pronounced in the “parabolic” pipe in Figure 19 where, remarkably, the largest positive pressure peak p m a x ( t ) , occurring just prior to the strongest contraction, is one order of magnitude larger than in its circular counterpart, in agreement with the data analysis above.

3.5.5. MOT-I ∪ OCL-5

Table A6, Figure 20 and Figure 21 combine M O T - I with the last occlusion scenario O C L - 5 from Table 1; O C L - 5 operates the higher occlusion values than the previous four scenarios, as o f 0.94 and o r 0.47 , and this motility protocol is only available in circular pipes since o f > o max . Table 5 reports that case C - C L O L O V I S V M A X reaches the highest retrograde amplitude and Reynolds value max ( Re ( | v ¯ m i n | ) ) 69.5 , twice as large as its counterparts O C L - 3 and O C L - 4 , and this occlusion protocol makes it possible to obtain the strongest retrograde/antegrade amplitude contrasts, as max ( R ( v ¯ ) ) 5.2 and max ( R ( v ) ) 4.3 , which are achieved in the low-viscous flow with wave speeds S L O W and VMAX, respectively. In contrast, circular periodic pipes show a relatively small increase in velocity contrasts compared to the previous scenario, O C L - 4 ; an antegrade maximum velocity amplitude v m a x still mainly dominates its retrograde counterpart, O C L - 4 alike, but this situation is reversed in the high-viscous flow at the two slowest wave speeds, with max ( R ( v ) ) = 1.03 obtained in the case H I V I S S L O W . The realized max ( Re ( | v ¯ m i n | ) ) 28.8 achieved in the L O V I S V M A X case is less than half of its closed pipe counterpart, but the obtained Re ( | v ¯ m i n | ) range exceeds all previous results in circular periodic pipes. It is interesting to note that despite a noticeable increase in the retrograde velocity amplitude, the maximum pressure intervals retain their O C L - 4 range; at the same time, like the retrograde velocity amplitude, the maximum mass change almost doubles its magnitude compared with O C L - 4 and reaches max | δ M | | 3.9 % | in the case C - C L O H I V I S S L O W again, while the same case retains the O C L - 4 range in periodic pipes, as max | δ M | | 1 % | .
Figure 20 and Figure 21 illustrate the case L O V I S S L O W in closed and periodic circular pipes, respectively. They show that the maximum of the retrograde amplitude | min t v ¯ | is about four times higher in the closed pipe than in the periodic one due to the larger value of retrograde peak | min t v m i n | in the closed case. Herein, while in periodic geometry these two highest amplitude retrogrades minimum are located just before the local and downstream occlusion peaks, the global retrograde minimum value min t v m i n is displayed directly at the time of highest occlusion in the closed tube. At the same time, the two local antegrade peaks of v ¯ m a x ( t ) , occurring during upstream and local relaxation, are larger and sharper in the periodic pipe than in the closed one, which contributes, in addition to the lower retrograde velocity amplitude, to explain the smaller velocity ratios reported in the periodic tubes. The intervals of the average pressure (first diagram) are very similar across the ten segments in the periodic pipes, but both limit values of these intervals increase monotonically in the closed pipe; in the central segment s . 6 (third diagram), the pressure minimum p m i n is about three times smaller in the closed tube than in the periodic tube, in agreement with its larger retrograde velocity amplitude during contraction, which should also contribute to a higher value of mass change during the motility cycle. However, the largest pressure intervals, both obtained in the case L O V I S VMAX, reduce their disparity between closed and periodic pipes according to Table 5. Overall, these two figures clearly illustrate that there is a notable difference in the antegrade and retrograde velocity distributions between the closed and periodic circular tubes when they become highly obstructed.

3.5.6. Motility Pattern M O T - I : Summary

Motility pattern M O T - I from Table 1 is combined with the five occlusion scenarios O C L - 1 ÷ O C L - 5 from Table 4; the first four scenarios are run in the two closed geometries ( P - C L O and C - C L O ) and in the two periodic geometries ( P - P E R and C - P E R ); O C L - 5 with o f > o max is restricted to the two circular configurations, C - C L O and C - P E R . The amplitude of the relative volume change max ( δ V ( t ) = V ( t ) V ( o = o n ) 1 ) varies from 6 % ( O C L - 1 ) to 16 % ( O C L - 5 ), while the relative mass variation reaches ± 5 % for the entire colonic motility cycle in the “parabolic” pipes due to the largest pressure oscillation in the case P - C L O H I V I S S L O W performing the most challenging protocols O C L - 3 & O C L - 4 . L B M - D B experienced no stability issues in these simulations using c s 2 = 1 3 , where the compressibility effect is not significant.
Giving motility protocol and geometry, the computational series consists from the six runs combining the two kinematic viscosities ν p h ( L O V I S & H I V I S ) with the three propagation speeds v w ( S L O W , F A S T & VMAX). Table A2, Table A3, Table A4, Table A5 and Table A6 report the pressure and velocity distributions in central segment 6 ( s . 6 ) for the five motility protocols, respectively; these results are resumed in Table 5. Figure 10, Figure 11, Figure 12, Figure 13, Figure 14, Figure 15, Figure 16, Figure 17, Figure 18, Figure 19, Figure 20 and Figure 21 illustrate these data, but also display the dynamics of segment-averaged pressure and longitudinal velocity distributions across the ten segments.
As a rule, the distribution shapes appear qualitatively similar within the same series, but smaller values of ν p h and higher values of v w make them more oscillating.
It is observed that the segment-averaged velocity intervals [ v ¯ m i n , v ¯ m a x ] and the segment-peak velocity intervals [ v m i n , v m a x ] scale nearly linearly with v w ; the ratios of their highest retrograde and antegrade amplitudes, respectively, R ( v ¯ ) = | v ¯ m i n | v ¯ m a x and R ( v ) = | v m i n | v m a x , then depend weakly on ν p h and v w . This property should make it possible to (approximately) predict the velocity intervals, and then the disparity between the antegrade and retrograde velocity amplitudes from a single combination ( ν p h , v w ) .
The largest retrograde velocity amplitudes, | v m i n | = | min t ( v m i n ( t ) ) | and | v ¯ m i n | = | min t ( v ¯ ( t ) ) | , commonly occur either at the highest local occlusion of a given segment or very shortly after it, when the downstream neighbor begins its contraction; these results are in qualitative agreement with the MRI estimates from Table 2 and Table 3.
Generally, Re ( | v ¯ m i n | ) increases progressively from the first (least occluded) scenario O C L - 1 to the last (most occluded) scenario O C L - 5 , providing the same fluid properties, geometry, propagating wave speed v w , and sound velocity c s . The greatest value of Reynolds number max ( Re ( | v ¯ m i n | ) ) is commonly obtained in the fastest low-viscous case L O V I S VMAX, where max ( R e ( | v ¯ m i n | ) ) [ 4.9 , 81.6 ] increases significantly from P - P E R O C L - 1 to C - C L O O C L - 4 at the same value c s 2 = 1 3 . Among the four first scenarios, O C L - 2 and O C L - 4 operate with higher relaxation and occlusion speeds than O C L - 1 and O C L - 3 due to the smaller values of o r and larger values o f o r in Table 4, the O C L - 2 & O C L - 4 then generally produce higher values of Re ( | v ¯ m i n | ) and larger pressure intervals than two of their counterparts, respectively. At the same time, O C L - 1 and O C L - 3 display significantly higher velocity contrasts than O C L - 2 and O C L - 4 , respectively, due to their relatively larger retrograde velocity amplitude; their values further increase to R ( v ¯ ) 5.2 ( L O V I S S L O W ) and R ( v ) 4.3 ( L O V I S VMAX) with the most obstructive scenario O C L - 5 in circular closed pipe ( C - C L O , Table A6).
Comparing the periodic and closed ends, the non-peripheral segments of the closed pipes produce antegrade velocity peaks due to the strong pressure oscillation during motility cycle in the last segment ( s . 10 ); in contrast, its upstream neighbor s . 9 displays the largest retrograde velocity amplitude. These extreme events, which are due to the last contraction, are almost absent in the periodic pipes where, on the other hand, the second retrograde peak amplitude occurs when the downstream neighbor of a given segment contracts to the maximum occlusion degree. The closed pipes produce larger averaged and peak retrograde velocity amplitudes, and therefore Re ( | v ¯ m i n | ) , but also larger pressure intervals than their periodic counterparts, with the three last occlusion scenarios O C L - 3 - O C L - 5 .
Concerning the dependence on the transverse shape, the “parabolic” (“P”) and circular (“C”) pipes, both closed and periodic, produce a similar range of velocity and pressure with the two least occluding scenarios, as O C L - 1 and O C L - 2 . However, O C L - 3 and O C L - 4 increase the degree of contraction to o f = o max , where three parabolic contours touch and form acute angles; in this non-smooth transverse form, the pressure peaks of p m a x ( t ) are located very shortly before local contraction o ( s . 6 ) = o f , in qualitative agreement with Figure 5 [9] (see Figure 16, Figure 17 and Figure 19). The pressure peaks become an order of magnitude larger in “parabolic” pipes, both closed and periodic, not only compared to less obstructive scenarios, but even relative to their circular counterparts, which contract and relax smooth bounding surfaces without fixed points (cf. Figure 18 and Figure 19). The largest observed pressure difference p m a x p m i n 2.2 mmHg then becomes reached in the case O C L - 4 H I V I S VMAX, in both closed and periodic “parabolic” tubes, and this pressure gap exceeds its obtained values in circular tubes, even when they deform under the heaviest obstruction scenario O C L - 5 . Consistently, O C L - 1 and O C L - 2 reach a slightly larger range of Re ( | v ¯ m i n | ) in circular pipes, while O C L - 3 and O C L - 4 replace them by “parabolic” tubes.
Recall that L B M - D B with the motility model M O T - I is believed to follow DCM dynamics in the MRI experiment [68], where the two numerical scenarios O C L - 1 ( o n 0.367 ) and O C L - 3 ( o n 0.539 ) both apply the relative MRI interval [ o r , n , o f , n ] = [ 0.2 , 0.6 ] . We observe that the experimental (MRI) and numerical ( L B M - D B ) mean velocity profiles display, in a given segment s . 6 , and particularly in closed pipes, a net global retrograde peak min t ( v ¯ ) , and it occurs either at the time of maximum occlusion or very shortly after this event. Quantitatively, O C L - 1 in Table A2 gives lower velocity contrasts and slightly lower values of Re ( | v ¯ m i n | ) than those of our estimate of the MRI graphical data in Table 3. Herein, L B M - D B offers max ( R ( v ¯ ) ) [ 1.39 , 2.12 ] on the four geometries in the case H I V I S S L O W , and thus LBM-DBt covers the H V O L MRI result R ( v ¯ ) = 2 , but the L B M - D B ratio of peak values R ( v ) [ 1.46 , 1.69 ] is smaller than our MRI estimate R ( v ) = | v m i n | v m a x = 3.25 . On the other hand, closed circular pipes, deforming under the O C L - 3 scenario, produce { S L O W , F A S T } peak ratios as R ( v ) = { 3.13 , 3.19 } ( H I V I S ) and R ( v ) = { 3.32 , 3.23 } ( L O V I S ). These results formally agree best with the MRI data, as R ( v ) = | v m i n | v m a x { 3.25 , 3.75 } ( H I V I S ) and R ( v ) = | v m i n | v m a x = { 3 , 3 } ( L O V I S ); however, it should be kept in mind that the L B M - D B peak values are the extreme values on a given segment, while the MRI peak values are measured on the tube axis. At the same time, the O C L - 3 offers Re ( | v ¯ m i n | ) = { 1.83 , 3.79 } ( H I V I S ) and Re ( | v ¯ m i n | ) = { 6.62 , 13.1 } ( L O V I S ) then exceed their MRI H V O L estimates of Re ( | v ¯ m i n | ) = { 0.48 , 1.92 } ( H I V I S ) and Re ( | v ¯ m i n | ) = { 1.8 , 4.9 } ( L O V I S ), respectively.
To conclude, we note that the amplitude of the MRI retrograde velocity min t | v ¯ m i n | in Table 3 decreases from low fluid volume ( L V O L , 60 % ) to high fluid volume ( H V O L , 80 % ), and we expect the same trend between H V O L and the full fluid, which is here modeled numerically. The least occluding motility protocol O C L - 1 is then the best coming close to the MRI data, but O C L - 1 is perhaps slightly less occluding than the motility scenario of an experimental device, while O C L - 3 , operating with the same relative interval, is far too obstructive.

3.6. MOT-II vs. MOT-I

Motility pattern M O T - I I from Table 1 is depicted in Figure 4 left. In an effort to mimic the DCM motility pattern [9], M O T - I I first reduces relaxation time from t r = t f to t r = 0.375 t f ; second, M O T - I I marks pause f after the forward step during the time of t f , p = t f ; third, M O T - I I obeys “an intestine law” such that a given segment begins to relax back to the neutral state when its downstream neighbor becomes most obstructed; fourth, M O T - I I imitates the larger elasticity effect [9] and then it returns back about 1.5 times longer than M O T - I ; overall, the M O T - I I motility cycle lasts about 15 % longer than that of M O T - I .
Motility pattern M O T - I I is now combined with the five occlusion configurations from Table 4; the results are reported in Appendix C.1, Figure A6, Figure A7, Figure A8, Figure A9, Figure A10, Figure A11, Figure A12 and Figure A13 and Table A7, Table A8, Table A9, Table A10 and Table A11; these results are discussed in this section. Table 5 summarizes the data; it indicates that the relative volume change δ V ( t ) of M O T - I I , from 12.5 % ( O C L - 1 ) to 30.12 % ( O C L - 5 ) per motility cycle, almost doubles in magnitude that of M O T - I .
Overall, M O T - I and M O T - I I show similar results, but they also manifest several distinctive features. Typically, the M O T - I I peak velocity profiles appear sharper; in particular, they manifest a sequence composed of two retrograde peaks v ¯ m i n ( t ) occurring during each occlusion event: the first begins at the time of pause f , and it is preceded by a large antegrade peak v ¯ m a x ( t ) of similar amplitude, while the second retrograde peak occurs towards the end of pause f , when the downstream neighbor becomes obstructed and rapidly pushes the flow backward, accompanied by the simultaneous pressure drop and a decrease in the magnitude of the antegrade peak upstream of the contracted haustra segment.
Like M O T - I , M O T - I I progressively increases from O C L - 1 to O C L - 5 , an amplitude of the segment-average retrograde velocity | min t v ¯ m i n ( t ) | ; typically, R e ( | v ¯ m i n | ) = | v ¯ m i n | l ν p h then grows by a factor of about 1.5–2.5 from M O T - I to M O T - I I in the same configuration. It could be suggested that an increase in reflux amplitude is due to more rapid relaxation and greater reduction in fluid volume during the motility cycle. So far, the case L O V I S V M A X reaches R e ( | v ¯ m i n | ) 90 in closed circular pipe under the most obstructing scenario O C L - 5 and, according to analysis in Appendix C.1, this value even underestimates the effective Reynolds number in an incompressible fluid due to the high numerical compressibility level of this fastest example (see below).
Similar to M O T - I , M O T - I I produces a higher retrograde amplitude | min t v ¯ m i n | (i) in closed pipes compared to periodic pipes and (ii), in circular pipes compared to “parabolic” pipes under O C L - 3 O C L - 5 . In closed pipes, M O T - I I increases the largest M O T - I ’s contrasts between retrograde and antegrade amplitudes, and O C L - 5 then reaches the highest values, as max ( R ( v ¯ ) ) 7 and max ( R ( v ) ) 4 in the case H I V I S S L O W (see line above the last in Table 5), where an average flow becomes retrograde from local contraction to final contraction due to the closure of the terminal segment, which then exhibits a high pressure peak above an averaged pressure level (see M O T - I I O C L - 5 , C - C L O in Figure A11).
At the same time, M O T - I I O C L - 5 provides the only example where the maximum of the antegrade mean velocity amplitude max t ( v ¯ ) exceeds its retrograde counterpart min t ( v ¯ ) over the entire series of six combinations ν p h v w ; this occurs in circular periodic pipes, where the largest contrasts max ( R ( v ¯ ) ) = 0.95 and max ( R ( v ) ) = 0.73 , obtained in the case L O V I S VMAX, are both lesser than one (see last line in Table 5 and C-PER in Figure A13).
The second motility pattern increases the magnitude of the relative mass change per motility cycle | max ( δ M ) | , which reaches | 17 % | in the stable case of c s 2 = 1 3 in the case H I V I S S L O W in the closed “parabolic” pipe (see M O T - I I O C L - 3 P - C L O in Table 5). An amplitude of | δ M | generally remains much smaller in the circular pipes on the one hand and in the periodic pipes on the other hand, with the same configuration.
Concerning pressure amplitude, the in vitro DCM in Figure 5 [9] and the in silico L B M - D B show the largest pressure oscillations just before and at the time of local contraction. “Parabolic” pipes produce more irregular flow patterns and larger pressure intervals than circular pipes with O C L - 3 and O C L - 4 , where M O T - I I O C L - 3 increases pressure oscillation by a factor of about three with respect to M O T - I and reaches p m a x p m i n 5 ÷ 6.2 mmHg in the stable case H I V I S V M A X of c s 2 = 1 3 , while the magnitude of the pressure drop remains 5–10 times lower in circular tubes.

3.7. LBM-DB: Stability vs. Compressibility

Our simulations of the segmented deformation are performed using the (third) reconstruction algorithm r - bc - corner from Section 2.2.4, which computes all incoming populations with boundary rules and thus avoids the equilibrium reconstruction [78] for corner links. The combination of the p - pressure B F L 4 rule in corners with the p - pressure + p - flow B F L 5 scheme otherwise is employed in all cases, but B F L 4 is switched to B F L 0 at the first iteration in “fresh” nodes, dropping its local post-collision correction K 4 n ^ q in Equation (11); n ^ q is then restored with the three local collision sub-iterations, without any streaming neither towards new-born nodes nor between them. L B M - D B then retains quite satisfactory robustness in piece-wise continuous motility protocols with respect to very small lumen aperture, fast occlusion speed, and rapid wave propagation, at least up to R e ( | v ¯ m i n | ) 82 ( M O T - I I O C L - 4 , L O V I S V M A X in Table A10) and R e 49 ( M O T - I I O C L - 5 , L O V I S F A S T  Table A11); these stable results are obtained using c s 2 = 1 3 in closed circular pipes.
In detail, L B M - D B poses no stability problems with M O T - I , but it begins to lose its robustness under M O T - I I in “parabolic” pipes subjected to occlusion scenarios O C L - 3 & O C L - 4 in the fastest flow regime L O V I S VMAX, most likely due to a rapid increase in fluid velocity and the large pressure oscillations in a very small opening bounded by overly curved surfaces with acute angles. The reported stable results are then obtained by increasing the Mach number by factors of 2 2 ( O C L - 3 ) and 2 3 ( O C L - 4 ) due to reduction in sound speed c s 2 in Equation (4a) (see Table A9 and Table A10). The circular pipes require a similar reduction of c s only in the most difficult case O C L - 5   L O V I S V M A X and only with closed ends (see C - C L O in Table A11).
Appendix C.1.3 examines the dependence of the results on c s in stable cases and concludes that even if the numerical compressibility effect has a relatively small impact on the shape of velocity distribution, it significantly increases the amplitude of the relative mass change | δ M | and reduces the magnitude of pressure and velocity fields. So, taking the example of O C L - 4   L O V I S VMAX, which now operates with a small value of c s in “parabolic” pipes, one would expect its pressure interval to exceed, in incompressible flow, the H I V I S   V M A X result of p m a x p m i n 5 ÷ 6.2 mmHg mentioned above.

3.8. LBM-DB: Numerical Performance

Finally, ParPac’s code structure and parallelization aspects are discussed in Appendix B.1. Regarding numerical efficiency, the ParPac code performed all colonic computations using the M R T collision operator on the DELL Laptop with intel CORE i7 10TH GEN using 4 CPU (ou cores, mpirun -np 4). Computational time in the typical case “C”: M O T - I : O C L - 1 lasts ≈ {1000–1310, 654, 444, 507} s on { 1 , 2 , 4 , 8 } CPUs (or cores), respectively, with 24,960 fluid grid nodes at the neutral occlusion for 21,760 temporal iterations. This number of steps is the same for all computations using the same (characteristic velocity) value of peristaltic speed V w 8 × 10 3 [ l . u ] . A processor count of 2–4 is then optimal at a given grid resolution, providing appropriate scalability without any specific optimization performed for this research study.

4. Conclusions

L B M - D B proposes to improve the boundary accuracy and quality of reconstruction in “fresh” fluid nodes using the most accurate and compact boundary schemes, which respect both the gradients of surrounding pressure and the curvature of flow. In this work, L B M - D B is applied to model the Navier–Stokes Newtonian fluid flow in a deforming channel and in a colonic tube, when its delimiting membrane deforms according to the prescribed intestinal law. The two intestinal motility patterns from Table 1, M R I (following [68]) and M O T - I I (following [9]), are simulated combining each of them with the five occlusion scenarios O C L - 1 ÷ O C L - 5 from Table 4 of increasing occlusion degree. These ten motility models are applied in the four geometries composed of closed and periodic colonic tubes of circular (“C”) and “parabolic” (“P”) apertures by combining the two viscosity values [68] with the three realistic propagation speeds [9,68]. In that, the dynamic “P”-opening aims to imitate the real triangular transverse shape [53] and its DCM representation [9]. We note that our additional results with intermediate neutral buffers of length l n = 0.2 cm, imitating the semilunar folds [9] located between the haustra units, do not appear to show any significant difference from those reported in this work. Next, the periodic and closed entry/exit boundary conditions are quite artificial in the colonic context but allow us to get an idea on the role of the inlet/outlet closure for all the scenarious modeled in the two geometries. Future work could expand on the physiological conditions directly informed by an experiment, e.g., extending the Dirichlet pressure L B M schemes [44,48,92].
The temporal evolution of pressure during the motility cycle is found in a qualitative agreement with the velocity profiles; the sequential sub-division of the pressure distributions into positive and negative branches with respect to the mean pressure of the colon appears to be in qualitative agreement with experimental and computational data [61]. It is worthwhile to notice that circular pipes which operate with a much smoother surface without any fixed points produce pressure peaks one order of magnitude smaller than the “parabolic” pipes at the high degree of occlusion. However, the largest relative pressure amplitude observed in the high-viscosity fastest flow ( H I V I S F A S T ) reaches the level of only 6 mmHg (with M O T - I I O C L - 4 in the “P” —pipe), which is below the manometric DCM measurement [9] of 25 mmHg monitored in a partially filled, non-Newtonian fluid flow subjected to active membrane suction [9].
On the other hand, in agreement with observation [62] “that the mixed flow of two similar amplitudes happens in an ascending colon”, our results show an almost simultaneous presence of both antegrade and retrograde velocity peaks in a given segment, occurring, respectively, shortly before and very soon after the moment of local contraction (and therefore downstream relaxation). This positioning of the segmented retrograde peaks agrees with the analysis of Table 3 on MRI velocity distributions [68]. In turn, the location of antegrade peaks is in agreement with the MRI velocity curve [64] monitored at the center of the contracting haustra unit.
Additionally, we observe that M O T - I I typically produces the multiple retrograde peaks of a similar amplitude during its pause interval after the highest contraction. The three retrograde velocity patterns around the moments of (i) the upstream, local, and downstream contractions are then almost the same, but they significantly decay in amplitude before and after these major occlusion events. The velocity profiles become particularly irregular in the “parabolic” highly occluded pipes, where the difference from their circular counterparts becomes clearly visible, especially in the closed case. In turn, the periodic pipes generally manifest the significant antegrade velocity peaks during the time of the upstream and local contractions, and their amplitude even exceeds the retrograde one for all six ν p h v w combinations under the strongest deformation M O T - I I O C L - 5 . Remarkably, we observe in all scenarios that the largest mean and peak velocity amplitudes, both antegrade and retrograde, scale almost linearly with peristaltic wave speed v w . Therefore, their respective ratios R ( v ¯ ) = | v ¯ m i n | v ¯ m a x , and especially R ( v ) = | v m i n | v m a x , are almost independent of v w and decay slightly as viscosity decreases.
In more detail, the two first scenarios, O C L - 1 and O C L - 2 , operate with the occlusion degrees in the range of the DCM experiments, as o f { 0.6 , 0.7 } . The observed velocity contrasts are then similar both in the closed and periodic pipes on the one hand and in the two transverse forms on the other hand. The largest velocity contrasts manifest themselves at the high-viscosity slowest flow ( H I V I S S L O W ) as R ( v ¯ ) 2 ( 1.7 ) and R ( v ) 1.7 ( 1.9 ) with M O T - I ( M O T - I I , respectively). However, these two characteristic values increase to R ( v ¯ ) 5.2 (7) and R ( v ) 2 (4, respectively) with the next three scenarios O C L - 3 O C L - 5 , where o f [ 0.86 , 0.94 ] approaches complete obstruction. In turn, the Reynolds number R e , which is based on the largest retrograde mean velocity amplitude in the central segment, regularly increases from O C L - 1 to O C L - 5 and from M O T - I to M O T - I I . Commonly, the low-viscosity fastest fluid flow ( L O V I S VMAX) shows the largest Reynolds range, as R e 82 with M O T - I O C L - 4 in the closed circular pipe, and this limit should extend above R e 90 with M O T - I I O C L - 5 at the same compressibility level.
Overall, L B M - D B demonstrates physically sound results over a surprisingly wide parameter range, well above experimental and numerical data [68], by performing computations with sound speed value c s = 1 / 3 [ l . u ] , where the Mach number is small, M a = V w / c s = 1.35 × 10 2 . L B M - D B then only loses its stability in the low-viscosity fastest flows with the most difficult M O T - I I protocols, as O C L - 3 O C L - 4 in “P”-pipes and O C L - 5 in the closed “C”-tube, most likely due to a large local velocity increase above V w . The remedy is found in a strong reduction in c s , which still demonstrates a qualitatively correct velocity response to the sequence of occlusion events. However, the pressure and velocity magnitudes then become lower than expected from an incompressible fluid. Further modeling on finer grids should increase the numerical values of kinematic viscosity and then enable stable computations of an incompressible fluid even with the fastest and most occlusive motility models in the “parabolic” pipes, which fit better the DCM geometry and produce more realistic pressure waves amplitudes than commonly employed cylindrical tubes.
Concerning enhancement of mixing and adsorption versus motility, our observations suggest that the chyme flow is composed of the segmented vortexes of variable amplitude, depending not only on the peristaltic wave speed, occlusion degree, and occluding speed, but also, in the event of strong contractions, on the entry/exit boundary conditions. It would be interesting to quantify the propulsive and retrograde contributions to mixing, e.g., by comparing the modeling of tracer dispersion following [4,5,8,27,28,29] with the estimates of either the level of stretching [24,25] or the degree of mixing efficiency [64,103].
From a hydrodynamic point of view, several studies seek to predict shear rates along and across a digestive flow for medical and pharmaceutical purposes, e.g., [11,15,68]; on the other hand, an adequate gut modeling requires to represent it by the non-Newtonian fluid. The LBM is advantageous in both cases, because it computes the shear stress components locally, by projecting the non-equilibrium onto the appropriate quadratic polynomials of the velocity components, e.g., [42,43,78,104,105]. Therefore, it is expected that the LBM could handle the highly viscous flow and the rheology of non-Newtonian fluids more easily and accurately than S P H [106]. Moreover, it is very likely that the numerical viscosity value is automatically increased locally with digestive rheology. For this reason, we do not attempt to stabilize some fastest and most occluding simulations found within the stability limits in this work. A natural extension of our work would then lead to the generalization of all L B M - D B steps for the space- and time-dependent viscosity distribution, with a subsequent extension for the LBM description of non-Newtonian rheology, such as the Bingham fluid law [107,108,109,110], the Ostwald de Waele “power law” relationship [111], the Carreau–Yasuda viscosity model [112,113], or their extensions to the Herschell–Bulkley constitutive model [114], which are mainly employed to represent intestinal flow, chyme, and soft feces, e.g., [3,5,15,24,62,115].
The L B M - D B does not require the conservation of the lumen volume; so far, its largest relative reduction reaches 31 % with the last motility protocol M O T - I I O C L - 5 . This numerical property can help in handling of physical situations in which lumen volume and fluid mass vary through a modeled colonic slice, e.g., due to the presence of air, a non-uniform water adsorption, an irregular mass transfer or a numerical modeling of the visco-elasticicity of the membrane. A subsequent involving stage would consist of replacing an heuristic “elastic relaxation motility” with an appropriate physical membrane description adapting the existing LBM models in elastic solids, e.g., [116,117,118,119,120].
Regarding better understanding of physiology and pathology, one may first think about an extension from the single peristaltic wave to the colinear and/or colliding waves following [28]. Also, L B M - D B can break easily the symmetry between the three “parabolic” sectors and attempt to mimic intrahaustral motility due to chaotic ripple activity [51]. Furthermore, while DCM [9,64], D C M D T [68] and currently L B M - D B neglect the longitudinal haustral contractions, which could be triggered by distension of the lumen due to intraluminal accumulation of colonic content, it should become possible to incorporate them into L B M - D B due to its intrinsic ability for mass and volume change. Finally, L B M - D B can be extended to process the 3D anatomical images of the colon without any major algorithmic changes.

Funding

This work is not supported by any funding.

Data Availability Statement

The data presented in this study are available on request from the corresponding author due to INRAE confidentiality policy.

Acknowledgments

The author is grateful to P. Klein (ITWM) for common work on the ParPac code between 1999 and 2002, to C. de Loubens (CNRS) for introducing me to the field of “gut modeling”, to A. Vikhansky (Siemens) for his suggestion of the deformable channel benchmark, to F. Deslandes (INRAE) for his advice on the update to modern OpenMPI, to V. Fromion (INRAE) for encouraging discussions, and to the Documentation Center (INRAE, Jouy-En-Josas) for its kind hospitality enabling this work to be carried out.

Conflicts of Interest

The author declares no conflict of interest.

Appendix A. Circular and “Parabolic” Membranes

This section completes Section 2.1. The two contours delimiting the 2D lumen cross-section of the colonic pipe, “parabolic” (“P”) in Figure 5, Figure 6, Figure 7 and Figure 8 and circular (“C”) in Figure 9, are dynamically constructed from a given occlusion degree o ( t ) . “C” defines the radius value r of the circular aperture providing a current open area value π r 2 s ( t ) with s ( t ) = 1 o ( t ) . “P” subdivides an open surface π r 2 s ( t ) into three sectors i = { a , b , c } , each delimited by the parabolic contour p ( i ) and the two linear segments which connect the center point O = ( 0 , 0 ) of an inscribed triangle to the two end vertices of p ( i ) ; e.g, if the triangle’s side “a” is opposite vertex A ( 0 , r ) , the sector’s area S ( a ) is bounded by parabola p ( a ) and two linear segments [ O , B ( 3 2 r , r 2 ) ] and [ O , C ( 3 2 r , r 2 ) ] . The three parabolic curves y = n x 2 + l are specified in Table A1 according to the sector’s area S ( i ) ; they are expressed in coordinate axes x = x cos [ ϕ ] + y sin [ ϕ ] , y = y cos [ ϕ ] x sin [ ϕ ] , which are, respectively, parallel and perpendicular to side i = { a , b , c } of the equilateral triangle inscribed in an outer cylinder of radius r. In this work, S ( i ) , then n ( i ) and l ( n ) , vary linearly in time together with o ( t ) and s ( t ) ; accordingly, r 2 varies linearly in cross-section “C”.
Table A1. Three parabolic contours y = n x 2 + l giving an outer cylinder radius r, r 2 = 9 r 2 , r 3 = 18 r 3 , occlusion degree o ( t ) and S ( i ) = π r 2 ( 1 - o ( t ) ) 3 when the three sectors have the same area S ( i ) .
Table A1. Three parabolic contours y = n x 2 + l giving an outer cylinder radius r, r 2 = 9 r 2 , r 3 = 18 r 3 , occlusion degree o ( t ) and S ( i ) = π r 2 ( 1 - o ( t ) ) 3 when the three sectors have the same area S ( i ) .
p(i) cos [ ϕ ] sin [ ϕ ] n ( S ) l ( n )
p ( a ) 10 12 3 S ( a ) r 2 r 3 r 4 ( 2 + 3 r n ( a ) )
p ( b ) 1 2 3 2 12 3 S ( b ) r 2 r 3 r 4 ( 2 3 r n ( b ) )
p ( c ) 1 2 3 2 12 3 S ( c ) r 2 r 3 r 4 ( 2 3 r n ( c ) )

Appendix B. Addition on LBM-DB in Section 2.2

Appendix B.1 extends Section 2.2.1 with respect to organization of the ParPac code. Appendix B.2 complements Section 2.2.3 for classification and implementation of the boundary rules. Appendix B.3 provides more details about Section 2.2.4. Appendix B.4 discusses numerical simulations in deformable channel flow.

Appendix B.1. The ParPac Code: Structure and Parallelization

The Parallel Particles C++ code (ParPac [93]) was developed during the period of 1999–2002. with d 2 Q 9 and d 3 Q 15 discrete velocity sets for modeling and vizualization of flow in complex structures. ParPac is based on the completely object-oriented design, where the mathematical optimization methods provide an optimal domain decomposition within the framework of the OpenMPI library (initially m p i 1.0 . 7 m p i c h 2 ). ParPac allows an almost automatic parallelization of a given application using its data exchange libraries, minimizing communication mechanisms. ParPac [93] was first applied to the modeling of single- and two-phase (air/oil) flows in generated fiber filters, but also for flow around filaments in melt spinning processes. Additionally, ParPac is equipped with a dynamic adaptation to ensure a well-balanced use of computing power in free-interface flow models. They were specified for metal casting in realistic 3D shapes [121] and for an academic study of Bingham fluid filling in expanding 2D cavities [107]; the free interface module was coupled with other CFD codes and applied for melt solidification. Later, ParPac was extended for the L B M solving of linear and non-linear, anisotropic advection-diffusion equations (see [98] and references herein). In the present work, ParPac is updated to the recent version of OpenMPI ( m p i c h 4.0 ) and extended for advanced boundary schemes [44,79]. On this basis, ParPac’s “Moving boundary module” [78] integrates novel reconstruction algorithms, particularly in corner geometries (Appendix B.2 and Appendix B.3). They are validated in modeling the deformable channel (Appendix B.4) and circular colonic contractions (Section 3 and Appendix C).
ParPac’s capability for both research and application purposes comes from the use of specific “point-class” lists, with the “point-class” containing all the necessary grid-collocated information, such as the local population solution and their moments. In detail, in the “Moving boundary module”, the solid surface is sub-divided into two sub-domains: (a) the static-solid  V s , which is permanent, and (b) the moving-solid V m ( t ) which moves and/or deforms but must not shift its contour further than towards the neighboring cell in a single time step. Accordingly, fluid grid nodes r V f ( t ) are sub-divided into three fluid lists:
1.
BulkList : r V l ( t ) has no solid directional neighbors;
2.
StaticBoundaryList : r b V b s ( t ) V b ( t ) has only static-solid directional neighbors;
3.
MovingBoundaryList : r b V b m ( t ) = V b ( t ) V b s ( t ) has at least one moving-solid neighbor.
These lists are dynamically updated according to the status of surrounding solid nodes. There is no need to allocate static-solid grid nodes V s . However, V m ( t ) grid nodes are allocated for simplicity, and they are sub-divided into two lists:
1.
SolidList : r V m s ( t ) has no fluid neighbors.
2.
SolidBoundaryList : r V m b ( t ) = V m ( t ) V m s ( t ) has at least one fluid neighbor.
Once the solid-fluid interface is updated, “fresh” nodes r n ( t + 1 ) R ( t + 1 ) are defined as those which move from V m b ( t ) to V b m ( t + 1 ) . ParPac triggers specific lists, as, for example, the three fluid lists that perform collision and streaming in c - s - s t e p (see Section 2.2.1). After that, the two fluid boundary lists construct their solution for the opposite links to the cut links via the pre-selected boundary rules in bc - s t e p (see Section 2.2.3); these rules may depend upon a number of available upstream fluid neighbors. Finally, the fresh nodes list reconstructs populations via the (re-)filling algorithms in f - s t e p (see Section 2.2.4). The neighborhood information required by a given triggered function must be pre-communicated via pre-defined buffers, either by a triggered list or a larger “parent” list, and then used for multiple steps. In order to decrease communication data, the propagation step is performed before bc - s t e p , allowing to implement the two-node M R boundary rules from Appendix B.2 through the one-neighbor communication.

Appendix B.2. Classification of the Boundary Rules

The boundary rule is said to be parametrized [44,45,48,79] when it shares the bulk parametrization property of the T R T operator, in our case here, when an incompressible steady-state solution of the momentum obeying N S E is fixed by Reynolds number R e and collision number Λ on a given grid. This property is not automatic, but it is ensured by B B and the following three principal, infinite-member classes of directional schemes [44,79,94]:
1.
The single-node class [79]: L I + = { L I 1 + , L I 3 + , L I 4 + } yields γ = γ ^ = 0 in Equation (10), with typically I 1 = 1 or β = 0 there, while no - nb links apply I 1 = 0 . L I + provides an exact linear (Couette) velocity profile in an arbitrarily grid-inclined channel using the Stokes equilibrium ( I n s e = 0 in Equation (4a)); this property is abbreviated as c - flow accuracy; this property extends to N S E equilibrium ( I n s e = 1 in Equation (4a)) in a straight channel. Three coefficients { α ^ , β , β ^ } are derived from the directional linearly accurate Taylor approximation of the Dirichlet velocity condition, and they are expressed by Equation (46) [79] via the free adjustable coefficients α q ( p ) and α q ( u ) from Equation (13).
(a)
Class L I = { L I k } sets α q ( p ) = 0 and defines its coefficients with respect to α q ( u ) ] 0 , δ q 1 ] with Equations (53)–(54) [79]. Additionally, class L I prescribes K k + = 0 in Equation (11), while its coefficient K k of the anti-symmetric post-collision correction is distinctive between the three sub-classes L I k , here k = { 1 , 3 , 4 } : (i) K 1 guarantees the parametrization; (ii) additionally, K 3 makes the Poiseuille Stokes flow exact in an arbitrarily inclined channel, abbreviated as p - flow ; (iii) K 4 holds parametrization and makes the Dirichlet velocity closure relation independent of the pressure gradient, abbreviated as p - pressure .
(b)
The local single-node class [79] E L I + = { E L I k + } is built with ideas from [94]; its coefficients are defined by Equations (53) and (54) [79]; E L I + does not need to have a fluid upstream neighbor because it sets β = 0 . However, in addition to L I k , E L I k + must prescribe the symmetric boundary value α q ( p ) e q + ( r b , t ) in Equation (12), as well as the symmetric post-collision terms K E L I + n ^ q + ( r b , t ) in Equation (11) with implementation [79]; the E L I + implementation in “fresh” nodes is postponed to future work because of these two additional approximations, which can be computed in principle with an initial estimate and then updated with the r - iter .
2.
The two-node class M R [48,78,79] ensures the combined accuracy p - p = p - flow p - pressure due to its five non-zero coefficients in Equation (10); it is most suitable with the two - nb   link offering I 1 = I 2 = 1 , while one - nb   link can apply it with I 2 = 1 , I 2 = 0 . Class A V M R , which is a recent extension [48,79] of M R , additionally makes the second-order pressure gradient q 2 e q + vanish from the Dirichlet velocity closure condition, abbreviated as p - p - p , and it demonstrates the most accurate N S E solutions [79]. A V M R operates with α q ( p ) = K + = 0 , and its coefficients are given in Table VI [79].
3.
Class L I 5 updates class M L I , originally described with the F ^ q correction given by Equation (5.6) [44] and Equations (50), (53) [81]. In the present context, L I 5 simply complements the L I 4 term of K 4 n ^ q in Equation (11) with parabolic correction 1 2 α ( u ) δ q 2 q 2 e q ( r b , t ^ ) ; L I 5 then satisfies both the p - pressure and p - flow criteria and therefore matches their combined accuracy of the p - p condition. Here, a three-point asymmetric finite-difference directional approximation of q 2 e q ( t ^ ) is computed with (1) the downstream Dirichlet value e q ( r b + δ q c q , t ^ ) , (2) the local solution e q ( r b , t ^ ) , and (3) the upstream solution e q ( r b c q , t ^ ) (see Equation (5.14) [44], Equation (53) [81]). In principle, L I 5 applies in corner links using the Dirichlet velocity values of their two directional solid neighbors [44,92], while a new node r n can compute e q ( r n ) with an extrapolated momentum value of ρ u ( r n ) on the first sub-step of the r - iter update, then with an effective momentum value in subsequent iterations.
Let us take an example of these boundary classes using the B F L coefficients: the original B F L 0 [100], B F L 4 (Equation (52) [79] with Table III) and B F L 5 (class M L I [44], Equation (53) [81]); they all share the same coefficients α ^ γ in Equation (10), then the same values α ( p ) and α ( u ) in Equations (13)) and the same Dirichlet term in Equation (12), but differ in their post-collision correction F ^ q in Equation (11). Giving a positive adjustable collision rate Λ = τ 1 2 from Equation (3), the four B F L k boundary classes read as
B F L k : γ = γ ^ = α ( p ) = 0 , δ [ 0 , 1 2 ] : α ^ = 2 δ , β = 1 α ^ , β ^ = 0 , α ( u ) = 2 , δ [ 1 2 , 1 ] : α ^ = 1 2 δ , β = 0 , β ^ = 1 α ^ , α ( u ) = 1 δ ; B F L 0 : F ^ q = 0 ; B F L 4 : F ^ q = K 4 n ^ q , K 4 = 2 + α ( u ) ( Λ 1 2 δ ) , B F L 5 : F ^ q = K 4 n ^ q + 1 2 α ( u ) δ q 2 q 2 e q ( r b ) , q 2 e q ( r b ) = 2 δ q + δ q ( e q ( r b + δ q c q ) e q ( r b ) δ q e q ( r b ) e q ( r b + δ q c q ) δ q ) .
Recall that δ q is the distance between boundary node r b and the solid surface at r b + δ q c q along the cut link c q , while in the opposite direction c q , we define δ q = 1 if r b + c q is a fluid node, otherwise r b is a corner node and δ q [ 0 , 1 ] is the distance to solid wall along an opposite cut link c q .

Appendix B.3. Reconstruction r-Step

The reconstruction r - s t e p of Section 2.2.4 applies the f - bc sub-step for all non-corner cut links, and then optionally involves sub-steps f - equil , f - extr , and f - bc - corner ; their implementation is discussed below.

Appendix B.3.1. Fill-Boundary f-bc Sub-Step

two - nb and one - nb c u t l i n k s c q define f q ( r n , t + 1 ) with bc - s t e p using Equations (10) and (11), where L I k schemes involve the following directional approximations [78] in the “fresh” nodes:
1.
β f q ( r n , t + 1 ) = β f ^ q ( r n c q , t ) due to streaming from r n c q to r n .
2.
β ^ f ^ q ( r n , t ) β ^ f q ( r n c q , t + 1 ) due to potential streaming from r n to r n c q .
3.
α ^ f ^ q ( r n , t ) α ^ ( 2 f ^ q ( r n c q , t ) f q ( r n c q , t + 1 ) ) where f q ( r n c q , t + 1 ) replaces f ^ q ( r n 2 c q , t ) in the two-point linear interpolation due to streaming step f ^ q ( r n 2 c q , t ) = f q ( r n c q , t + 1 ) .
4.
n ^ q ( r n , t ) n ^ q ( r n c q , t ) in Equation (11) with the second-order accuracy. Since n ^ q ( r n c q , t ) is not available in no - nb   links , n ^ q ( r n , t + 1 ) = 0 applies there as an initial guess and then, optionally, r - iter updates this approximation.

Appendix B.3.2. Preliminary Velocity Approximation a-vel

In the “fresh” nodes, the filling with f - equil needs to approximate velocity u ( r n , t ^ ) by first carrying out the a - vel sub-step. f - bc basically needs u ( r n , t ^ ) only on the first iterative sub-step for parabolic correction 1 2 α ( u ) q 2 e q when using class L I 5 , or else in the term of α ( p ) e q + in Equation (12) when using class E L I + with α ( p ) I n s e 0 . The approximation of velocity u ( r n , t ^ ) can be computed as follows:
1.
in the presence of corner links, u ( r n , t ^ ) averages their linear velocity approximations which are computed from the two Dirichlet values u b ( r n + δ ± q , t ^ ) .
2.
otherwise, u ( r n , t ^ ) averages the first-order accurate linear interpolations computed along the cut links between the Dirichlet value u b ( r n + δ q c q , t ^ ) and the known value u ( r n c q , t ^ ) .

Appendix B.3.3. Fill-Extrapolation Sub-Step f-extr

In a few cases, an averaged extrapolation performs over all “old” neighbors of a “fresh” node r n :
1.
{ f q ( r n , t + 1 ) } for tangential links and no - nb   links with the r - extr sub-step.
2.
{ f ^ q ( r n , t ) } with the r - bc - corner sub-step in no - nb   links as an initial guess.
A formally more accurate approximation, such as an average of directional linear extrapolations, has not been not observed to improve the whole accuracy of the scheme; however, these extrapolations considerably grow the number of communications in parallel implementation.

Appendix B.3.4. Fill-Equilibrium Sub-Step f-equil

When the first two filling sub-steps f - stream and f - bc are completed, the equilibrium filling f - equil applies using Equation (4) with ρ ^ = ρ 0 in the following form in the set T e :
f q ( r n , t + 1 ) = e q ( ρ , u ( r n , t ^ ) ) , q T e , e q ( ρ , u ) = ρ e q ( 1 , 0 ) + e q ( 0 , u ) ,
where the current local density value ρ is determined from its definition: ρ = q T e f q ( r n , t + 1 ) + q T e f q ( r n , t + 1 ) by solving the following linear equation:
ρ ( 1 q T e e q ( 1 , 0 ) ) = q T e e q ( 0 , u ) + q T e f q ( r n , t + 1 ) .
Here, set T e includes an immobile population f 0 ( r n , t + 1 ) with any one of the three reconstruction sub-steps: (i) r - equil , (ii) r - extr , and (iii) r - bc - corner ; set T e also contains tangential links and no - nb links with r - equil , where u = u ( r n , t ^ ) was first set by a preliminary approximation a - vel . On the other hand, r - extr and r - bc - corner , but also r - equil on the following sub-iterations of r - iter , derive u = ρ 0 1 j from momentum j ( r n , t ^ ) of the moving populations, already reconstructed in a given “fresh” node.

Appendix B.4. Deformable Channel Flow: Simulation Results

We consider the 2D straight channel ( x , y ) where the inlet and outlet cross-sections extend continuously along the x-axis according to an applied positive parabolic profile u b ( y ) = u m a x y 2 h 2 1 x ; the maximum stretching speed max y u b ( y ) = u m a x 1 x then occurs on the two horizontal walls y = ± h , while the bisection points with the central axis y = 0 remain stationary. The Stokes and N S E channel solutions should develop linear pressure gradient p = 2 ν u m a x h 2 1 x and an associated parabolic velocity profile u t h = u b ( y ) ; an example of a numerical solution is shown in Figure A1. Note that the magnitude of an acute angle between the outlet transverse shape and the horizontal wall continuously diminishes, making this benchmark stimulating for the L B M - D B application to circular contractions in “parabolic” (triangular) pipes. All figures display the deviation from an exact solution for all fluid points together and in the last three channel configurations in Figure A1.
Figure A1. The set of fluid grid points ( x , y [ h , h ] ) is displayed at t = { 0 , 300 , 500 , 700 } steps (from left to right and top to bottom). The entry and exit are initially vertical and distanced by L = 49 [ l . u ] ; they are subjected to permanent stretching according to the parabolic velocity profile u x t h ( y ) = u m a x y 2 h 2 1 x . Data [ l . u ]: h = 10 , ν = 1 , u m a x = 0.05 .
Figure A1. The set of fluid grid points ( x , y [ h , h ] ) is displayed at t = { 0 , 300 , 500 , 700 } steps (from left to right and top to bottom). The entry and exit are initially vertical and distanced by L = 49 [ l . u ] ; they are subjected to permanent stretching according to the parabolic velocity profile u x t h ( y ) = u m a x y 2 h 2 1 x . Data [ l . u ]: h = 10 , ν = 1 , u m a x = 0.05 .
Fluids 10 00022 g0a1
Figure A2, Figure A3, Figure A4 and Figure A5 display velocity and pressure error estimates in the four diagrams for each combination of the boundary schemes, as err ( u x ) = u x ( x , y ) u x t h ( y ) versus y (a) and x (b), err ( p ) = p ( x , y ) p t h ( x ) versus x (c), and the set of two estimates { L 2 ( u x ) , L 2 ( p ) } (d), which are computed as
L 2 ( ϕ ) = r V f ( t ) ( ϕ ϕ t h ) 2 r V f ( t ) ϕ t h 2 , ϕ ( r , t ) .
In this, p t h ( x ) = p 0 ( t ) + p ( x x 0 ) is traced through a numerical value of p 0 ( t ) at a fixed point ( x 0 , y 0 ) in the center of the initial fluid domain.
Figure A2 compares the B B (left) and B F L 0 (right) applying them for all cut link. The B B displays velocity fluctuations up to six nodes from the horizontal walls (first diagram), in both entry and exit moving fronts, and they are associated with the destabilizing pressure oscillation on the outlet interface. B F L 0 effectively improves the solution, concentrating stronger velocity deviations (i) near the horizontal walls and (ii) on the moving exit front; B F L 0 reduces an amplitude of err ( u x ) by a factor of two to seven and the pressure discrepancy on the moving front by an order of magnitude. The last two diagrams show that B B produces velocity error estimate which oscillates around a time-independent value of L 2 ( u x ) 2 × 10 2 while the pressure estimate matches L 2 ( p ) [ 4 , 7 ] × 10 3 until the B B loses its approximate solution. At the same time, B F L 0 monotonously increases both error estimates over time and locates velocity error in an interval one order smaller, as L 2 ( u x ) [ 2 , 4 ] × 10 3 when L 2 ( p ) [ 1 , 7 ] × 10 3 .
Figure A2. Velocity and pressure deviations from the exact solution are output for all fluid grid nodes in a deformable channel when t = { 300 , 500 , 700 } steps . (top) err ( u x ) versus y; (second from top) err ( u x ) versus x; (second from bottom) err ( p ) versus x; (bottom) { L 2 ( u x ) , L 2 ( p ) } . This figure: all cut links apply B B (left) and B F L 0 (right).
Figure A2. Velocity and pressure deviations from the exact solution are output for all fluid grid nodes in a deformable channel when t = { 300 , 500 , 700 } steps . (top) err ( u x ) versus y; (second from top) err ( u x ) versus x; (second from bottom) err ( p ) versus x; (bottom) { L 2 ( u x ) , L 2 ( p ) } . This figure: all cut links apply B B (left) and B F L 0 (right).
Fluids 10 00022 g0a2
Figure A3 applies the most accurate schemes, A V M R for two - nb links and B F L 5 for one - nb links , while no - nb links apply either B B (left) or B F L 4 (right). Although B B only applies in the corner links, it spoils the whole solution and produces even larger velocity and, above all, pressure errors than in Figure A2 (left). However, even if B F L 4 replaces B B in the corner links alone, the velocity error near the horizontal walls reduces by a factor of three to six, and the magnitude of pressure deviation on the moving front diminishes by more than an order. Overall, combination A V M R B F L 5 B F L 4 (right) produces similar error distributions with B F L 0 in Figure A2 (right), but it further decreases velocity error across the channel and pressure discrepancy on the deforming front. Both error estimates L 2 ( u x ) and L 2 ( p ) now increase monotonously over time in both simulations; B F L 4 in the corners (right) reduces their magnitudes by one order compared to B B (left), and also diminishes both error estimates compared to B F L 0 applied in all links (Figure A2, right).
Figure A3. Similar as in Figure A2. (left): A V M R - B F L 5 - B B . (right): A V M R - B F L 5 - B F L 4 . (left,right): two - nb   links apply A V M R ; one - nb   links apply B F L 5 ; no - nb   links apply B B (left) and B F L 4 (right).
Figure A3. Similar as in Figure A2. (left): A V M R - B F L 5 - B B . (right): A V M R - B F L 5 - B F L 4 . (left,right): two - nb   links apply A V M R ; one - nb   links apply B F L 5 ; no - nb   links apply B B (left) and B F L 4 (right).
Fluids 10 00022 g0a3
Figure A4 keeps B F L 4 in the corner links, otherwise it applies p - pressure B F L 4 (left) and p - p B F L 5 (right). Consistent with its parabolic ( p - flow ) velocity accuracy, B F L 5 slightly diminishes the velocity error but shows almost the same pressure distribution due to an equivalent level of pressure accuracy and the same processing of corner links; the corresponding L 2 error estimates are very similar in both cases, with a slightly larger magnitude of L 2 ( u x ) , but a slightly smaller magnitude of L 2 ( p ) , in the case B F L 4 (left). In turn, A V M R B F L 5 B F L 4 in Figure A3 (right) gains slightly compared to B F L 5 B F L 5 B F L 4 in Figure A4 (right) due to advanced pressure accuracy of A V M R compared to B F L 5 .
Figure A4. Similar to Figure A2. (left): B F L 4 B F L 4 B F L 4 where all cut links apply B F L 4 . (right): B F L 5 B F L 5 B F L 4 , where the no - nb   links apply B F L 4 , all others apply B F L 5 .
Figure A4. Similar to Figure A2. (left): B F L 4 B F L 4 B F L 4 where all cut links apply B F L 4 . (right): B F L 5 B F L 5 B F L 4 , where the no - nb   links apply B F L 4 , all others apply B F L 5 .
Fluids 10 00022 g0a4
Figure A5 compares reconstruction algorithms r - equil (left) and r - bc - corner (right) using combination A V M R B F L 5 B F L 0 in the two cases. Recall that r - equil applies an equilibrium reconstruction not only in tangential links but also in corner links. r - equil displays a very high-velocity oscillation in the most acute angles on the moving exit, and it does not decrease using the iterative update r - iter . Furthermore, when replacing B F L 0 with B F L 4 , r - equil destabilizes due to its inaccurate reconstruction of the non-equilibrium in the “fresh” nodes. On the other hand, r - equil shows a smaller pressure deviation from the linear profile on the deforming exit front. Consistently, r - equil manifests large oscillations in L 2 ( u x ) but a smaller magnitude pressure-error estimate L 2 ( p ) . r - bc - corner can reach the r - equil pressure accuracy level by either replacing B F L 0 with the B F L 4 in corners (see in Figure A3 on the right) or reducing the free-tunable collision freedom Λ in Equation (3) towards zero.
Finally, by replacing the r - bc - corner treatment of the corner links in the new fluid nodes by an averaged extrapolation r - extr , the results of A V M R B F L 4 B F L 0 remain very similar to those in Figure A5, because the B F L 0 rule applied here in the corner links is based on the extrapolated population solution. On the other hand, based on the results above, r - bc - corner improves its accuracy by replacing B F L 0 by B F L 4 in corner links, where B F L 4 benefits from the iteratively updated post-collision correction. All in all, these benchmark simulations demonstrate that the treatment of the velocity boundary condition on the corner links can determine total accuracy.
Figure A5. Similar to Figure A2. The reconstruction algorithm is r - equil (left) and r - bc - corner (right). A V M R B F L 4 B F L 0 is applied in two - nb , one - nb and no - nb   links , respectively.
Figure A5. Similar to Figure A2. The reconstruction algorithm is r - equil (left) and r - bc - corner (right). A V M R B F L 4 B F L 0 is applied in two - nb , one - nb and no - nb   links , respectively.
Fluids 10 00022 g0a5

Appendix C. Complementary Numerical Results

Appendix C.1 describes results of the numerical experiments carried out by combining motility pattern M O T - I I from Table 1 with the five occlusion scenarios of Table 4, and compares them with the M O T - I results from Section 3.5; Section 3.6 summarizes this comparison. Appendix C.2 quantifies M O T - I and M O T - I I results in Table A2, Table A3, Table A4, Table A5 and Table A6, and Table A7, Table A8, Table A9, Table A10 and Table A11, respectively; these results are summarized in Table 5.
Table A2. M O T - I O C L - 1 . Data: [ v ¯ m i n , v ¯ m a x ] v w and [ v m i n , v m a x ] v w are multiplied by 10.
Table A2. M O T - I O C L - 1 . Data: [ v ¯ m i n , v ¯ m a x ] v w and [ v m i n , v m a x ] v w are multiplied by 10.
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w × 10 R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w × 10 R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 4.31 , 3.11 ] 1.39 [ 10.0 , 6.37 ] 1.57 0.41 [ 1.50 × 10 3 , 2.49 × 10 3 ] 1.80 × 10 3
F A S T [ 4.31 , 3.59 ] 1.20 [ 10.0 , 6.70 ] 1.50 0.83 [ 7.79 × 10 3 , 5.82 × 10 3 ] 1.19 × 10 3
V M A X [ 4.19 , 4.23 ] 0.99 [ 9.18 , 5.97 ] 1.54 2.01 [ 5.62 × 10 2 , 4.09 × 10 2 ] 7.84 × 10 4
L O V I S S L O W [ 4.25 , 4.02 ] 1.06 [ 9.64 , 6.36 ] 1.52 1.39 [ 2.12 × 10 3 , 1.39 × 10 3 ] 9.13 × 10 4
F A S T [ 4.12 , 4.37 ] 0.94 [ 9.08 , 5.77 ] 1.57 2.69 [ 9.40 × 10 3 , 7.27 × 10 3 ] 7.09 × 10 4
V M A X [ 3.99 , 4.87 ] 0.82 [ 8.41 , 6.10 ] 1.38 6.52 [ 6.64 × 10 2 , 5.35 × 10 2 ] 5.50 × 10 4
C - C L O H I V I S S L O W [ 4.48 , 3.24 ] 1.39 [ 9.38 , 6.41 ] 1.46 0.43 [ 1.71 × 10 3 , 1.71 × 10 3 ] 7.54 × 10 4
F A S T [ 4.46 , 3.83 ] 1.17 [ 9.34 , 6.72 ] 1.39 0.86 [ 8.52 × 10 3 , 4.71 × 10 3 ] 2.28 × 10 4
V M A X [ 4.44 , 4.35 ] 1.02 [ 8.86 , 5.99 ] 1.48 2.13 [ 6.20 × 10 2 , 4.43 × 10 2 ] 1.88 × 10 4
L O V I S S L O W [ 4.42 , 4.18 ] 1.06 [ 8.94 , 6.37 ] 1.40 1.44 [ 2.33 × 10 3 , 1.52 × 10 3 ] 3.50 × 10 5
F A S T [ 4.42 , 4.48 ] 0.99 [ 8.74 , 5.78 ] 1.51 2.89 [ 1.03 × 10 2 , 7.74 × 10 3 ] 2.94 × 10 4
V M A X [ 4.37 , 4.97 ] 0.88 [ 8.09 , 6.32 ] 1.28 7.13 [ 6.96 × 10 2 , 5.35 × 10 2 ] 5.12 × 10 4
P - P E R H I V I S S L O W [ 2.93 , 1.56 ] 1.88 [ 7.07 , 4.29 ] 1.65 0.28 [ 8.83 × 10 4 , 2.35 × 10 3 ] 1.38 × 10 3
F A S T [ 2.99 , 1.63 ] 1.83 [ 6.76 , 4.29 ] 1.58 0.58 [ 4.87 × 10 3 , 5.56 × 10 3 ] 9.26 × 10 4
V M A X [ 3.03 , 1.88 ] 1.61 [ 6.57 , 4.29 ] 1.53 1.45 [ 4.81 × 10 2 , 3.53 × 10 2 ] 5.99 × 10 4
L O V I S S L O W [ 3.02 , 1.71 ] 1.77 [ 6.65 , 4.29 ] 1.55 0.99 [ 1.63 × 10 3 , 1.16 × 10 3 ] 7.08 × 10 4
F A S T [ 3.03 , 1.99 ] 1.52 [ 6.52 , 4.27 ] 1.53 1.98 [ 8.63 × 10 3 , 6.35 × 10 3 ] 5.42 × 10 4
V M A X [ 2.99 , 2.55 ] 1.17 [ 6.36 , 4.02 ] 1.58 4.88 [ 7.00 × 10 2 , 5.33 × 10 2 ] 5.62 × 10 4
C - P E R H I V I S S L O W [ 3.12 , 1.47 ] 2.12 [ 6.87 , 4.08 ] 1.69 0.30 [ 7.83 × 10 4 , 1.61 × 10 3 ] 5.00 × 10 4
F A S T [ 3.18 , 1.55 ] 2.05 [ 6.62 , 4.08 ] 1.62 0.61 [ 5.12 × 10 3 , 3.71 × 10 3 ] 9.72 × 10 5
V M A X [ 3.22 , 2.00 ] 1.61 [ 6.48 , 4.07 ] 1.59 1.54 [ 4.94 × 10 2 , 4.48 × 10 2 ] 2.47 × 10 4
L O V I S S L O W [ 3.21 , 1.83 ] 1.75 [ 6.54 , 4.08 ] 1.60 10.5 × 10 1 [ 1.67 × 10 3 , 1.38 × 10 3 ] 1.18 × 10 4
F A S T [ 3.22 , 2.11 ] 1.53 [ 6.45 , 4.04 ] 1.59 2.10 [ 8.78 × 10 3 , 8.37 × 10 3 ] 3.35 × 10 4
V M A X [ 3.26 , 2.67 ] 1.22 [ 6.39 , 3.86 ] 1.65 5.32 [ 6.82 × 10 2 , 6.94 × 10 2 ] 3.83 × 10 4
Table A3. M O T - I O C L - 2 . Data: [ v ¯ m i n , v ¯ m a x ] v w is multiplied by 10.
Table A3. M O T - I O C L - 2 . Data: [ v ¯ m i n , v ¯ m a x ] v w is multiplied by 10.
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w × 10 R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 5.88 , 3.98 ] 1.48 [ 1.38 , 1.52 ] 9.06 × 10 1 0.56 [ 2.71 × 10 3 , 2.04 × 10 2 ] 1.20 × 10 3
F A S T [ 6.37 , 5.09 ] 1.25 [ 1.49 , 1.55 ] 9.63 × 10 1 1.22 [ 1.36 × 10 2 , 4.46 × 10 2 ] 3.04 × 10 3
V M A X [ 6.99 , 6.20 ] 1.13 [ 1.37 , 1.51 ] 9.06 × 10 1 3.35 [ 1.07 × 10 1 , 1.15 × 10 1 ] 4.10 × 10 3
L O V I S S L O W [ 6.52 , 5.67 ] 1.15 [ 1.46 , 1.54 ] 9.49 × 10 1 2.13 [ 3.98 × 10 3 , 6.75 × 10 3 ] 3.77 × 10 3
F A S T [ 7.33 , 6.48 ] 1.13 [ 1.27 , 1.47 ] 8.63 × 10 1 4.78 [ 1.79 × 10 2 , 1.32 × 10 2 ] 4.32 × 10 3
V M A X [ 7.83 , 7.28 ] 1.07 [ 1.24 , 1.31 ] 9.46 × 10 1 12.8 [ 1.38 × 10 1 , 1.15 × 10 1 ] 5.53 × 10 3
C - C L O H I V I S S L O W [ 6.14 , 4.50 ] 1.36 [ 1.30 , 1.43 ] 9.11 × 10 1 0.59 [ 2.92 × 10 3 , 3.20 × 10 3 ] 1.31 × 10 3
F A S T [ 6.36 , 5.65 ] 1.13 [ 1.35 , 1.43 ] 9.42 × 10 1 1.22 [ 1.57 × 10 2 , 8.03 × 10 3 ] 4.91 × 10 4
V M A X [ 7.45 , 6.70 ] 1.11 [ 1.23 , 1.38 ] 8.91 × 10 1 3.58 [ 1.17 × 10 1 , 8.44 × 10 2 ] 1.46 × 10 3
L O V I S S L O W [ 6.91 , 6.42 ] 1.08 [ 1.29 , 1.41 ] 9.13 × 10 1 2.26 [ 4.38 × 10 3 , 2.77 × 10 3 ] 1.21 × 10 3
F A S T [ 7.82 , 6.87 ] 1.14 [ 1.24 , 1.34 ] 9.28 × 10 1 5.11 [ 1.99 × 10 2 , 1.55 × 10 2 ] 1.53 × 10 3
V M A X [ 8.34 , 7.97 ] 1.05 [ 1.26 , 1.20 ] 1.05 13.6 [ 1.62 × 10 1 , 1.31 × 10 1 ] 7.04 × 10 4
P - P E R H I V I S S L O W [ 4.25 , 3.92 ] 1.09 [ 1.03 , 1.65 ] 6.25 × 10 1 0.41 [ 3.01 × 10 3 , 2.03 × 10 2 ] 2.88 × 10 3
F A S T [ 4.50 , 3.94 ] 1.14 [ 1.01 , 1.68 ] 6.04 × 10 1 0.87 [ 8.27 × 10 3 , 4.44 × 10 2 ] 4.23 × 10 3
V M A X [ 4.72 , 4.08 ] 1.16 [ 0.96 , 1.66 ] 5.82 × 10 1 2.26 [ 8.76 × 10 2 , 1.15 × 10 1 ] 4.98 × 10 3
L O V I S S L O W [ 4.65 , 4.02 ] 1.16 [ 1.00 , 1.68 ] 5.97 × 10 1 1.52 [ 2.81 × 10 3 , 6.73 × 10 3 ] 4.75 × 10 3
F A S T [ 4.75 , 4.11 ] 1.16 [ 0.93 , 1.62 ] 5.74 × 10 1 3.10 [ 1.63 × 10 2 , 1.35 × 10 2 ] 5.11 × 10 3
V M A X [ 4.69 , 4.97 ] 0.94 [ 0.98 , 1.48 ] 6.59 × 10 1 7.66 [ 1.50 × 10 1 , 1.18 × 10 1 ] 5.20 × 10 3
C - P E R H I V I S S L O W [ 4.80 , 3.76 ] 1.28 [ 0.10 , 1.53 ] 6.51 × 10 1 0.46 [ 1.95 × 10 3 , 3.16 × 10 3 ] 2.43 × 10 4
F A S T [ 5.06 , 3.91 ] 1.30 [ 1.00 , 1.54 ] 6.52 × 10 1 0.97 [ 9.64 × 10 3 , 7.10 × 10 3 ] 1.14 × 10 3
V M A X [ 5.22 , 3.99 ] 1.31 [ 0.90 , 1.51 ] 5.98 × 10 1 2.51 [ 1.06 × 10 1 , 6.93 × 10 2 ] 1.93 × 10 3
L O V I S S L O W [ 5.17 , 3.94 ] 1.31 [ 0.93 , 1.53 ] 6.07 × 10 1 1.69 [ 3.46 × 10 3 , 2.14 × 10 3 ] 1.71 × 10 3
F A S T [ 5.24 , 4.38 ] 1.20 [ 0.88 , 1.48 ] 5.99 × 10 1 3.42 [ 1.92 × 10 2 , 1.32 × 10 2 ] 1.99 × 10 3
V M A X [ 5.06 , 5.15 ] 0.98 [ 0.89 , 1.39 ] 6.40 × 10 1 8.26 [ 1.68 × 10 1 , 1.27 × 10 1 ] 1.09 × 10 3
Table A4. M O T - I O C L - 3 .
Table A4. M O T - I O C L - 3 .
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 1.78 , 0.67 ] 2.66 [ 4.23 , 2.06 ] 2.05 1.71 [ 4.75 × 10 2 , 1.24 × 10 1 ] 5.45 × 10 2
F A S T [ 1.88 , 0.80 ] 2.34 [ 4.62 , 2.25 ] 2.05 3.61 [ 1.14 × 10 1 , 2.96 × 10 1 ] 4.85 × 10 2
V M A X [ 1.89 , 1.16 ] 1.63 [ 4.69 , 2.35 ] 2.00 9.07 [ 4.16 × 10 1 , 1.12 ] 4.61 × 10 2
L O V I S S L O W [ 1.92 , 1.00 ] 1.93 [ 4.75 , 2.33 ] 2.04 6.26 [ 2.09 × 10 2 , 5.23 × 10 2 ] 4.67 × 10 2
F A S T [ 1.94 , 1.27 ] 1.53 [ 4.63 , 2.34 ] 1.98 12.7 [ 5.56 × 10 2 , 1.59 × 10 1 ] 4.55 × 10 2
V M A X [ 1.96 , 1.46 ] 1.34 [ 4.31 , 2.73 ] 1.58 32.0 [ 2.95 × 10 1 , 7.61 × 10 1 ] 3.97 × 10 2
C - C L O H I V I S S L O W [ 1.90 , 0.61 ] 3.11 [ 3.79 , 1.21 ] 3.13 1.83 [ 5.90 × 10 3 , 8.13 × 10 3 ] 1.60 × 10 2
F A S T [ 1.97 , 0.74 ] 2.68 [ 4.13 , 1.29 ] 3.19 3.79 [ 1.70 × 10 2 , 1.68 × 10 2 ] 1.04 × 10 2
V M A X [ 2.03 , 0.95 ] 2.15 [ 4.12 , 1.23 ] 3.36 9.75 [ 1.34 × 10 1 , 9.67 × 10 2 ] 6.22 × 10 3
L O V I S S L O W [ 2.03 , 0.87 ] 2.34 [ 4.18 , 1.26 ] 3.32 6.62 [ 4.89 × 10 3 , 3.27 × 10 3 ] 7.67 × 10 3
F A S T [ 2.01 , 1.01 ] 2.00 [ 4.00 , 1.24 ] 3.23 13.1 [ 2.27 × 10 2 , 1.70 × 10 2 ] 5.38 × 10 3
V M A X [ 1.91 , 1.19 ] 1.61 [ 3.52 , 1.41 ] 2.51 31.2 [ 1.72 × 10 1 , 1.40 × 10 1 ] 6.21 × 10 3
P - P E R H I V I S S L O W [ 0.83 , 0.52 ] 1.59 [ 2.50 , 2.20 ] 1.14 0.80 [ 4.73 × 10 2 , 1.22 × 10 1 ] 3.89 × 10 2
F A S T [ 1.04 , 0.59 ] 1.76 [ 3.01 , 2.38 ] 1.26 1.99 [ 1.08 × 10 1 , 2.94 × 10 1 ] 3.71 × 10 2
V M A X [ 1.28 , 0.72 ] 1.78 [ 3.28 , 2.49 ] 1.32 6.12 [ 3.72 × 10 1 , 1.12 ] 3.85 × 10 2
L O V I S S L O W [ 1.19 , 0.67 ] 1.79 [ 3.26 , 2.46 ] 1.33 3.89 [ 1.93 × 10 2 , 5.20 × 10 2 ] 3.75 × 10 2
F A S T [ 1.32 , 0.75 ] 1.75 [ 3.15 , 2.49 ] 1.26 8.61 [ 4.84 × 10 2 , 1.58 × 10 1 ] 3.97 × 10 2
V M A X [ 1.35 , 0.83 ] 1.63 [ 3.13 , 2.96 ] 1.06 22.0 [ 2.21 × 10 1 , 7.31 × 10 1 ] 4.36 × 10 2
C - P E R H I V I S S L O W [ 0.92 , 0.44 ] 2.09 [ 2.42 , 2.00 ] 1.21 0.88 [ 4.10 × 10 3 , 6.00 × 10 3 ] 7.18 × 10 3
F A S T [ 1.02 , 0.48 ] 2.13 [ 2.39 , 2.02 ] 1.18 1.95 [ 1.03 × 10 2 , 1.28 × 10 2 ] 4.01 × 10 3
V M A X [ 1.09 , 0.57 ] 1.92 [ 2.27 , 2.05 ] 1.11 5.23 [ 1.11 × 10 1 , 5.43 × 10 2 ] 1.72 × 10 3
L O V I S S L O W [ 1.07 , 0.52 ] 2.06 [ 2.33 , 2.05 ] 1.14 3.50 [ 3.47 × 10 3 , 2.33 × 10 3 ] 2.49 × 10 3
F A S T [ 1.09 , 0.60 ] 1.81 [ 2.20 , 2.02 ] 1.09 7.11 [ 2.08 × 10 2 , 9.13 × 10 3 ] 1.32 × 10 3
V M A X [ 1.04 , 0.66 ] 1.57 [ 1.90 , 1.77 ] 1.08 16.9 [ 1.74 × 10 1 , 8.83 × 10 2 ] 2.58 × 10 3
Table A5. M O T - I O C L - 4 .
Table A5. M O T - I O C L - 4 .
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 1.32 , 7.18 × 10 1 ] 1.83 [ 3.17 , 2.90 ] 1.09 1.26 [ 6.82 × 10 2 , 1.19 × 10 1 ] 3.73 × 10 3
F A S T [ 1.45 , 1.09 ] 1.33 [ 3.65 , 3.02 ] 1.21 2.78 [ 1.99 × 10 1 , 3.22 × 10 1 ] 2.97 × 10 3
V M A X [ 1.65 , 1.35 ] 1.22 [ 4.11 , 3.27 ] 1.26 7.91 [ 6.99 × 10 1 , 1.49 ] 2.81 × 10 2
L O V I S S L O W [ 1.60 , 1.28 ] 1.25 [ 4.00 , 3.17 ] 1.26 5.22 [ 3.64 × 10 2 , 6.64 × 10 2 ] 1.74 × 10 2
F A S T [ 1.69 , 1.37 ] 1.23 [ 4.09 , 3.25 ] 1.26 11.0 [ 8.83 × 10 2 , 2.22 × 10 1 ] 3.44 × 10 2
V M A X [ 2.05 , 1.69 ] 1.22 [ 4.58 , 3.25 ] 1.41 33.5 [ 5.90 × 10 1 , 1.27 ] 4.91 × 10 2
C - C L O H I V I S S L O W [ 1.47 , 6.32 × 10 1 ] 2.33 [ 3.01 , 3.55 ] 0.85 1.41 [ 6.65 × 10 3 , 9.03 × 10 3 ] 2.21 × 10 2
F A S T [ 1.62 , 8.32 × 10 1 ] 1.95 [ 3.47 , 3.60 ] 0.96 3.11 [ 2.52 × 10 2 , 1.90 × 10 2 ] 1.44 × 10 2
V M A X [ 1.81 , 1.18 ] 1.53 [ 3.66 , 3.45 ] 1.06 8.67 [ 2.03 × 10 1 , 1.37 × 10 1 ] 9.09 × 10 3
L O V I S S L O W [ 1.75 , 1.06 ] 1.64 [ 3.64 , 3.54 ] 1.03 5.70 [ 7.55 × 10 3 , 4.49 × 10 3 ] 1.08 × 10 2
F A S T [ 1.83 , 1.24 ] 1.47 [ 3.65 , 3.36 ] 1.09 11.9 [ 3.35 × 10 2 , 2.44 × 10 2 ] 8.45 × 10 3
V M A X [ 1.91 , 1.50 ] 1.28 [ 3.53 , 3.00 ] 1.18 31.1 [ 2.60 × 10 1 , 2.12 × 10 1 ] 1.49 × 10 2
P - P E R H I V I S S L O W [ 0.79 , 6.42 × 10 1 ] 1.24 [ 2.41 , 3.49 ] 0.69 0.76 [ 6.78 × 10 2 , 1.18 × 10 1 ] 2.41 × 10 2
F A S T [ 1.07 , 6.86 × 10 1 ] 1.55 [ 2.97 , 3.52 ] 0.85 2.05 [ 2.01 × 10 1 , 3.19 × 10 1 ] 5.14 × 10 3
V M A X [ 1.42 , 7.88 × 10 1 ] 1.80 [ 3.28 , 3.61 ] 0.91 6.82 [ 7.21 × 10 1 , 1.49 ] 3.17 × 10 2
L O V I S S L O W [ 1.29 , 7.40 × 10 1 ] 1.74 [ 3.23 , 3.57 ] 0.91 4.20 [ 3.70 × 10 2 , 6.63 × 10 2 ] 2.28 × 10 2
F A S T [ 1.51 , 8.48 × 10 1 ] 1.78 [ 3.66 , 3.60 ] 1.02 9.85 [ 9.28 × 10 2 , 2.22 × 10 1 ] 3.63 × 10 2
V M A X [ 1.76 , 1.06 ] 1.66 [ 4.41 , 3.46 ] 1.28 28.8 [ 5.39 × 10 1 , 1.28 ] 3.79 × 10 2
C - P E R H I V I S S L O W [ 0.94 , 6.54 × 10 1 ] 1.44 [ 2.56 , 3.88 ] 0.66 0.90 [ 6.87 × 10 3 , 9.01 × 10 3 ] 1.15 × 10 2
F A S T [ 1.10 , 6.45 × 10 1 ] 1.70 [ 2.44 , 3.92 ] 0.62 2.11 [ 1.71 × 10 2 , 1.92 × 10 2 ] 6.68 × 10 3
V M A X [ 1.21 , 0.72 ] 1.69 [ 2.22 , 3.88 ] 0.57 5.82 [ 1.59 × 10 1 , 7.40 × 10 2 ] 3.91 × 10 3
L O V I S S L O W [ 1.18 , 6.94 × 10 1 ] 1.70 [ 2.29 , 3.92 ] 0.59 3.86 [ 4.90 × 10 3 , 3.23 × 10 3 ] 4.71 × 10 3
F A S T [ 1.22 , 7.75 × 10 1 ] 1.57 [ 2.33 , 3.77 ] 0.62 7.95 [ 3.01 × 10 2 , 1.50 × 10 2 ] 3.71 × 10 3
V M A X [ 1.32 , 8.54 × 10 1 ] 1.55 [ 2.47 , 3.07 ] 0.80 21.6 [ 2.59 × 10 1 , 1.41 × 10 1 ] 8.02 × 10 3
Table A6. M O T - I O C L - 5 .
Table A6. M O T - I O C L - 5 .
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
C - C L O H I V I S S L O W [ 3.65 , 7.36 × 10 1 ] 4.96 [ 6.81 , 3.40 ] 2.00 3.51 [ 1.52 × 10 2 , 2.32 × 10 2 ] 3.89 × 10 2
F A S T [ 4.21 , 8.57 × 10 1 ] 4.91 [ 7.55 , 3.21 ] 2.35 8.08 [ 4.28 × 10 2 , 4.85 × 10 2 ] 2.19 × 10 2
V M A X [ 4.43 , 9.69 × 10 1 ] 4.57 [ 8.19 , 2.77 ] 2.95 21.3 [ 2.35 × 10 1 , 1.48 × 10 1 ] 9.69 × 10 3
L O V I S S L O W [ 4.42 , 8.55 × 10 1 ] 5.18 [ 7.94 , 2.99 ] 2.66 14.4 [ 9.96 × 10 3 , 7.67 × 10 3 ] 1.36 × 10 2
F A S T [ 4.38 , 1.08 ] 4.04 [ 8.35 , 2.55 ] 3.28 28.6 [ 3.56 × 10 2 , 2.29 × 10 2 ] 8.01 × 10 3
V M A X [ 4.26 , 1.51 ] 2.81 [ 8.83 , 2.07 ] 4.26 69.5 [ 2.52 × 10 1 , 2.02 × 10 1 ] 1.80 × 10 2
C - P E R H I V I S S L O W [ 1.05 , 6.75 × 10 1 ] 1.56 [ 5.13 , 4.96 ] 1.03 1.01 [ 8.93 × 10 3 , 2.18 × 10 2 ] 1.21 × 10 2
F A S T [ 1.02 , 7.04 × 10 1 ] 1.45 [ 5.05 , 4.99 ] 1.01 1.96 [ 2.11 × 10 2 , 4.60 × 10 2 ] 5.63 × 10 3
V M A X [ 1.36 , 8.46 × 10 1 ] 1.61 [ 4.52 , 4.91 ] 0.92 6.54 [ 1.24 × 10 1 , 1.28 × 10 1 ] 3.02 × 10 3
L O V I S S L O W [ 1.19 , 7.72 × 10 1 ] 1.54 [ 4.81 , 4.97 ] 0.97 3.87 [ 4.01 × 10 3 , 7.13 × 10 3 ] 3.43 × 10 3
F A S T [ 1.50 , 9.04 × 10 1 ] 1.66 [ 4.29 , 4.85 ] 0.88 9.78 [ 2.43 × 10 2 , 1.63 × 10 2 ] 3.50 × 10 3
V M A X [ 1.76 , 9.95 × 10 1 ] 1.77 [ 3.62 , 4.26 ] 0.85 28.8 [ 2.29 × 10 1 , 1.07 × 10 1 ] 1.12 × 10 2
Table A7. M O T - I I O C L - 1 . Data: [ v ¯ m i n , v ¯ m a x ] v w and [ v m i n , v m a x ] v w  are multiplied by 10.
Table A7. M O T - I I O C L - 1 . Data: [ v ¯ m i n , v ¯ m a x ] v w and [ v m i n , v m a x ] v w  are multiplied by 10.
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w × 10 R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w × 10 R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 7.44 , 4.27 ] 1.74 [ 15.6 , 9.00 ] 1.73 0.71 [ 3.39 × 10 3 , 4.94 × 10 3 ] 9.93 × 10 4
F A S T [ 7.63 , 4.86 ] 1.57 [ 16.4 , 9.62 ] 1.70 1.46 [ 1.54 × 10 2 , 1.77 × 10 2 ] 2.18 × 10 3
V M A X [ 8.33 , 5.42 ] 1.54 [ 15.9 , 8.90 ] 1.79 4.00 [ 1.02 × 10 1 , 1.12 × 10 1 ] 3.31 × 10 3
L O V I S S L O W [ 7.92 , 5.09 ] 1.56 [ 16.2 , 9.32 ] 1.74 2.59 [ 4.01 × 10 3 , 4.41 × 10 3 ] 2.83 × 10 3
F A S T [ 8.70 , 5.68 ] 1.53 [ 15.6 , 8.60 ] 1.81 5.68 [ 1.64 × 10 2 , 1.79 × 10 2 ] 3.76 × 10 3
V M A X [ 10.0 , 6.26 ] 1.60 [ 15.1 , 8.01 ] 1.88 16.3 [ 1.10 × 10 1 , 1.11 × 10 1 ] 5.96 × 10 3
C - C L O H I V I S S L O W [ 7.13 , 4.49 ] 1.59 [ 14.6 , 9.15 ] 1.60 0.68 [ 3.49 × 10 3 , 3.98 × 10 3 ] 8.71 × 10 4
F A S T [ 7.27 , 5.04 ] 1.44 [ 15.1 , 9.71 ] 1.56 1.40 [ 1.53 × 10 2 , 1.43 × 10 2 ] 6.32 × 10 5
V M A X [ 8.27 , 5.62 ] 1.47 [ 14.8 , 8.98 ] 1.65 3.97 [ 9.98 × 10 2 , 8.50 × 10 2 ] 4.14 × 10 4
L O V I S S L O W [ 7.81 , 5.28 ] 1.48 [ 15.0 , 9.40 ] 1.60 2.55 [ 3.93 × 10 3 , 3.43 × 10 3 ] 2.56 × 10 4
F A S T [ 8.67 , 5.87 ] 1.48 [ 14.5 , 8.71 ] 1.67 5.66 [ 1.61 × 10 2 , 1.35 × 10 2 ] 5.08 × 10 4
V M A X [ 9.91 , 6.43 ] 1.54 [ 14.2 , 8.16 ] 1.75 16.2 [ 1.12 × 10 1 , 8.83 × 10 2 ] 7.48 × 10 4
P - P E R H I V I S S L O W [ 4.12 , 2.97 ] 1.39 [ 8.39 , 7.44 ] 1.13 0.40 [ 3.53 × 10 3 , 4.75 × 10 3 ] 7.45 × 10 4
F A S T [ 4.43 , 3.17 ] 1.40 [ 8.62 , 7.50 ] 1.15 0.85 [ 1.56 × 10 2 , 1.70 × 10 2 ] 1.45 × 10 3
V M A X [ 4.67 , 3.31 ] 1.41 [ 7.60 , 7.13 ] 1.06 2.24 [ 9.85 × 10 2 , 1.08 × 10 1 ] 2.26 × 10 3
L O V I S S L O W [ 4.60 , 3.26 ] 1.41 [ 8.15 , 7.39 1.10 1.50 [ 3.94 × 10 3 , 4.24 × 10 3 ] 1.90 × 10 3
F A S T [ 4.67 , 3.39 ] 1.38 [ 7.62 , 6.82 ] 1.12 3.05 [ 1.57 × 10 2 , 1.74 × 10 2 ] 2.60 × 10 3
V M A X [ 4.58 , 3.58 ] 1.28 [ 7.00 , 6.30 ] 1.11 7.47 [ 9.98 × 10 2 , 1.14 × 10 1 ] 3.86 × 10 3
C - P E R H I V I S S L O W [ 4.18 , 3.02 ] 1.38 [ 8.26 , 7.04 ] 1.17 0.40 [ 3.56 × 10 3 , 3.83 × 10 3 ] 8.05 × 10 4
F A S T [ 4.45 , 3.19 ] 1.39 [ 8.41 , 7.04 ] 1.19 0.85 [ 1.52 × 10 2 , 1.37 × 10 2 ] 2.95 × 10 4
V M A X [ 4.58 , 3.41 ] 1.34 [ 7.46 , 6.50 ] 1.15 2.20 [ 9.57 × 10 2 , 8.50 × 10 2 ] 1.00 × 10 4
L O V I S S L O W [ 4.56 , 3.27 ] 1.39 [ 7.89 , 6.82 ] 1.16 1.49 [ 3.83 × 10 3 , 3.31 × 10 3 ] 4.77 × 10 5
F A S T [ 4.56 , 3.54 ] 1.29 [ 7.40 , 6.17 ] 1.20 2.97 [ 1.52 × 10 2 , 1.39 × 10 2 ] 2.07 × 10 4
V M A X [ 4.49 , 3.79 ] 1.19 [ 7.54 , 5.85 ] 1.29 7.33 [ 9.80 × 10 2 , 8.99 × 10 2 ] 2.92 × 10 4
Table A8. M O T - I I O C L - 2 .
Table A8. M O T - I I O C L - 2 .
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 1.31 , 0.92 ] 1.42 [ 3.08 , 2.87 ] 1.07 1.26 [ 8.69 × 10 3 , 1.15 × 10 2 ] 4.46 × 10 4
F A S T [ 1.54 , 1.04 ] 1.48 [ 3.78 , 2.95 ] 1.28 2.96 [ 3.81 × 10 2 , 4.71 × 10 2 ] 2.98 × 10 3
V M A X [ 1.72 , 1.46 ] 1.18 [ 3.94 , 2.57 ] 1.53 8.28 [ 2.49 × 10 1 , 0.30 ] 5.86 × 10 3
L O V I S S L O W [ 1.64 , 1.29 ] 1.27 [ 4.06 , 2.77 ] 1.47 5.34 [ 9.93 × 10 3 , 1.2 × 10 2 ] 8.64 × 10 3
F A S T [ 1.78 , 1.54 ] 1.16 [ 3.79 , 2.35 ] 1.61 11.6 [ 4.00 × 10 2 , 5.00 × 10 2 ] 7.03 × 10 3
V M A X [ 2.01 , 1.67 ] 1.21 [ 3.31 , 2.49 ] 1.33 32.9 [ 2.68 × 10 1 , 0.35 ] 1.47 × 10 2
C - C L O H I V I S S L O W [ 1.38 , 1.00 ] 1.38 [ 3.02 , 2.71 ] 1.11 1.32 [ 9.42 × 10 3 , 9.59 × 10 3 ] 1.21 × 10 2
F A S T [ 1.56 , 1.17 ] 1.33 [ 3.49 , 2.61 ] 1.34 3.00 [ 4.20 × 10 2 , 3.54 × 10 2 ] 8.57 × 10 3
V M A X [ 1.74 , 1.45 ] 1.20 [ 3.48 , 2.36 ] 1.48 8.35 [ 2.77 × 10 1 , 2.34 × 10 1 ] 7.63 × 10 3
L O V I S S L O W [ 1.66 , 1.34 ] 1.24 [ 3.61 , 2.38 ] 1.52 5.41 [ 1.10 × 10 2 , 9.01 × 10 3 ] 1.02 × 10 2
F A S T [ 1.84 , 1.49 ] 1.23 [ 3.38 , 2.26 ] 1.50 12.0 [ 4.41 × 10 2 , 3.85 × 10 2 ] 7.92 × 10 3
V M A X [ 2.22 , 1.43 ] 1.55 [ 3.30 , 2.14 ] 1.54 36.3 [ 2.91 × 10 1 , 2.70 × 10 1 ] 9.75 × 10 3
P - P E R H I V I S S L O W [ 0.79 , 0.78 ] 1.01 [ 1.79 , 2.59 ] 0.69 0.76 [ 7.99 × 10 3 , 1.5 × 10 2 ] 3.17 × 10 3
F A S T [ 0.94 , 0.83 ] 1.13 [ 2.06 , 2.64 ] 0.78 1.80 [ 3.82 × 10 2 , 4.83 × 10 2 ] 5.22 × 10 3
V M A X [ 1.10 , 0.87 ] 1.27 [ 2.03 , 2.50 ] 0.81 5.30 [ 2.54 × 10 1 , 2.92 × 10 1 ] 7.14 × 10 3
L O V I S S L O W [ 1.04 , 0.84 ] 1.23 [ 2.09 , 2.62 ] 0.80 3.38 [ 1.01 × 10 2 , 1.18 × 10 2 ] 6.29 × 10 3
F A S T [ 1.15 , 0.92 ] 1.26 [ 1.96 , 2.34 ] 0.84 7.53 [ 4.02 × 10 2 , 4.70 × 10 2 ] 8.04 × 10 3
V M A X [ 1.21 , 1.01 ] 1.19 [ 1.95 , 1.78 ] 1.10 19.7 [ 2.15 × 10 1 , 3.29 × 10 1 ] 1.32 × 10 2
C - P E R H I V I S S L O W [ 0.93 , 0.87 ] 1.07 [ 1.81 , 2.46 ] 0.74 0.89 [ 9.15 × 10 3 , 9.93 × 10 3 ] 7.83 × 10 3
F A S T [ 1.07 , 0.90 ] 1.20 [ 2.01 , 2.47 ] 0.82 2.06 [ 4.22 × 10 2 , 3.59 × 10 2 ] 5.34 × 10 3
V M A X [ 1.25 , 0.88 ] 1.41 [ 1.92 , 2.19 ] 0.87 5.98 [ 2.72 × 10 1 , 2.21 × 10 1 ] 4.37 × 10 3
L O V I S S L O W [ 1.18 , 0.86 ] 1.37 [ 2.00 , 2.36 ] 0.85 3.84 [ 1.09 × 10 2 , 8.80 × 10 3 ] 4.53 × 10 3
F A S T [ 1.28 , 0.93 ] 1.39 [ 1.85 , 2.01 ] 0.92 8.38 [ 4.32 × 10 2 , 3.58 × 10 2 ] 4.47 × 10 3
V M A X [ 1.23 , 0.99 ] 1.24 [ 1.74 , 1.72 ] 1.01 20.1 [ 2.78 × 10 1 , 2.71 × 10 1 ] 7.04 × 10 3
Table A9. M O T - I I O C L - 3 ; P - C L O and P - P E R replace c s 2 = 1 3 by c s 2 = 1 3 × 4 2 in the case LOVIS − VMAX★★.
Table A9. M O T - I I O C L - 3 ; P - C L O and P - P E R replace c s 2 = 1 3 by c s 2 = 1 3 × 4 2 in the case LOVIS − VMAX★★.
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 3.20 , 7.30 × 10 1 ] 4.38 [ 7.50 , 3.41 ] 2.20 3.07 [ 1.92 × 10 1 , 3.29 × 10 1 ] 1.69 × 10 1
F A S T [ 3.33 , 9.18 × 10 1 ] 3.63 [ 7.92 , 4.09 ] 1.94 6.40 [ 5.46 × 10 1 , 9.42 × 10 1 ] 1.47 × 10 1
V M A X [ 3.60 , 1.64 ] 2.19 [ 8.24 , 5.92 ] 1.39 17.3 [ 2.04 , 3.21 ] 1.20 × 10 1
L O V I S S L O W [ 3.47 , 1.39 ] 2.50 [ 7.80 , 4.79 ] 1.63 11.3 [ 1.05 × 10 1 , 1.70 × 10 1 ] 1.33 × 10 1
F A S T [ 3.91 , 1.91 ] 2.05 [ 8.43 , 7.51 ] 1.12 25.5 [ 2.75 × 10 1 , 4.05 × 10 1 ] 1.10 × 10 1
V M A X   [ 3.36 , 1.19 ] 2.81 [ 7.29 , 3.17 ] 2.30 54.8 [ 4.40 × 10 1 , 6.22 × 10 1 ] 1.56 × 10 1
C - C L O H I V I S S L O W [ 2.97 , 8.22 × 10 1 ] 3.61 [ 5.69 , 1.92 ] 2.97 2.85 [ 1.31 × 10 2 , 1.32 × 10 2 ] 3.41 × 10 2
F A S T [ 3.08 , 9.69 × 10 1 ] 3.18 [ 5.94 , 2.00 ] 2.97 5.92 [ 6.06 × 10 2 , 3.51 × 10 2 ] 2.02 × 10 2
V M A X [ 3.51 , 1.13 ] 3.11 [ 6.98 , 1.90 ] 3.68 16.8 [ 4.06 × 10 1 , 1.96 × 10 1 ] 1.24 × 10 2
L O V I S S L O W [ 3.15 , 1.04 ] 3.03 [ 6.67 , 1.90 ] 3.50 10.3 [ 1.61 × 10 2 , 8.31 × 10 3 ] 1.82 × 10 2
F A S T [ 3.62 , 1.19 ] 3.03 [ 7.12 , 1.79 ] 3.97 23.6 [ 6.58 × 10 2 , 3.11 × 10 2 ] 1.17 × 10 2
V M A X [ 3.53 , 1.34 ] 2.64 [ 6.37 , 1.72 ] 3.70 57.7 [ 4.24 × 10 1 , 2.02 × 10 1 ] 1.24 × 10 2
P - P E R H I V I S S L O W [ 1.15 , 7.61 × 10 1 ] 1.51 [ 2.64 , 3.63 ] 7.28 × 10 1 1.10 [ 1.76 × 10 1 , 3.15 × 10 1 ] 1.28 × 10 1
F A S T [ 1.83 , 1.10 ] 1.67 [ 3.84 , 4.56 ] 8.41 × 10 1 3.52 [ 5.13 × 10 1 , 9.24 × 10 1 ] 1.27 × 10 1
V M A X [ 3.32 , 2.04 ] 1.63 [ 5.57 , 6.01 ] 9.26 × 10 1 16.0 [ 2.10 , 3.04 ] 1.20 × 10 1
L O V I S S L O W [ 1.90 , 1.26 ] 1.51 [ 4.00 , 4.94 ] 8.10 × 10 1 6.21 [ 9.83 × 10 2 , 1.64 × 10 1 ] 1.21 × 10 1
F A S T [ 3.81 , 2.42 ] 1.58 [ 6.01 , 7.38 ] 8.15 × 10 1 24.9 [ 2.95 × 10 1 , 3.57 × 10 1 ] 1.19 × 10 1
V M A X   [ 1.66 , 1.25 ] 1.33 [ 3.27 , 4.06 ] 8.05 × 10 1 27.1 [ 4.22 × 10 1 , 5.86 × 10 1 ] 1.56 × 10 1
C - P E R H I V I S S L O W [ 1.05 , 7.43 × 10 1 ] 1.41 [ 2.00 , 2.94 ] 6.82 × 10 1 1.00 [ 1.27 × 10 2 , 8.48 × 10 3 ] 1.08 × 10 2
F A S T [ 1.23 , 8.56 × 10 1 ] 1.43 [ 2.35 , 3.01 ] 7.80 × 10 1 2.35 [ 5.95 × 10 2 , 3.03 × 10 2 ] 5.93 × 10 3
V M A X [ 1.43 , 9.42 × 10 1 ] 1.52 [ 2.58 , 2.97 ] 8.68 × 10 1 6.86 [ 3.94 × 10 1 , 1.66 × 10 1 ] 3.32 × 10 3
L O V I S S L O W [ 1.38 , 9.12 × 10 1 ] 1.51 [ 2.59 , 3.01 ] 8.60 × 10 1 4.50 [ 1.59 × 10 2 , 7.12 × 10 3 ] 5.12 × 10 3
F A S T [ 1.49 , 9.60 × 10 1 ] 1.55 [ 2.55 , 2.89 ] 8.82 × 10 1 9.71 [ 6.30 × 10 2 , 2.57 × 10 2 ] 3.12 × 10 3
V M A X [ 1.72 , 1.07 ] 1.61 [ 2.38 , 2.45 ] 9.71 × 10 1 28.1 [ 4.13 × 10 1 , 1.93 × 10 1 ] 5.03 × 10 3
Table A10. M O T - I I O C L - 4 . P - C L O L O V I S F A S T   , P - P E R L O V I S V M A X   , and P - C L O L O V I S V M A X   replace c s 2 = 1 3 by c s 2 = 1 3 × 4 n , n = { 1 , 2 , 3 } , respectively. C - C L O L O V I S V M A X   ( c s 2 = 1 3 × 4 3 ) and C - P E R L O V I S V M A X   ( c s 2 = 1 3 × 4 2 ) are given for comparison with their stable counterparts of V M A X ( c s 2 = 1 3 ) .
Table A10. M O T - I I O C L - 4 . P - C L O L O V I S F A S T   , P - P E R L O V I S V M A X   , and P - C L O L O V I S V M A X   replace c s 2 = 1 3 by c s 2 = 1 3 × 4 n , n = { 1 , 2 , 3 } , respectively. C - C L O L O V I S V M A X   ( c s 2 = 1 3 × 4 3 ) and C - P E R L O V I S V M A X   ( c s 2 = 1 3 × 4 2 ) are given for comparison with their stable counterparts of V M A X ( c s 2 = 1 3 ) .
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
P - C L O H I V I S S L O W [ 3.39 , 1.04 ] 3.27 [ 7.90 , 4.41 ] 1.79 3.26 [ 1.88 × 10 1 , 2.83 × 10 1 ] 2.33 × 10 1
F A S T [ 3.91 , 1.38 ] 2.83 [ 9.13 , 5.23 ] 1.74 7.50 [ 5.7 × 10 1 , 8.52 × 10 1 ] 2.20 × 10 1
V M A X [ 4.88 , 2.31 ] 2.11 [ 11.5 , 7.64 ] 1.50 23.4 [ 2.24 , 3.94 ] 2.21 × 10 1
L O V I S S L O W [ 4.41 , 3.05 ] 1.44 [ 10.6 , 6.23 ] 1.70 14.4 [ 1.08 × 10 1 , 1.59 × 10 1 ] 1.93 × 10 1
F A S T   [ 4.15 , 1.78 ] 2.34 [ 9.28 , 5.95 ] 1.56 27.1 [ 1.92 × 10 1 , 3.44 × 10 1 ] 2.49 × 10 1
V M A X   [ 2.57 , 1.17 ] 2.19 [ 5.56 , 3.82 ] 1.45 41.9 [ 1.97 × 10 1 , 2.45 × 10 1 ] 4.36 × 10 1
C - C L O H I V I S S L O W [ 3.28 , 1.64 ] 2.00 [ 6.27 , 5.10 ] 1.23 3.14 [ 1.81 × 10 2 , 1.55 × 10 2 ] 5.06 × 10 2
F A S T [ 3.57 , 1.94 ] 1.84 [ 7.08 , 5.55 ] 1.28 6.85 [ 8.62 × 10 2 , 4.81 × 10 2 ] 3.53 × 10 2
V M A X [ 4.75 , 2.82 ] 1.69 [ 9.41 , 5.25 ] 1.79 22.8 [ 5.76 × 10 1 , 2.98 × 10 1 ] 2.87 × 10 2
L O V I S S L O W [ 4.37 , 2.33 ] 1.88 [ 8.77 , 5.50 ] 1.59 14.3 [ 2.27 × 10 2 , 1.20 × 10 2 ] 3.53 × 10 2
F A S T [ 4.98 , 3.20 ] 1.56 [ 9.69 , 5.06 ] 1.92 32.5 [ 9.30 × 10 2 , 4.79 × 10 2 ] 3.03 × 10 2
V M A X [ 5.00 , 4.01 ] 1.25 [ 8.58 , 5.83 ] 1.47 81.6 [ 6.19 × 10 1 , 3.47 × 10 1 ] 6.13 × 10 2
cf. V M A X V M A X   [ 2.75 , 1.42 ] 1.94 [ 4.91 , 4.06 ] 1.21 44.8 [ 8.78 × 10 2 , 6.14 × 10 2 ] 2.50 × 10 1
P - P E R H I V I S S L O W [ 1.48 , 1.43 ] 1.04 [ 3.48 , 5.32 ] 6.54 × 10 1 1.42 [ 1.78 × 10 1 , 2.75 × 10 1 ] 1.95 × 10 1
F A S T [ 2.31 , 1.59 ] 1.46 [ 5.46 , 5.66 ] 9.65 × 10 1 4.44 [ 5.66 × 10 1 , 8.43 × 10 1 ] 1.95 × 10 1
V M A X [ 3.74 , 1.95 ] 1.92 [ 8.16 , 8.02 ] 1.02 18.0 [ 2.21 , 3.81 ] 2.19 × 10 1
L O V I S S L O W [ 3.17 , 1.80 ] 1.77 [ 7.26 , 6.11 ] 1.19 10.4 [ 1.12 × 10 1 , 1.71 × 10 1 ] 2.05 × 10 1
F A S T [ 4.06 , 2.03 ] 2.00 [ 8.50 , 9.44 ] 9.01 × 10 1 26.5 [ 2.79 × 10 1 , 5.41 × 10 1 ] 2.30 × 10 1
V M A X   [ 2.13 , 1.86 ] 1.15 [ 4.62 , 5.11 ] 9.05 × 10 1 34.7 [ 4.19 × 10 1 , 6.82 × 10 1 ] 2.64 × 10 1
C - P E R H I V I S S L O W [ 1.55 , 1.55 ] 1.00 [ 2.95 , 4.93 ] 5.98 × 10 1 1.49 [ 1.73 × 10 2 , 1.27 × 10 2 ] 2.64 × 10 2
F A S T [ 1.99 , 1.70 ] 1.17 [ 3.78 , 5.08 ] 7.43 × 10 1 3.82 [ 8.64 × 10 2 , 5.30 × 10 2 ] 1.67 × 10 2
V M A X [ 2.46 , 1.88 ] 1.31 [ 4.32 , 4.99 ] 8.64 × 10 1 11.8 [ 5.85 × 10 1 , 3.17 × 10 1 ] 1.25 × 10 2
L O V I S S L O W [ 2.28 , 1.83 ] 1.25 [ 4.20 , 5.11 ] 8.23 × 10 1 7.45 [ 2.30 × 10 2 , 1.30 × 10 2 ] 1.34 × 10 2
F A S T [ 2.60 , 1.89 ] 1.38 [ 4.29 , 4.80 ] 8.94 × 10 1 17.0 [ 9.34 × 10 2 , 4.99 × 10 2 ] 1.28 × 10 2
V M A X [ 3.06 , 1.98 ] 1.55 [ 4.21 , 4.04 ] 1.04 50.0 [ 5.54 × 10 1 , 2.65 × 10 1 ] 2.19 × 10 2
cf. V M A X V M A X   [ 2.31 , 2.00 ] 1.15 [ 4.19 , 4.37 ] 9.59 × 10 1 37.7 [ 1.34 × 10 1 , 1.05 × 10 1 ] 8.06 × 10 2
Table A11. M O T - I I O C L - 5 ; C - C L O L O V I S   V M A X   applies c s 2 = 1 3 × 4 3 .
Table A11. M O T - I I O C L - 5 ; C - C L O L O V I S   V M A X   applies c s 2 = 1 3 × 4 3 .
G - C a s e ν ph v w [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max Re ( | v ¯ min | ) [ p min , p max ] [ mmHg ] δ M
C - C L O H I V I S S L O W [ 6.61 , 9.45 × 10 1 ] 7.00 [ 12.8 , 3.19 ] 4.02 6.35 [ 5.17 × 10 2 , 5.00 × 10 2 ] 1.75 × 10 1
F A S T [ 7.08 , 1.13 ] 6.28 [ 13.7 , 3.51 ] 3.91 13.6 [ 1.40 × 10 1 , 9.64 × 10 2 ] 1.30 × 10 1
V M A X [ 7.28 , 1.34 ] 5.42 [ 13.9 , 3.69 ] 3.75 34.9 [ 8.91 × 10 1 , 3.08 × 10 1 ] 1.36 × 10 1
L O V I S S L O W [ 7.13 , 1.21 ] 5.89 [ 13.8 , 3.68 ] 3.73 23.3 [ 3.51 × 10 2 , 1.43 × 10 2 ] 1.25 × 10 1
F A S T [ 7.46 , 1.54 ] 4.85 [ 14.0 , 3.64 ] 3.84 48.7 [ 1.44 × 10 1 , 4.44 × 10 2 ] 1.52 × 10 1
V M A X   [ 5.44 , 1.05 ] 5.17 [ 11.8 , 3.40 ] 3.46 88.8 [ 1.58 × 10 1 , 5.88 × 10 2 ] 7.07 × 10 1
C - P E R H I V I S S L O W [ 1.31 , 1.77 ] 7.41 × 10 1 [ 2.57 , 6.21 ] 4.13 × 10 1 1.26 [ 1.92 × 10 2 , 2.48 × 10 2 ] 1.68 × 10 2
F A S T [ 1.42 , 1.90 ] 7.48 × 10 1 [ 2.76 , 6.35 ] 4.35 × 10 1 2.73 [ 1.16 × 10 1 , 5.31 × 10 2 ] 5.98 × 10 3
V M A X [ 1.85 , 2.05 ] 9.05 × 10 1 [ 3.61 , 6.21 ] 5.82 × 10 1 8.90 [ 8.39 × 10 1 , 2.04 × 10 1 ] 4.85 × 10 3
L O V I S S L O W [ 1.67 , 1.98 ] 8.44 × 10 1 [ 3.26 , 6.32 ] 5.17 × 10 1 5.47 [ 3.22 × 10 2 , 8.65 × 10 3 ] 3.59 × 10 3
F A S T [ 1.98 , 2.11 ] 9.38 × 10 1 [ 3.83 , 6.07 ] 6.31 × 10 1 12.9 [ 1.36 × 10 1 , 3.12 × 10 2 ] 7.97 × 10 3
V M A X [ 2.23 , 2.33 ] 9.53 × 10 1 [ 4.01 , 5.50 ] 7.28 × 10 1 36.3 [ 8.65 × 10 1 , 2.48 × 10 1 ] 3.48 × 10 2

Appendix C.1. Motility Pattern MOT-II

Recall that when M O T - I I is operated with wave speed v w = 2 cm / s in Figure 4 left, then (a) segment n [ 2 , 10 ] reaches the most relaxed position o ( t ) = o r when t = ( 0.375 + ( n 2 ) ) s ; (b) it reaches the highest occlusion degree t f = 1 s later, (c) and performs there pause f during t f , p = 1 s ; (d) when segment n starts its back relaxation, its downstream neighbor n + 1 becomes most occluded; (e) the motility cycle over one segment lasts t w = 7 s against t w = 5 s with M O T - I ; and (f) the output is performed every 0.1 s . Velocity speeds S L O W and F A S T lengthen all these time intervals by factors of 5 and 5 / 2 , respectively. Figure A6, Figure A7, Figure A8, Figure A9, Figure A10, Figure A11, Figure A12 and Figure A13 illustrate the M O T - I I results with five occlusion scenarios and three wave speeds.
Figure A6. Case P - C L O : M O T - I I : O C L - 1 : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Figure A6. Case P - C L O : M O T - I I : O C L - 1 : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g0a6
Figure A7. Case P - P E R : M O T - I I : O C L - 2 : H I V I S : S L O W . Legende (right-middle) to (left-top) and (left-middle).
Figure A7. Case P - P E R : M O T - I I : O C L - 2 : H I V I S : S L O W . Legende (right-middle) to (left-top) and (left-middle).
Fluids 10 00022 g0a7
Figure A8. P - P E R : M O T - I I : O C L - 3 : H I V I S : V M A X . Legende (left-top) to (left-middle) and (right-middle).
Figure A8. P - P E R : M O T - I I : O C L - 3 : H I V I S : V M A X . Legende (left-top) to (left-middle) and (right-middle).
Fluids 10 00022 g0a8
Figure A9. Case P - P E R : M O T - I I : O C L - 3 : L O V I S : V M A X   (case ( L - V )   in Table 5 with c s 2 = 1 3 × 1 4 2 ). Legende (left-top) to (left-middle) and (right-middle).
Figure A9. Case P - P E R : M O T - I I : O C L - 3 : L O V I S : V M A X   (case ( L - V )   in Table 5 with c s 2 = 1 3 × 1 4 2 ). Legende (left-top) to (left-middle) and (right-middle).
Fluids 10 00022 g0a9
Figure A10. Case C - C L O : M O T - I I : O C L - 4 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Figure A10. Case C - C L O : M O T - I I : O C L - 4 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g0a10
Figure A11. Case P - C L O : M O T - I I : O C L - 4 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Figure A11. Case P - C L O : M O T - I I : O C L - 4 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g0a11
Figure A12. Case C - C L O : M O T - I I : O C L - 5 : L O V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Figure A12. Case C - C L O : M O T - I I : O C L - 5 : L O V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g0a12
Figure A13. Case C - P E R : M O T - I I : O C L - 5 : L O V I S : S L O W . Legende (left-top) to (left-middle) and (right-middle).
Figure A13. Case C - P E R : M O T - I I : O C L - 5 : L O V I S : S L O W . Legende (left-top) to (left-middle) and (right-middle).
Fluids 10 00022 g0a13

Appendix C.1.1. The MOT-II ∪ OCL-1

Table A7 addresses M O T - I I subjected to occlusion scenario O C L - 1 to compare it with motility M O T - I O C L - 1 from Table A2. According to the summary in Table 5, the largest velocity contrasts increase slightly from M O T - I to M O T - I I in closed pipes; in contrast, they decrease significantly in periodic pipes where, however, the largest retrograde amplitudes continue to dominate their antegrade counterparts for all six ν p h v w experiments. max ( R e ( | v ¯ m i n | ) ) increases in all cases from M O T - I to M O T - I I due to larger retrograde amplitudes by factors of about 2 and 1.5 , respectively, in closed and periodic pipes; the relative volume change min t δ V ( t ) also increases in magnitude with the same rate, as from 6.7 % to 12.7 % in periodic pipes, and very similar in closed tubes.
Figure A6 displays the case H I V I S S L O W in the “parabolic” closed pipe ( P - C L O ) to compare with its M O T - I counterpart in Figure 10. The segment-averaged profiles of pressure p ¯ ( t ) and velocity v ¯ ( t ) show several similar features with the two motility patterns: (i) in the last segment, p ¯ ( t ) oscillates before and after subsequent contractions across the colon, and it increases its amplitude with time; (ii) therefore, there is a permanent positive pressure gradient across the segments, which induces a retrograde staircase flow of increasing amplitude from the first to the last segment. In this, M O T - I I approximately doubles the magnitude of the average pressure and mean velocity oscillations that are displayed by M O T - I in Figure 10.
Locally, in s . 6 (third and fourth diagrams), v ¯ ( t ) changes direction from antegrade to retrograde during its local contraction and keeps this reflux direction until the last contraction in closed segment s . 10 . Peak pressure and velocity intervals, such as [ p m i n , p m a x ] and [ v ¯ m i n , v ¯ m a x ] , almost double their amplitudes compared to M O T - I , in agreement with the increase in the average amplitudes of two fields observed above. In contrast with M O T - I , M O T - I I maintains the values of mean and peak (retrograde) velocity amplitudes | v ¯ ( t ) | and | v m i n ( t ) | during the first half of the pause f interval of the local contraction (here, 1 2 t f , p = 2.5 s ); then, the negative-valued v m i n ( t ) increases rapidly until the downstream contraction, where it fails again together with p m i n . Hence, unlike with M O T - I where there is no pause between the contraction and the return to the neutral state, pause f induces staircase profiles v ¯ ( t ) and v m i n ( t ) , decreasing amplitudes downstream. P - C L O and C - C L O maintain distribution shapes across all six ν p h v w runs, where a decrease in retrograde velocity amplitude during pause f , due to the increase in the local pressure with downstream contraction, becomes more visible in faster flows.

Appendix C.1.2. MOT-II ∪ OCL-2

This second combination is addressed in Table A8 to compare it with M O T - I O C L - 2 in Table A3. These results show that M O T - I I only slightly increases the mean velocity contrast max ( R ( v ¯ ) ) compared to M O T - I O C L - 2 but significantly increases max ( R ( v ) ) whose values exceed now one in all regimes in the closed pipes, and with the fastest wave speed V M A X in periodic tubes. According to the summary of results in Table 5, the values of max ( R e ( | v ¯ m i n | ) ) , and hence the highest retrograde amplitudes, increase by a factor of about 2.5 in all four geometries compared to M O T - I . In turn, compared to the previous M O T - I I O C L - 1 scenario, the pressure intervals and velocity amplitudes increase noticeably due to faster occlusion speeds, but R ( v ¯ ) ) and R ( v ) remain smaller, similar to their M O T - I counterparts in Table A2 ( O C L - 1 ) and Table A3 ( O C L - 2 ).
Figure A7 and Figure 14 allow to compare M O T - I I to M O T - I in the case H I V I S S L O W in “parabolic” periodic pipes ( P - P E R ). These two motility models display qualitatively similar sequences of peak values with respect to occlusion events, but M O T - I I makes them sharper and of greater amplitude. Here, the two commonly observed antegrade peaks v m a x ( t ) occur just before the upstream and local contractions, while the retrograde peaks v m i n ( t ) have less disparate amplitudes and are clearly identified during sequential events of highest occlusion, where v ¯ ( t ) regularly changes sign with M O T - I I . The antegrade velocity peaks drop during the upstream and local pause f time interval, while the retrograde peaks only slightly reduce their amplitudes. O C L - 2 shows the same range of antegrade and retrograde velocity peak amplitude, max t v m a x ( t ) and | min t v m i n ( t ) | with both motility patterns, M O T - I & M O T - I I , consistent with the above analysis of their contrasts in periodic pipes.

Appendix C.1.3. MOT-II ∪ OCL-3

Table A9 addresses this third motility protocol M O T - I I to compare it with M O T - I O C L - 3 from Table A4. Recall that o f reaches, with O C L - 3 , the maximum available value o max in non-overlapping parabolic contours. Table 5 indicates that M O T - I I increases max ( R e ( | v ¯ m i n | ) ) by a factor of about 1.8 over M O T - I in all four geometries; the velocity contrasts also increase to max ( R ( v ¯ ) ) 4.4 and max ( R ( v ¯ ) ) 3.6 . These values are reached in the case H I V I S S L O W with closed pipes, “parabolic” and circular, respectively.
Figure A8 and Figure 17 allow us to compare M O T - I I and M O T - I in the case H I V I S V M A X in a “parabolic” periodic tube ( P - P E R ). Note that M O T - I I produces twice as large a relative volume change; M O T - I I averaged pressure and velocity distributions become sharper across the ten segments, with multiple retrograde velocity peaks during sequential pause f events; the amplitudes of pressure and velocity peaks then exceed M O T - I by a factor of two–three, and average velocity v ¯ ( t ) displays rapid oscillations due to the superposition of antegrade and retrograde sharp profiles. Finally, due to the strong oscillations of pressure and velocity extreme values at the highest contraction o f = o max , the “parabolic” pipes face a relatively high mass change, up to max ( | δ M | ) | 17 % | , while it stays below max ( | δ M | ) | 3.4 % | in circular pipes, as has already been observed with M O T - I O C L - 3 .
Figure A9 displays the case P - P E R L O V I S V M A X   ( ( L - V )   in Table 5), which is computed with a small sound speed in Equation (4a), as c s 2 = 1 3 × 1 4 2 , due to the first stability problem encountered. In central segment s . 6 , the retrograde and antegrade velocity profiles retain their H I V I S shapes displayed in Figure A8 but reduce the disparity of H I V I S retrograde to antegrade amplitudes by the factor of 1.5 ÷ 2 , while the pressure amplitude reduces by the significant factor of about 7. In contrast, C - P E R , where V M A X remains stable with c s 2 = 1 3 in both H I V I S and L O V I S fluids, retains very close characteristic values between these two viscosity values in Table A9.
A similar difference in pressure and velocity amplitudes in Table A9 between the H I V I S and L O V I S V M A X   results obtained with two distinct sound speeds in the case P - P E R   F A S T leads us to conclude that a decrease in c s regularly leads to a reduction in the amplitude of these two hydrodynamic fields. Hence, at a given compressibility level using c s 2 = 1 3 , one would expect to obtain a higher retrograde velocity amplitude and the Reynolds number than reported here using small values of c s , as R e ( | v ¯ m i n | ) 27.1 with ( L - V )   in Table 5.

Appendix C.1.4. MOT-II ∪ OCL-4

This fourth motility protocol is addressed in Table A10. Recall that O C L - 3 and O C L - 4 operate with o f = o max but O C L - 4 prescribes a smaller value for o r , and then this occlusion scenario produces faster relaxation and contraction speeds. Like M O T - I , M O T - I I displays similar pressure and velocity shapes with O C L - 3 and O C L - 4 at the same sound velocity, but O C L - 4 increases peak interval [ v m i n , v m a x ] , R e ( | v ¯ m i n | ) and the amplitude of relative mass change | δ M | , especially in the fastest regimes; at the same time, O C L - 4 reduces velocity contrasts R ( v ¯ ) and R ( v ) due to its higher antegrade amplitudes.
Figure A10 and Figure 18 allow to compare M O T - I I to M O T - I in the case H I V I S V M A X in closed circular pipe ( C - C L O ); they show that the averaged pressure and velocity amplitudes increase by a factor of two. On the other hand, in s . 6 , the highest retrograde, both averaged and peak velocity amplitudes, which occur rapidly during the highest contraction events with M O T - I , are shifted towards the end of the pause f interval, when the downstream neighbor becomes obstructed and it quickly repels the flow upstream due to a very strong local pressure drop; accordingly, peak antegrade distribution v m a x ( t ) decreases its amplitude.
Figure A11 also addresses the case H I V I S V M A X but in the closed “parabolic” pipe ( P - C L O ), where the averaged pressure amplitude further increases by a factor of two compared to the closed circular pipe in Figure A10, and by a factor of 2 ÷ 3 compared to the M O T - I counterpart. Retrograde peak distribution v m i n ( t ) doubles its amplitude in s . 6 ; the most significant reflux then happens very quickly as soon as the highest degree of occlusion is reached, particularly upstream, locally and downstream, while its secondary peaks occur sequentially halfway through the time interval of pause f . Here, the magnitude of the relative mass change | δ M | increases with respect to M O T - I , but it remains much smaller in the circular pipe of Figure A10 than in the “parabolic” pipe of Figure A11.
According to Table A10, the case L O V I S V M A X is stabilized in “parabolic” pipes with even smaller values of sound speed, as c s 2 = 1 3 × 1 2 6 ( P - C L O ) and c s 2 = 1 3 × 1 2 4 ( P - P E R ). On the other hand, Table A10 allows us to compare the stable results of L O V I S V M A X in circular pipes between c s 2 = 1 3 and c s 2 = 1 3 × 4 n , with n = 3 ( C - C L O ) and n = 2 ( C - P E R ), respectively. These results confirm that the higher compressibility due to a smaller value of c s produces (a) reduction in the pressure interval [ p m i n , p m a x ] by factors of [ 7 , 5 ] ( n = 6 ) and [ 4 , 2 ] ( n = 4 ), respectively; (b) reduction in R e ( | v ¯ m i n | ) by a factor of about 1.8 ( n = 6 ) and 1.3 ( n = 4 ); (c) increase in | δ M | by a factor of 4 ( n = { 4 , 6 } ), while peak values [ v m i n , v m a x ] and R ( v ) remain almost independent of sound speed value. This test confirms that even if the main features of a flow pattern stay relatively correct with very small values of sound velocity, the obtained pressure intervals, mean velocity contrasts, and the R e ( | v ¯ m i n | ) range are smaller than would be expected from incompressible fluid.

Appendix C.1.5. The MOT-II ∪ OCL-5

The fifth, most occluding scenario is covered in Table A11; here, since o f 0.938 > o max , O C L - 5 is limited to circular pipes. According to Table 5, | min t δ V ( t ) | increases from approximately |−16%| to |−31%| between M O T - I and M O T - I I , and | δ M | grows even faster. Furthermore, the contrast of averaged velocities R ( v ¯ ) increases for all six runs in closed pipes ( C - C L O ) compared with M O T - I O C L - 5 in Table A6, and especially with S L O W propagation speed where R ( v ¯ ) reaches the highest values on all simulations, as R ( v ¯ ) { 6 , 7 } in L O V I S and H I V I S fluids, respectively. In contrast, all velocity ratios decrease in periodic pipes ( C - P E R ) by a factor of about two, and then they become less than one, meaning that the highest antegrade velocity amplitude dominates the retrograde one. This is explained by the reduction in the amplitude of retrograde velocity, where max ( R e ( | v ¯ m i n | ) ) then reduces from 50 ( O C L - 4 ) to about 36 ( O C L - 5 ).
Figure A11 addresses the case L O V I S S L O W in closed circular pipe ( C - C L O ) to compare it with M O T - I in Figure 20. The two figures show that the two largest antegrade peaks v m a x ( t ) due to upstream and local relaxations retain their amplitudes. In contrast, mean-velocity distribution v ¯ ( t ) is now maintained retrograde in s . 6 from an upstream contraction to a final contraction; the five peaks of v m i n ( t ) (Figure 20) become replaced by the five “fingers” which lengthen and shorten during the sequential intervals of pause f ; their amplitude increases spectacularly from | v m i n | 8 v w ( M O T - I ) to | v m i n | 14 v w ( M O T - I I ); | v ¯ | and R e ( | v ¯ m i n | ) then increase by a factor of approximately 1.6 and | δ M | grows from | 1.36 % | to | 12.5 % | .
Figure A13 shows that the same protocol L O V I S S L O W produces a quite different flow pattern in a circular periodic pipe ( C - P E R ). Here, unlike in closed pipes, the two largest antegrade peaks v m a x ( t ) due to upstream and local relaxations increase their amplitudes from max ( v m a x ) 5 v w ( M O T - I ) to max ( v m a x ) 7 v w ( M O T - I I ); on the other hand, the retrograde amplitudes decrease from max ( | v m i n | ) 5 v w ( M O T - I ) to max ( | v m i n | ) 3 v w ( M O T - I I ); this explains R ( v ) and R ( v ¯ ) reduction in Table A11 by a factor of two. Furthermore, for the first time, max t ( v ¯ ) exceeds its retrograde counterpart | min t ( v ¯ ) | for all six ν p h v w configurations, offering R ( v ¯ ) < 1 in a whole series of computations.
Finally, the fastest case C-CLO LOVIS VMAX is the first configuration where c s 2 = 1 3 loses robustness in the circular pipes. It stabilizes with c s 2 = 1 3 × 4 3 , where it reaches R e ( | v ¯ m i n | ) = max ( R e ( | v ¯ m i n | ) ) 89 , and the value of the Reynolds number should become even higher in an incompressible flow according to numerical experiments discussed previously in Table A10.

Appendix C.2. LBM-DB Results in Table A2Table A11

Table A2, Table A3, Table A4, Table A5 and Table A6 provide the results of the five occlusion scenarios from Table 4 applied with the M O T - I motility pattern in Section 3.5; Table A7, Table A8, Table A9, Table A10, Table A11 address similar experiments for M O T - I I from Appendix C.1. These results are compared in Section 3.6.

References

  1. Corsetti, M.; Costa, M.; Bassotti, G.; Bharucha, A.E.; Borrelli, O.; Dinning, P.; Di Lorenzo, C.; Huizinga, J.D.; Jimenez, M.; Rao, S.; et al. First translational consensus on terminology and definitions of colonic motility in animals and humans studied by manometric and other techniques. Nat. Rev. Gastroenterol. Hepatol. 2019, 16, 559–579. [Google Scholar] [CrossRef] [PubMed]
  2. Schulze, K.S. The imaging and modelling of the physical processes involved in digestion and absorption. Acta Physiol. 2014, 213, 394–405. [Google Scholar] [CrossRef] [PubMed]
  3. Palmada, N.; Cater, J.E.; Cheng, L.K.; Suresh, V. Experimental and computational studies of peristaltic flow in a duodenal model. Fluids 2022, 7, 40. [Google Scholar] [CrossRef]
  4. Zhang, Y.; Wu, P.; Jeantet, R.; Dupont, D.; Delaplace, G.; Chen, X.D.; Xiao, J. How motility can enhance mass transfer and absorption in the duodenum: Taking the structure of the villi into account. Chem. Eng. Sci. 2020, 213, 115406. [Google Scholar] [CrossRef]
  5. Lim, Y.F.; de Loubens, C.; Love, R.J.; Lentle, R.G.; Janssen, P.W.M. Flow and mixing by small intestine villi. Food Funct. 2015, 6, 1787. [Google Scholar] [CrossRef] [PubMed]
  6. Cremer, J.; Arnoldini, M.; Hwa, T. Effect of water flow and chemical environment on microbiota growth and composition in the human colon. Proc. Natl. Acad. Sci. USA 2017, 114, 6438–6443. [Google Scholar] [CrossRef]
  7. Pal, A.; Indireshkumar, K.; Schwizer, W.; Abrahamsson, B.; Fried, M.; Brasseur, J.G. Gastric flow and mixing studied using computer simulation. Proc. R. Soc. B Biol. Sci. 2004, 271, 2587–2594. [Google Scholar] [CrossRef] [PubMed]
  8. Zha, J.; Zou, S.; Xao, J.; Liu, X.; Delaplace, G.; Jeantet, R.; Dupont, D.; Wu, P.; Chen, X.D.; Xiao, J. The role of circular folds in mixing intensification in the small intestine: A numerical study. Chem. Eng. Sci. 2021, 229, 116079. [Google Scholar] [CrossRef]
  9. Stamatopoulos, K.; Karandikar, S.; Goldstein, M.; O’Farrell, C.; Marciani, L.; Sulaiman, S.; Hoad, C.L.; Simmons, M.J.H.; Batchelor, H.K. Dynamic colon model (DCM): A cine-MRI informed biorelevant in vitro model of the human proximal large intestine characterized by positron imaging techniques. Pharmaceutics 2020, 12, 659. [Google Scholar] [CrossRef]
  10. O’Farrell, C.; Stamatopoulos, K.; Simmons, M.; Batchelor, H. In vitro models to evaluate ingestible devices: Present status and current trends. Adv. Drug Deliver. Rev. 2021, 178, 113924. [Google Scholar] [CrossRef]
  11. Schütt, M.; Stamatopoulos, K.; Batchelor, H.K.; Simmons, M.J.H.; Alexadis, A. Modelling and Simulation of the Drug Release from a Solid Dosage Form in the Human Ascending Colon: The Influence of Different Motility Patterns and Fluid Viscosities. Pharmaceutics 2021, 13, 859. [Google Scholar] [CrossRef] [PubMed]
  12. Wei, L.Y.; Yu, Y.-J.; Fei, F.; Zheng, M.-Y.; Zhang, S.W. High-resolution colonic manometry and its clinical application in patients with colonic dysmotility: A review. World J. Clin. Cases 2019, 7, 2675–2686. [Google Scholar]
  13. Huizinga, J.D.; Hussain, A.; Chen, J.-H. Interstitial cells of Cajal and human colon motility in health and disease. Am. J. Physiol. Gastrointest. Liver Physiol. 2021, 321, G552–G575. [Google Scholar] [CrossRef] [PubMed]
  14. Spencer, N.J.; Kyloh, M.; Wattchow, D.A.; Thomas, A.; Sia, T.C.; Brookes, S.J.; Nicholas, S.J. Characterization of motor patterns in isolated human colon: Are there differences in patients with slow-transit constipation? Am. J. Physiol. Gastrointest. Liver Physiol. 2012, 302, G34–G43. [Google Scholar] [CrossRef] [PubMed]
  15. Ahmad, F.; de Loubens, C.; Dubreuil, A.; Faucheron, J.L.; Tanguy, S. Towards an assessment of rectal function by coupling X-ray defecography and fluid mechanical modelling. Annu. Int. Conf. IEEE Eng. Med. Biol. Soc. 2022, 2022, 4962–4965. [Google Scholar]
  16. Trzpis, M.; Sun, G.; Chen, J.-H.; Huizinga, J.-D.; Broens, P. Novel insights into physiological mechanisms underlying fecal continence. Am. J. Physiol. Gastrointest. Liver Physiol. 2023, 324, G1–G9. [Google Scholar] [CrossRef] [PubMed]
  17. Dinning, P.G.; Wiklendt, L.; Maslen, L.; Gibbins, I.; Patton, V.; Arkwright, J.W.; Lubowski, D.Z.; O’Grady, G.; Bampton, P.A.; Brookes, S.J.; et al. Quantification of in vivo colonic motor patterns in healthy humans before and after a meal revealed by high-resolution fiber-optic manometry. Neurogastroenterol. Motil. 2014, 26, 1443–1457. [Google Scholar] [CrossRef] [PubMed]
  18. Kuizenga, M.H.; Sia, T.C.; Dodds, K.N.; Wiklendt, L.; Arkwright, J.W.; Thomas, A.; Brookes, S.J.; Spencer, N.J.; Wattchow, D.A.; Dinning, P.G.; et al. Neurally mediated propagating discrete clustered contractions superimposed on myogenic ripples in ex vivo segments of human ileum. Am. J. Physiol. Gastrointest. Liver Physiol. 2015, 308, G1–G11. [Google Scholar] [CrossRef]
  19. Lin, A.Y.; Du, P.; Dinning, P.G.; Arkwright, J.W.; Kamp, J.P.; Cheng, L.K.; Bissett, I.P.; O’Grady, G. High-resolution anatomic correlation of cyclic motor patterns in the human colon: Evidence of a rectosigmoid brake. Am. J. Physiol. Gastrointest. Liver Physiol. 2017, 312, G508–G515. [Google Scholar] [CrossRef]
  20. Dinning, P.G. A new understanding of the physiology and pathophysiology of colonic motility? Neurogastroenterol. Motil. 2018, 30, e13395. [Google Scholar] [CrossRef] [PubMed]
  21. Milkova, N.; Parsons, S.P.; Ratcliffe, E.; Huizinga, J.D. On the nature of high-amplitude propagating pressure waves in the human colon. Am. J. Physiol. Gastrointest. Liver Physiol. 2020, 318, G646–G660. [Google Scholar] [CrossRef] [PubMed]
  22. Li, Y.; Xu, R.; Xiu, H.; Kong, F. Development of a small intestinal simulator to assess the intestinal mixing and transit as affected by digesta viscosity. Innovat. Food Sci. Emerg. Technol. 2022, 82, 103202. [Google Scholar] [CrossRef]
  23. Huizinga, J.D.; Parsons, S.P.; Chen, J.-H.; Pawelka, A.; Pistilli, M.; Li, C.; Yu, Y.; Ye, P.; Liu, Q.; Tong, M.; et al. Motor patterns of the small intestine explained by phase-amplitude coupling of two pacemaker activities: The critical importance of propagation velocity. Am. J. Physiol. Cell Physiol. 2015, 309, C403–C414. [Google Scholar] [CrossRef] [PubMed]
  24. Fullard, L.; Lammers, W.; Wakec, G.C.; Ferrua, M.J. Propagating longitudinal contractions in the ileum of the rabbit—Efficiency of advective mixing. Food Funct. 2014, 5, 2731–2742. [Google Scholar] [CrossRef] [PubMed]
  25. Fullard, L.; Lammers, W.J.; Ferrua, M.J. Advective mixing due to longitudinal and segmental contractions in the ileum of the rabbit. J. Food Eng. 2015, 160, 1–10. [Google Scholar] [CrossRef]
  26. Wang, Y.; Brasseur, J.G.; Banco, G.G.; Webb, A.G.; Ailiani, A.C.; Neuberger, T. A multi-scale Lattice Boltzmann model of macro-to-micro scale transport, with applications to gut function. Philos. Trans. R. Soc. A Math. Phys. Eng. Sci. 2010, 368, 2683–2880. [Google Scholar]
  27. Wang, Y.; Brasseur, J.G. Three-dimensional mechanisms of macro-to-micro-scale transport and absorption enhancement by gut villi motions. Phys. Rev. E 2017, 95, 062412. [Google Scholar] [CrossRef] [PubMed]
  28. Agbesi, R.A.; Chevalier, N. Flow and mixing induced by single, colinear, and colliding contractile waves in the intestine. Phys. Rev. Fluids 2022, 7, 043101. [Google Scholar] [CrossRef]
  29. Palmada, N.; Cater, J.E.; Cheng, L.K.; Suresh, V. Anatomically realistic computational model of flow and mixing in the human duodenum. Phys. Fluids 2023, 35, 011907. [Google Scholar] [CrossRef]
  30. Pistilli, M.J. Mechanisms Underlying Rhythmic Activities in the Gastrointestinal Tract. Master’s Thesis, McMaster University, Hamilton, ON, Canada, 2012. [Google Scholar]
  31. Pritchard, S.E.; Paul, J.; Major, G.; Marciani, L.; Gowland, P.A.; Spiller, R.C.; Hoad, C.L. Assessment of motion of colonic contents in the human colon using MRI tagging. Neurogastroenterol. Motil. 2017, 29, e13091. [Google Scholar] [CrossRef]
  32. de Loubens, C.; Lentle, R.G.; Love, R.J.; Hulls, C.; Janssen, P.W. Fluid mechanical consequences of pendular activity, segmentation and pyloric outflow in the proximal duodenum of the rat and the guinea pig. J. R. Soc. Interface 2013, 10, 20130027. [Google Scholar] [CrossRef]
  33. de Loubens, C.; Lentle, R.G.; Hulls, C.; Janssen, P.W.; Love, R.J.; Chamber, J.P. Characterisation of Mixing in the Proximal Duodenum of the Rat during Longitudinal Contractions and Comparison with a Fluid Mechanical Model Based on Spatiotemporal Motility Data. PLoS ONE 2014, 9, e95000. [Google Scholar] [CrossRef]
  34. Canon, W.B. The movements of the intestines studied by means of the Röntgen rays. J. Med. Res. 1902, 7, 72–75. [Google Scholar] [CrossRef]
  35. Womack, W.A.; Kvietys, P.R.; Granger, D.N. Villious Motility. In Handbook of Physiology, 2nd ed.; Schultz, S.G., Wood, J.D., Rauner, B.B., Eds.; American Physiological Society: Bethesda, MD, USA, 1989; pp. 975–986. [Google Scholar]
  36. Codutti, A.; Cremer, J.; Alim, K. Changing Flows Balance Nutrient Absorption and Bacterial Growth along the Gut. Phys. Rev. Lett. 2022, 129, 138101. [Google Scholar] [CrossRef] [PubMed]
  37. Lammers, W.J. Spatial and temporal coupling between slow waves and pendular contractions. Am. J. Physiol.-Gastrointest. Liver Physiol. 2005, 289, G898–G903. [Google Scholar] [CrossRef] [PubMed]
  38. Lee, H.-T.; Hennig, G.W.; Fleming, N.W.; Keef, K.D.; Spencer, N.J.; Ward, S.M.; Sanders, K.M.; Smith, T.K. The mechanism and spread of pacemarker activity through myenteric interstitial cells of Cajal in human small intestine. Gastroentorology 2007, 132, 1852–1865. [Google Scholar] [CrossRef] [PubMed]
  39. Melepattu, M.P.; de Loubens, C. Steady streaming flow induced by active biological microstructures; application to small intestine villi. Phys. Fluids 2022, 34, 061905. [Google Scholar] [CrossRef]
  40. Higuera, F.; Succi, S.; Benzi, R. Lattice Gas Dynamics with Enhanced Collisions. Europhys. Lett. 1989, 9, 345. [Google Scholar] [CrossRef]
  41. Qian, Y.; d’Humières, D.; Lallemand, P. Lattice BGK models for Navier-Stokes equation, Europhys. Lett. 1992, 17, 479. [Google Scholar]
  42. Lallemand, P.; Luo, L.-S. Theory of the lattice Boltzmann method: Dispersion, dissipation, isotropy, Galilean invariance, and stability. Phys. Rev. E. 2000, 61, 6546. [Google Scholar] [CrossRef]
  43. d’Humières, D.; Ginzburg, I.; Krafczyk, M.; Lallemand, P.; Luo, L.-S. Multiple-relaxation-time lattice Boltzmann models in three dimensions. Philos. Trans. R. Soc. A 2002, 360, 437. [Google Scholar] [CrossRef] [PubMed]
  44. Ginzburg, I.; Verhaeghe, F.; d’Humières, D. Two-relaxation-time Lattice Boltzmann scheme: About parametrization, velocity, pressure and mixed boundary conditions. Commun. Comput. Phys. 2008, 3, 427–478. [Google Scholar]
  45. d’Humières, D.; Ginzburg, I. Viscosity independent numerical errors for Lattice Boltzmann models: From recurrence equations to “magic” collision numbers. Comput. Math. Appl. 2009, 58, 823. [Google Scholar] [CrossRef]
  46. Silva, G.; Semiao, V. Truncation errors and the rotational in-variance of three-dimensional lattice models in the lattice Boltzmann method. J. Comput. Phys. 2014, 269, 259. [Google Scholar] [CrossRef]
  47. Chai, Z.; Shi, B. Multiple-relaxation-time lattice Boltzmann method for the Navier-Stokes and nonlinear convection-diffusion equations: Modeling, analysis, and elements. Phys. Rev. E 2020, 102, 023306. [Google Scholar] [CrossRef] [PubMed]
  48. Ginzburg, I. Steady-state Two-relaxation-time Lattice Boltzmann formulation for transport and flow, closed with the compact multireflection boundary and interface-conjugate schemes. J. Comp. Sci. 2021, 54, 101215. [Google Scholar] [CrossRef]
  49. Christensen, J. Gross and Microscopic Anatomy of the Large Intestine. Chapter 2 in The Large Intestine: Physiology, Pathophysiology, and Disease; Raven Press: New York, NY, USA, 1991. [Google Scholar]
  50. Quan, X.; Yang, Z.; Xue, M.; Chen, J.-H.; Huizinga, J.D. Relationships between motor patterns and intraluminal pressure in the 3-taeniated proximal colon of the rabit. Nat. Sci. Rep. 2017, 7, 42293. [Google Scholar]
  51. Chen, J.-H.; Yang, Z.; Yu, Z.Y.; Huizinga, J.D. Haustral boundary contractions in the proximal 3-taeniated colon of the rabbit. Am. J. Physiol. 2016, 310, G181–G192. [Google Scholar] [CrossRef]
  52. High-definition spatiotemporal mapping of contractile activity in the isolated proximal colon of the rabbit. J. Comp. Physiol. B 2008, 178, 257–268. [CrossRef]
  53. Huizinga, J.D.; Pervez, M.; Nirmalathasan, S.; Chen, J.-H. Characterization of haustral activity in the human colon. Am. J. Physiol. Gastrointest. Liver Physiol. 2021, 320, G1067–G1080. [Google Scholar] [CrossRef] [PubMed]
  54. Blausen.com staff. Medical gallery of Blausen Medical 2014. WikiJ. Med. 2014, 1. ISSN 2002-4436. [Google Scholar] [CrossRef]
  55. Spencer, N.J.; Dinning, P.G.; Brookes, S.J.; Costa, M. Insights into the mechanisms underlying colonic motor patterns. J. Physiol. 2016, 594, 4099–4116. [Google Scholar] [CrossRef] [PubMed]
  56. Chen, J.H.; Yu, Y.; Yang, Z.; Yu, W.Z.; Chen, W.L.; Yu, H.; Kim, M.J.-M.; Huang, M.; Tan, S.; Luo, H.; et al. Intraluminal pressure patterns in the human colon assessed by high-resolution manometry. Sci. Rep. 2017, 7, 41436. [Google Scholar] [CrossRef] [PubMed]
  57. Pervez, M.; Ratcliffe, E.; Parsons, S.P.; Chen, J.-H.; Huizinga, J.D. The cyclic motor patterns in the human colon. Neurogasterinterol. Motil. 2020, 32, e13807. [Google Scholar] [CrossRef] [PubMed]
  58. Jaung, R.; Varghese, C.; Lin, A.Y.; Paskaranandavadivel, N.; Du, P.; Rowbotham, D.; Dinning, P.; O’Grady, G.; Bissett, I. High Resolution Colonic Manometry Pressure Profiles Are Similar in Asymptomatic Diverticulosis and Controls. Digest. Dis. Sci. 2021, 66, 832–842. [Google Scholar] [CrossRef] [PubMed]
  59. Cremer, J.; Segota, I.; Yang, C.-Y.; Arnoldini, M.; Sauls, J.T.; Zhang, Z.; Gutierrez, E.; Groisman, A.; Hwa, T. Effect of flow and peristaltic mixing on bacterial growth in a gut-like channel. Proc. Natl. Acad. Sci. USA 2016, 113, 11414–11419. [Google Scholar] [CrossRef] [PubMed]
  60. Labarthe, S.; Polizzi, B.; Phan, T.; Goudon, T.; Ribot, T.; Laroche, T. A mathematical model to investigate the key drivers of the biogeography of the colon microbiota. J. Theor. Biol. 2019, 462, 552–581. [Google Scholar] [CrossRef] [PubMed]
  61. Arkwright, J.W.; Dickson, A.; Maunder, S.A.; Blenman, N.G.; Lim, J.; O’Grady, G.; Archer, R.; Costa, M.; Spencer, N.J.; Brookes, S.J.; et al. The effect of luminal content and rate of occlusion on the interpretation of colonic manometry. Neurogastroenterol. Motil. Off. J. Eur. Gastrointest. Motil. Soc. 2013, 25, e52–e59. [Google Scholar] [CrossRef]
  62. Kamaltdinov, M.R.; Kuchumov, A.G.; Sadeghy, K. Chyme flow numerical modelling in the colon in normal conditions and in disorders. Ser. Biomech. 2022, 36, 105–112. [Google Scholar]
  63. Stamatopoulos, K.; Batchelor, H.K.; Simmons, M.J. Dissolution profile of theophylline modified release tablets, using a biorelevant Dynamic Colon Model (DCM). Eur. J. Pharm. Biopharm. 2016, 108, 9–17. [Google Scholar] [CrossRef] [PubMed]
  64. O’Farrell, C.; Hoad, C.L.; Stamatopoulos, K.; Marciani, L.; Sulaiman, S.; Simmons, M.J.H.; Batchelor, H.K. Luminal Fluid Motion Inside an In Vitro Dissolution Model of the Human Ascending Colon Assessed Using Magnetic Resonance Imaging. Pharmaceutics 2021, 13, 1545. [Google Scholar] [CrossRef] [PubMed]
  65. Sinnot, M.D.; Cleary, P.W.; Arkwright, J.W.; Dinning, P.G. Investigating the relationships between peristaltic contraction and fluid transport in the human colon using Smoothed Particle Hydrodynamics. Comput. Biol. Med. 2012, 42, 492–503. [Google Scholar] [CrossRef] [PubMed]
  66. Alexiadis, A.; Stamatopoulos, K.; Wen, W.; Batchelor, H.K.; Batchelor, H.K.; Bakalis, S.; Barigou, M.; Simmons, M.J.H. Using discrete multi-physics for detailed exploration of hydrodynamics in an in vitro colon system. Comput. Biol. Med. 2017, 81, 188–198. [Google Scholar] [CrossRef] [PubMed]
  67. Schütt, M.; Stamatopoulos, K.; Simmons, M.J.H.; Batchelor, H.K.; Alexadis, A. Modelling and simulation of the hydrodynamics and mixing profiles in the human proximal colon using Discrete Multiphysics. Comput. Biol. Med. 2020, 121, 103819. [Google Scholar] [CrossRef] [PubMed]
  68. Schütt, M.; O’Farrell, C.; Stamatopoulos, K.; Hoad, C.L.; Marciani, L.; Sulaiman, S.; Simmons, M.J.H.; Batchelor, H.K.; Alexadis, A. Simulating the Hydrodynamic Conditions of the Human Ascending Colon: A Digital Twin of the Dynamic Colon Model. Pharmaceutics 2022, 14, 1402. [Google Scholar] [CrossRef] [PubMed]
  69. Lätt, J.; Malaspinas, O.; Kontaxakis, D.; Parmigiani, A.; Lagrava, D.; Brogi, F.; Belgacem, M.B.; Thorimbert, Y.; Leclaire, S.; Li, S.; et al. Palabos: Parallel Lattice Boltzmann Solver. Comput. Math. Appl. 2021, 81, 334–350. [Google Scholar] [CrossRef]
  70. Ginzburg, I.; Silva, G.; Talon, L. Analysis and improvement of Brinkman lattice Boltzmann schemes: Bulk, boundary, interface. Similarity and distinctness with finite elements in heterogeneous porous media. Phys. Rev. E. 2015, 91, 023307-32. [Google Scholar] [CrossRef] [PubMed]
  71. Ouared, R.; Chopard, B.; Stahl, B.; Rfenacht, D.; Yilmaz, H.; Courbebaisse, G. Thrombosis modeling in intracranial aneurysms: A lattice Boltzmann numerical algorithm. Comput. Phys. Commun. 2008, 179, 128–131. [Google Scholar] [CrossRef]
  72. Zhang, Y. Hemodynamic Investigation and Thrombosis Modeling of Intracranial Aneurysms. Ph. D. Dissertation, INSA, Lyon, France, 2015. [Google Scholar]
  73. Zhang, Y.; Wang, Y.; Kao, E.; Flórez-Valencia, L.; Courbebaisse, G. Towards optimal flow diverter porosity for the treatment of intracranial aneurysm. J. Biomech. 2019, 82, 20–27. [Google Scholar] [CrossRef] [PubMed]
  74. Noël, R.; Ge, F.; Zhang, Y.; Navarro, L.; Courbebaisse, G. Lattice Boltzmann Method for Modelling of Biological Phenomena. In Proceedings of the 25th European Signal Processing Conference (EUSIPCO), Kos Island, Greece, 28 August–2 September 2017. [Google Scholar]
  75. Krüger, T. Computer Simulation Study of Collective Phenomena in Dense Suspensions of Red Blood Cells Under Shear; Springer Spektrum|Springer Fachmedien Wiesbaden GmbH: Wiesbaden, Germany, 2012. [Google Scholar]
  76. Chopard, B.; de Sousa, D.R.; Latt, J.; Dubois, F.; Yourassowsky, C.; Van Antwerpen, P.; Eker, O.; Vanhamme, L.; Perez-Morga, D.; Courbebaisse, G. A physical description of the adhesion and aggregation of platelets. R. Soc. Open Sci. 2017, 4, 170219. [Google Scholar] [CrossRef]
  77. Kotsalos, C. Digital Blood. Ph. D. Dissertation, Université de Genève, Geneva, Switzerland, 2020. [Google Scholar]
  78. Ginzburg, I.; d’Humières, D. Multi-reflection boundary conditions for lattice Boltzmann models. Phys. Rev. E 2003, 68, 066614-1-30. [Google Scholar]
  79. Ginzburg, I.; Silva, G.; Marson, F.; Chopard, B.; Latt, J. Unified directional parabolic-accurate lattice Boltzmann boundary schemes for grid-rotated narrow gaps and curved walls in creeping and inertial fluid flows. Phys. Rev. E 2023, 107, 025303-49. [Google Scholar] [CrossRef] [PubMed]
  80. Ginzbourg, I.; d’Humières, D. Local Second-Order Boundary method for Lattice Boltzmann models. J. Stat. Phys. 1996, 84, 927. [Google Scholar] [CrossRef]
  81. Khirevich, S.; Ginzburg, I.; Tallarek, U. Coarse- and fine-grid numerical behavior of MRT/TRT Lattice-Boltzmann schemes in regular and random sphere packings. J. Comput. Phys. 2014, 281, 708. [Google Scholar] [CrossRef]
  82. Ginzburg, I.; Silva, G. Mass-balance and locality vs accuracy with the new boundary and interface-conjugate approaches in advection-diffusion lattice Boltzmann method. Phys. Fluids 2021, 33, 057104. [Google Scholar] [CrossRef]
  83. Yu, H.; Chen, X.; Wang, Z.; Deep, D.; Lima, E.; Zhao, Y.; Teague, S.D. Mass-conserved volumetric lattice Boltzmann method for complex flows with fully moving boundaries. Phys. Rev. E 2014, 89, 063304. [Google Scholar] [CrossRef]
  84. Feng, Z.G.; Michaelides, E.E. The immersed boundary-lattice Boltzmann method for solving fluid–particles interaction problems. J. Comput. Phys. 2004, 195, 602. [Google Scholar] [CrossRef]
  85. Pepona, M.; Favier, J. A coupled Immersed Boundary—Lattice Boltzmann method for incompressible flows through moving porous media. J. Comput. Phys. 2016, 321, 1170. [Google Scholar] [CrossRef]
  86. Seta, T.; Rojas, R.; Hayashi, K.; Tomiyama, A. Implicit-correction-based immersed boundary-lattice Boltzmann method with two-relaxation time. Phys. Rev. E 2014, 89, 023307. [Google Scholar] [CrossRef]
  87. Wang, S.; Yan, H.; Cai, Y.; Pan, G.; Zhang, G.; Song, D. Accurate simulations of moving flexible objects with an improved immersed boundary-lattice Boltzmann method. Phys. Fluids 2024, 36, 11. [Google Scholar]
  88. Thorimbert, Y.; Marson, F.; Parmigiani, A.; Chopard, B.; Latt, J. Lattice Boltzmann simulation of dense rigid spherical particle suspensions using immersed boundary method. Comput. Fluids 2018, 166, 286–294. [Google Scholar] [CrossRef]
  89. Puleo, S.; Pasta, S.; Scardulla, F.; D’Acquisto, L. Fluid-Solid Interaction Analysis for Developing In-Situ Strainand Flow Sensors for Prosthetic Valve Monitoring. Sensors 2024, 24, 5040. [Google Scholar] [CrossRef]
  90. Chen, L.; Yu, Y.; Hou, G. Sharp-interface immersed boundary lattice Boltzmann method with reduced spurious-pressure oscillations for moving boundaries. Phys. Rev. E 2013, 87, 053306. [Google Scholar] [CrossRef]
  91. Fringand, T.; Cheylan, I.; Lenoir, M.; Mace, L.; Favier, J. A stable and explicit fluid–structure interaction solver based on lattice-Boltzmann and immersed boundary methods. Comput. Methods Appl. Mech. Eng. 2024, 421, 116777. [Google Scholar] [CrossRef]
  92. Ginzburg, I.; Verhaeghe, F.; d’Humières, D. Study of the simple hydrodynamic solutions with the two-relaxation-times Lattice Boltzmann scheme. Commun. Comput. Phys. 2008, 3, 519–581. [Google Scholar]
  93. Ginzburg, I.; Klein, P.; Reinel-Bitzer, D.; Steiner, K.; Wiegman, A. Annual Report 2000; Fraunhofer-Institut für Techno- und Wirtschaftsmathematik (ITWM): Kaiserslautern, Germany, 2000; pp. 42–51. [Google Scholar]
  94. Marson, F.; Thorimbert, Y.; Chopard, B.; Ginzburg, I.; Latt, J. Enhanced single-node Lattice Boltzmann boundary condition for fluid flows. Phys. Rev. E. 2021, 103, 053308. [Google Scholar] [CrossRef] [PubMed]
  95. Chen, L.; Yu, Y.; Lu, J.; Hou, G. A comparative study of lattice Boltzmann methods using bounce-back schemes and immersed boundary ones for flow acoustic problems. Int. J. Numer. Meth. Fluids 2014, 74, 439–467. [Google Scholar] [CrossRef]
  96. Tao, S.; Hu, J.; Guo, Z. An investigation on momentum exchange methods and refilling algorithms for lattice Boltzmann simulation of particulate flows. Comput. Fluids 2016, 133, 1–14. [Google Scholar] [CrossRef]
  97. Talon, L.; Bauer, D.; Gland, N.; Youssef, S.; Auradou, H.; Ginzburg, I. Assessment of the two relaxation time Lattice-Boltzmann scheme to simulate Stokes flow in porous media. Water Resour. Res. 2012, 48, W04526. [Google Scholar]
  98. Ginzburg, I. Multiple anisotropic collisions for advection–diffusion Lattice Boltzmann schemes. Adv. Water Res. 2013, 51, 381–404. [Google Scholar] [CrossRef]
  99. Ginzbourg, I.; Adler, P.M. Boundary flow condition analysis for the three-dimensional lattice Boltzmann model. J. Phys. II France 1994, 4, 191–214. [Google Scholar] [CrossRef]
  100. Bouzidi, M.; Firdaouss, M.; Lallemand, P. Momentum transfer of a Boltzmann-lattice fluid with boundaries. Phys. Fluids 2001, 13, 3452–3459. [Google Scholar] [CrossRef]
  101. Yu, D.; Mei, R.; Luo, L.-S.; Shyy, W. Viscous flow computations with the method of lattice Boltzmann equation. Prog. Aerosp. Sci. 2003, 39, 329. [Google Scholar] [CrossRef]
  102. Zhao, W.; Huang, J.; Yong, W.-A. Boundary conditions for kinetic theory based models I: Lattice Boltzmann models. SIAM J. Multisc. Model. Simul. 2019, 17, 854–872. [Google Scholar] [CrossRef]
  103. El Omari, K.; Younes, E.; Burghelea, T.; Castelain, C.; Moguen, Y.; Le Guer, Y. Active chaotic mixing in a channel with rotating arc-walls. Phys. Rev. Fluids 2022, 6, 024502. [Google Scholar] [CrossRef]
  104. Aharonov, E.; Rothman, D.H. Non-Newtonian flow (through porous media): A Lattice-Boltzmann method. Geophsy. Res. Lett. 1993, 20, 679–682. [Google Scholar] [CrossRef]
  105. Giraud, L.; d’Humières, D.; Lallemand, P. A lattice Boltzmann model for Jeffreys viscoelastic fluid. EPL (Europhys. Lett.) 1998, 42, 625. [Google Scholar] [CrossRef]
  106. Duque-Daza, C.; Alexiadis, A. A Simplified Framework for Modelling Viscoelastic Fluids in Discrete Multiphysics. Chem. Eng. 2021, 5, 61. [Google Scholar] [CrossRef]
  107. Ginzburg, I.; Steiner, K. A free-surface Lattice Boltzmann method for modelling the filling of expanding cavities by Bingham Fluids. Philos. Trans. R. Soc. Lond. A 2002, 360, 453–466. [Google Scholar] [CrossRef] [PubMed]
  108. Talon, L.; Bauer, D. On the determination of a generalized Darcy equation for yield-stress fluid in porous media using a Lattice-Boltzmann TRT scheme. Eur. Phys. J. E 2013, 36, 139. [Google Scholar] [CrossRef] [PubMed]
  109. Thibaud, C.; Laurent, T. Generalization of Darcy’s law for Bingham fluids in porous media: From flow-field statistics to the flow-rate regimes. Phys. Rev. E 2015, 91, 023011. [Google Scholar]
  110. Ferrari, M.A.; Franco, A.T. Exploring the periodic behavior of the lid-driven cavity flow filled with a Bingham fluid. J. Non-Newtonian Fluid Mech. 2023, 36, 105030. [Google Scholar] [CrossRef]
  111. Bresolin, C.S.; Fiorot, G. On the stability and accuracy of TRT Lattice-Boltzmann method for non-Newtonian Ostwald-de Waele fluid flows. Comput. Fluids 2024, 282, 106388. [Google Scholar] [CrossRef]
  112. Boyd, J.; Buick, J.M.; Green, S. Analysis of the Casson and Carreau-Yasuda Non-Newtonian Blood Models in Steady and Oscillatory Flows Using the Lattice Boltzmann Method. Phys. Fluids 2007, 19, 093103. [Google Scholar] [CrossRef]
  113. Xie, C.; Balhoff, M.T. Lattice Boltzmann modeling of the apparent viscosity of thinning-elastic fluids in porous media. Transp. Porous Media 2021, 137, 63–86. [Google Scholar] [CrossRef]
  114. Batôt, G.; Talon, L.; Peysson, Y.; Fleury, M.; Bauer, D. Analytical and numerical investigation of the advective and dispersive transport in Herschel–Bulkley fluids by means of a Lattice–Boltzmann Two-Relaxation-Time scheme. Chem. Eng. Sci. 2016, 141, 271–281. [Google Scholar] [CrossRef]
  115. Patarin, J.; Bleses, D.; Magnin, A.; Guerin, S.; Malbert, C.-H. Rheological characterization of gastric juices from bread with different amylose/amylopectin ratios. Food Digest. Res. Curr. Opin. 2015, 6, 2–9. [Google Scholar] [CrossRef]
  116. Murayama, T.; Yoshino, M.; Hirata, T. Three-Dimensional Lattice Boltzmann Simulation of Two-Phase Flow Containing a Deformable Body with a Viscoelastic Membrane. Commun. Comput. Phys. 2011, 9, 1397–1413. [Google Scholar] [CrossRef]
  117. Murthy, J.S.N.; Kolluru, P.K.; Kumaran, V.; Ansumali, S. Lattice Boltzmann Method for Wave Propagation in Elastic Solids. Comm. Comput. Phys. 2017, 23, 1223–1240. [Google Scholar]
  118. Maquart, T.; Noël, R.; Courbebaisse, G.; Navarro, L. Towards a Lattice Boltzmann Method for Solids—Application to Static Equilibrium of Isotropic Materials. Appl. Sci. 2022, 12, 4627. [Google Scholar] [CrossRef]
  119. Boolakee, O.; Geier, M.; De Lorenzis, L. Dirichlet and Neumann boundary conditions for a lattice Boltzmann scheme for linear elastic solids on arbitrary domains. Comput. Methods Appl. Mech. Eng. 2023, 415, 116225. [Google Scholar] [CrossRef]
  120. Müller, H.; Panov, J.; Müller, R. Improvement of the stability with a TRT-collision scheme in a Lattice Boltzmann method for elastic solids. Proc. Appl. Math. Mech (PAMM) 2024. [Google Scholar] [CrossRef]
  121. Ginzburg, I.; Steiner, K. Lattice Boltzmann model for free-surface flow and its application to filling process in casting. J. Comp. Phys. 2003, 185, 61–99. [Google Scholar] [CrossRef]
Figure 1. The large intestine. Reproduced from [54].
Figure 1. The large intestine. Reproduced from [54].
Fluids 10 00022 g001
Figure 2. “The 3D model of the biomechanical Dynamic Colon Model of human proximal colon with focus on caecum—ascending region”, reproduced from Figure S2a [9]; segments 1 and 10 are adjacent to the Caecum and the Hepatic flexures, respectively.
Figure 2. “The 3D model of the biomechanical Dynamic Colon Model of human proximal colon with focus on caecum—ascending region”, reproduced from Figure S2a [9]; segments 1 and 10 are adjacent to the Caecum and the Hepatic flexures, respectively.
Fluids 10 00022 g002
Figure 3. “PC cine-MRI image of the DCM at the midpoint (cross-section) of segment 6”, reproduced from Figure 7C [64].
Figure 3. “PC cine-MRI image of the DCM at the midpoint (cross-section) of segment 6”, reproduced from Figure 7C [64].
Fluids 10 00022 g003
Figure 4. Relative occlusion variation o ( t ) / o n 1 and the relative total volume variation δ V ( t ) = V ( t ) V ( o = o n ) 1 are displayed for 10 segments together; the motility cycle is followed by equilibration towards the steady state at the neutral occlusion degree. (left): M O T - I S L O W from Table 1 with [ o r , n , o f , n ] = [ 0.2 , 0.6 ] ; (right): M O T - I I V M A X from Table 1 with [ o r , n , o f , n ] = [ 0.6 , 0.4 ] .
Figure 4. Relative occlusion variation o ( t ) / o n 1 and the relative total volume variation δ V ( t ) = V ( t ) V ( o = o n ) 1 are displayed for 10 segments together; the motility cycle is followed by equilibration towards the steady state at the neutral occlusion degree. (left): M O T - I S L O W from Table 1 with [ o r , n , o f , n ] = [ 0.2 , 0.6 ] ; (right): M O T - I I V M A X from Table 1 with [ o r , n , o f , n ] = [ 0.6 , 0.4 ] .
Fluids 10 00022 g004
Figure 5. (top-left) O C L - 1 in “P” pipe; (top-right) o ( t ) = o r 0.293 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f = o Δ , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 1 in Table A2.
Figure 5. (top-left) O C L - 1 in “P” pipe; (top-right) o ( t ) = o r 0.293 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f = o Δ , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 1 in Table A2.
Fluids 10 00022 g005
Figure 6. (top-left) O C L - 2 in “P” pipe; (top-right) o ( t ) = o r 0.201 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f 0.704 , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 2 in Table A3.
Figure 6. (top-left) O C L - 2 in “P” pipe; (top-right) o ( t ) = o r 0.201 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f 0.704 , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 2 in Table A3.
Fluids 10 00022 g006
Figure 7. (top-left) O C L - 3 in “P” pipe; (top-right) o ( t ) = o r 0.431 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f = o max 0.862 , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 3 in Table A4.
Figure 7. (top-left) O C L - 3 in “P” pipe; (top-right) o ( t ) = o r 0.431 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f = o max 0.862 , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 3 in Table A4.
Fluids 10 00022 g007
Figure 8. (top-left) O C L - 4 in “P” pipe; (top-right) o ( t ) = o r 0.246 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f = o max 0.862 , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 4 in Table A5.
Figure 8. (top-left) O C L - 4 in “P” pipe; (top-right) o ( t ) = o r 0.246 , t = 25 s ; (bottom-left) o ( t ) [ o n , o f ] , t = 28 s ; (bottom-right) o ( t ) = o f = o max 0.862 , t = 30 s . Velocity field P - P E R H I V I S S L O W with M O T - I O C L - 4 in Table A5.
Fluids 10 00022 g008
Figure 9. “C” pipe with [ o r , n , o f , n ] = [ 1 5 , 3 5 ] [68]. (left) O C L - 1 , o f = o Δ ; (middle) O C L - 3 , o f = o max ; (right) O C L - 5 , o f = 2 5 ( 4 3 3 π ) .
Figure 9. “C” pipe with [ o r , n , o f , n ] = [ 1 5 , 3 5 ] [68]. (left) O C L - 1 , o f = o Δ ; (middle) O C L - 3 , o f = o max ; (right) O C L - 5 , o f = 2 5 ( 4 3 3 π ) .
Fluids 10 00022 g009
Figure 10. Case P - C L O : M O T - I : O C L - 1 : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Figure 10. Case P - C L O : M O T - I : O C L - 1 : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g010
Figure 11. Case C - C L O : M O T - I : O C L - 1 : L O V I S : F A S T . Legende (left-middle) to (left-top) and (right-middle).
Figure 11. Case C - C L O : M O T - I : O C L - 1 : L O V I S : F A S T . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g011
Figure 12. Case C - C L O : O C L - 1 : M O T - I : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Figure 12. Case C - C L O : O C L - 1 : M O T - I : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g012
Figure 13. Case P - P E R : M O T - I : O C L - 1 : L O V I S : F A S T . Legende (left-middle) to (left-top) and (right-middle).
Figure 13. Case P - P E R : M O T - I : O C L - 1 : L O V I S : F A S T . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g013
Figure 14. Case P - P E R : M O T - I : O C L - 2 : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Figure 14. Case P - P E R : M O T - I : O C L - 2 : H I V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g014
Figure 15. Case P - P E R : M O T - I : O C L - 2 : L O V I S : V M A X . Legende (right-middle) to (left-top) and (left-middle).
Figure 15. Case P - P E R : M O T - I : O C L - 2 : L O V I S : V M A X . Legende (right-middle) to (left-top) and (left-middle).
Fluids 10 00022 g015
Figure 16. Case P - C L O : M O T - I : O C L - 3 : H I V I S : V M A X . Legende (left-top) to (left-middle) and (right-middle).
Figure 16. Case P - C L O : M O T - I : O C L - 3 : H I V I S : V M A X . Legende (left-top) to (left-middle) and (right-middle).
Fluids 10 00022 g016
Figure 17. Case P - P E R : M O T - I : O C L - 3 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Figure 17. Case P - P E R : M O T - I : O C L - 3 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g017
Figure 18. Case C - C L O : M O T - I : O C L - 4 : H I V I S : V M A X .Legende (left-middle) to (left-top) and (right-middle).
Figure 18. Case C - C L O : M O T - I : O C L - 4 : H I V I S : V M A X .Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g018
Figure 19. Case P - C L O : M O T - I : O C L - 4 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Figure 19. Case P - C L O : M O T - I : O C L - 4 : H I V I S : V M A X . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g019
Figure 20. Case C - C L O : M O T - I : O C L - 5 : L O V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Figure 20. Case C - C L O : M O T - I : O C L - 5 : L O V I S : S L O W . Legende (left-middle) to (left-top) and (right-middle).
Fluids 10 00022 g020
Figure 21. Case C - P E R : M O T - I : O C L - 5 : L O V I S : S L O W . Legende (left-top) to (left-middle) and (right-middle).
Figure 21. Case C - P E R : M O T - I : O C L - 5 : L O V I S : S L O W . Legende (left-top) to (left-middle) and (right-middle).
Fluids 10 00022 g021
Table 1. Motility pattern is composed of 5 steps: (a) relax t r [ s ] ; (b) pause r t r , p [ s ] ; (c) forward t f [ s ] ; (d) pause f t f , p [ s ] , and (e) back t b [ s ] . The motility time of one segment t w sums all five components. The motility cycle over N s segments lasts T N s = t w + ( N s 2 ) t f , hereafter with N s = 10 .
Table 1. Motility pattern is composed of 5 steps: (a) relax t r [ s ] ; (b) pause r t r , p [ s ] ; (c) forward t f [ s ] ; (d) pause f t f , p [ s ] , and (e) back t b [ s ] . The motility time of one segment t w sums all five components. The motility cycle over N s segments lasts T N s = t w + ( N s 2 ) t f , hereafter with N s = 10 .
MotilityTravel Speed v w   ( PPW ) v w t r [ s ] t r , p [ s ] t f [ s ] t f , p [ s ] t b [ s ] t w [ s ] T 10 [ s ]
MOT-ISLOW v w = 0.4 cm / s [68]5050152565
FAST v w = 0.8 cm / s [68] 2.5 0 2.5 0 7.5 12.5 32.5
VMAX v w = 2 cm / s 10103513
MOT-IISLOW v w = 0.4 cm / s 1.875 055 23.125 3575
FAST v w = 0.8 cm / s 0.9375 0 2.5 2.5 11.5625 17.5 37.5
VMAX v w = 2 cm / s , Figure 5 [9] 0.375 011 4.625 715
Table 2. Estimates according to graphical output [68] of the MRI stream-wise velocity profiles in s . 6 : (a) [ v ¯ m i n , v ¯ m a x ] ; (b) [ v m i n , v m a x ] ; (c) t m i n [ 0 , t w ] ; (d) Re ( | v ¯ m i n | ) = | v ¯ m i n | l ν p h [ l = 2 cm ]. Motility pattern is imitated by the M O T - I in Table 1.
Table 2. Estimates according to graphical output [68] of the MRI stream-wise velocity profiles in s . 6 : (a) [ v ¯ m i n , v ¯ m a x ] ; (b) [ v m i n , v m a x ] ; (c) t m i n [ 0 , t w ] ; (d) Re ( | v ¯ m i n | ) = | v ¯ m i n | l ν p h [ l = 2 cm ]. Motility pattern is imitated by the M O T - I in Table 1.
M R I ν ph MRI [68] [ v ¯ min , v ¯ max ] [ cm / s ] [ v min , v max ] [ cm / s ] t min [ s ] Re ( | v ¯ min | )
L V O L - S L O W : low fluid volume of 60%, v w = 0.4 c m / s
M R I - 1 L O V I S Figure 7 [68] [ 0.1 , 0.1 ] [ 0.55 , 0.25 ] 13 0.327
M R I - 2 H I V I S Figure 7 [68] [ 0.5 , 0.1 ] [ 0.8 , 0.125 ] 12 0.48
H V O L - S L O W : high fluid volume of 80%, v w = 0.4 c m / s
M R I - 3 L O V I S Figure 9 [68] [ 0.55 , 0.25 ] [ 0.75 , 0.25 ] 10 1.8
M R I - 4 H I V I S Figure 9 [68] [ 0.5 , 0.25 ] [ 0.8 , 0.25 ] 15 0.48
L V O L - F A S T : low fluid volume of 60%, v w = 0.8 c m / s
M R I - 5 L O V I S Figure 10 [68] [ 1 , 0.25 ] [ 1.5 , 0.4 ] 7.96 6.53
M R I - 6 H I V I S Figure 10 [68] [ 1 , 0.2 ] [ 1.5 , 0.2 ] 6.82 1.92
H V O L - F A S T : high fluid volume of 80%, v w = 0.8 c m / s
M R I - 7 L O V I S Figure 11 [68] [ 0.75 , 0.4 ] [ 1.5 , 0.5 ] 7.96 4.9
M R I - 8 H I V I S Figure 11 [68] [ 1 , 0.4 ] ] 1.5 , 0.4 ] 6.82 1.92
Table 3. The relative values according to Table 2.
Table 3. The relative values according to Table 2.
MRI ν ph MRI [68] [ v ¯ min , v ¯ max ] v w R ( v ¯ ) = | v ¯ min | v ¯ max [ v min , v max ] v w R ( v ) = | v min | v max t min ( t f + t r ) t b Re ( | v ¯ min | )
L V O L - S L O W : low fluid volume of 60%, v w = 0.4 c m / s , t r = t f 5 s , t b = 15 s
M R I - 1 L O V I S Figure 7 [68] [ 0.25 , 0.25 ] 1 [ 1.375 , 0.625 ] 2.2 0.2 0.327
M R I - 2 H I V I S Figure 7 [68] [ 1.25 , 0.25 ] 5 [ 2 , 0.31 ] 6.45 0.1 ( 3 ) 0.48
L V O L - F A S T : low fluid volume of 60%, v w = 0.8 c m / s , t r = t f = 2.5 s , t b = 7.5 s
M R I - 5 L O V I S Figure 10 [68] [ 1.25 , 0.3125 ] 4 [ 1.875 , 0.5 ] 3.75 0.4 6.53
M R I - 6 H I V I S Figure 10 [68] [ 1.25 , 0.25 ] 5 [ 1.875 , 0.25 ] 7.5 0.24 1.92
H V O L - S L O W : high fluid volume of 80%, v w = 0.4 c m / s , t r = t f 5 s , t b = 15 s
M R I - 3 L O V I S Figure 9 [68] [ 1.375 , 0.625 ] 2.2 [ 1.875 , 0.625 ] 30 1.8
M R I - 4 H I V I S Figure 9 [68] [ 1.25 , 0.625 ] 2 [ 2 , 0.625 ] 3.25 0 . ( 3 ) 0.48
H V O L - F A S T : high fluid volume of 80%, v w = 0.8 c m / s , t r = t f = 2.5 s , t b = 7.5 s
M R I - 7 L O V I S Figure 11 [68] [ 0.9375 , 0.5 ] 1.675 [ 1.875 , 0.625 ] 3 0.4 4.9
M R I - 8 H I V I S Figure 11 [68] [ 1.25 , 0.5 ] 2.5 [ 1.875 , 0.5 ] 3.75 0.24 1.92
Table 4. The six occlusion scenarios are ordered according to o f . (a) The occlusion case; (b) the forward occlusion o f ; (c) the neutral occlusion o n ; (d) the relax occlusion o r ; (e) the relative interval with respect to o n : [ o r , n = o r o n 1 , o f , n = o f o n 1 ] ; (f) the relative neutral occlusion with respect to one of an inscribed triangle: o n , = o n o Δ 1 , with o Δ = 1 s Δ 0.586503 , s Δ = 3 3 4 π 0.413497 .
Table 4. The six occlusion scenarios are ordered according to o f . (a) The occlusion case; (b) the forward occlusion o f ; (c) the neutral occlusion o n ; (d) the relax occlusion o r ; (e) the relative interval with respect to o n : [ o r , n = o r o n 1 , o f , n = o f o n 1 ] ; (f) the relative neutral occlusion with respect to one of an inscribed triangle: o n , = o n o Δ 1 , with o Δ = 1 s Δ 0.586503 , s Δ = 3 3 4 π 0.413497 .
OCL-Case o f o n o r [ o r , n , o f , n ] o n , Construct
O C L - 1 o Δ = 1 3 3 4 π 5 8 ( 1 3 3 4 π ) 1 2 ( 1 3 3 4 π ) [ 1 5 , 3 5 ] 3 8 o f = o Δ is fixed, [ o r , n , o f , n ] from [68]
0.587 0.367 0.293 [ 0.2 , 0.6 ] 0.375
O C L - 2 42 19 π 30 19 π 12 19 π [ 3 5 , 2 5 ] −1+ 120 76 π 57 3 the “P” tube follows Figure S1 from [9]
0.704 0.503 0.201 [ 0.6 , 0.4 ] 0.143
O C L - 3 o max = 1 3 4 π 5 8 ( 1 5 3 4 π ) 1 2 ( 1 3 4 π ) [ 1 5 , 3 5 ] 19 3 + 12 π 8 ( 3 3 4 π ) o f = o max is fixed, [ o r , n , o f , n ] from [68]
0.862 0.539 0.431 [ 0.2 , 0.6 ] 0.08
O C L - 4 o max = 1 3 4 π 5 7 ( 1 3 4 π ) 2 7 ( 1 3 4 π ) [ 3 5 , 2 5 ] 8 ( 2 3 + π ) 7 ( 3 3 4 π ) o f = o max is fixed, [ o r , n , o f , n ] from the O C L - 2
0.862 0.616 0.246 [ 0.6 , 0.4 ] 0.05
O C L - 5 2 5 ( 4 3 3 π ) o Δ = 1 3 3 4 π 1 5 ( 4 3 3 π ) [ 1 5 , 3 5 ] 0 o f is larger than in Figure 5 from [9], [ o r , n , o f , n ] from [68]
0.938 0.587 0.469 [ 0.2 , 0.6 ] 0
Table 5. The extreme values over the six runs over { ν p h , v w } in Table A2, Table A3, Table A4, Table A5, Table A6, Table A7, Table A8, Table A9, Table A10 and Table A11 with ν p h = { H I V I S , L O V I S } , v w = { S L O W , F A S T , V M A X } ; (L-V) ★★ and (L-V) ★★★ are computed with c s 2 = 1 3 × 4 2 and c s 2 = 1 3 × 4 4 , respectively; otherwise, c s 2 = 1 3 . Note: min ( ) = min ν p h , v w ( ) , max ( ) = max ν p h , v w ( ) , the particular combination ( ν p h , v w ) where this extreme value is reached is abbreviated in parentheses; the maximum values and largest intervals by geometry and motility pattern are in bold in the columns 4–7.
Table 5. The extreme values over the six runs over { ν p h , v w } in Table A2, Table A3, Table A4, Table A5, Table A6, Table A7, Table A8, Table A9, Table A10 and Table A11 with ν p h = { H I V I S , L O V I S } , v w = { S L O W , F A S T , V M A X } ; (L-V) ★★ and (L-V) ★★★ are computed with c s 2 = 1 3 × 4 2 and c s 2 = 1 3 × 4 4 , respectively; otherwise, c s 2 = 1 3 . Note: min ( ) = min ν p h , v w ( ) , max ( ) = max ν p h , v w ( ) , the particular combination ( ν p h , v w ) where this extreme value is reached is abbreviated in parentheses; the maximum values and largest intervals by geometry and motility pattern are in bold in the columns 4–7.
Occl.Geom. min ( Re ( | v ¯ min | ) ) max ( R ( v ¯ ) ) max ( R ( v ) ) max ( Re ( | v ¯ min | ) ) max ( [ p min , p max ] [ mmHg ] ) max | ( δ M ) | min t δ V ( t ) Table
Motility pattern M O T - I
O C L - 1 P - C L O 0.41 ( H - S ) 1.39 ( H - S ) 1.57 ( H - S ) 6.52 ( L - V ) [ 6.64 , 5.35 ] × 10 2 ( L - V ) | 1.80 × 10 3 | ( H - S ) 6.68 × 10 2 Table A2
C C L O 0.43 ( H - S ) 1.39 ( H - S ) 1.51 ( L - F ) 7.13 ( L - V ) [ 6.96 , 5.35 ] × 10 2 ( L - V ) | 7.54 × 10 4 | ( H - S ) 6.68 × 10 2
P - P E R 0.28 ( H - S ) 1.88 ( H - S ) 1.65 ( H - S ) 4.88 ( L - V ) [ 7.00 , 5.33 ] × 10 2 ( L - V ) | 1.38 × 10 3 | ( H - S ) 6.55 × 10 2
C - P E R 0.30 ( H - S ) 2.12 ( H - S ) 1.69 ( H - S ) 5.32 ( L - V ) [ 6.82 , 6.94 ] × 10 2 ( L - V ) | 5.04 × 10 4 | ( H - S ) 6.55 × 10 2
O C L - 2 P - C L O 0.56 ( H - S ) 1.48 ( H - S ) 0.95 ( L - S ) 12.8 ( L - V ) [ 1.38 , 1.15 ] × 10 1 ( L - V ) | 4.32 × 10 3 | ( L - F ) 7.77 × 10 2 Table A3
C C L O 0.59 ( H - S ) 1.36 ( H - S ) 1.05 ( L - V ) 13.6 ( L - V ) [ 1.62 , 1.31 ] × 10 1 ( L - V ) | 1.53 × 10 3 | ( L - F ) 7.77 × 10 2
P - P E R 0.40 ( H - S ) 1.16 ( H - V ) / ( L - F ) 0.66 ( L - V ) 7.66 ( L - V ) [ 1.50 , 1.18 ] × 10 1 ( L - V ) | 5.30 × 10 3 | ( L - V ) 7.63 × 10 2
C - P E R 0.46 ( H - S ) 1.31 ( H V ) / ( L S ) 0.65 ( H - F ) 8.26 ( L - V ) [ 1.68 , 1.27 ] × 10 1 ( L - V ) | 1.99 × 10 3 | ( L - V ) 7.63 × 10 2
O C L - 3 P - C L O 1.71 ( H - S ) 2.66 ( H - S ) 2.05 ( H S ) / ( H F ) 32.0 ( L - V ) [ 0.42 , 1.12 ] ( H - V ) | 5.45 × 10 2 | ( H - S ) 1.35 × 10 1 Table A4
C C L O 1.83 ( H - S ) 3.11 ( H - S ) 3.36 ( H - V ) 31.2 ( L - V ) [ 1.72 , 1.40 ] × 10 1 ( L - V ) | 1.60 × 10 2 | ( H - S ) 1.35 × 10 1
P - P E R 0.80 ( H - S ) 1.79 ( L - S ) 1.33 ( L - S ) 22.0 ( L - V ) [ 0.37 , 1.12 ] ( H - V ) | 4.36 × 10 2 | ( L - V ) 1.32 × 10 1
C - P E R 0.88 ( H - S ) 2.13 ( H - F ) 1.21 ( H - S ) 16.9 ( L - V ) [ 1.74 , 0.88 ] × 10 1 ( L - V ) | 7.18 × 10 3 | ( H - S ) 1.32 × 10 1
O C L - 4 P - C L O 1.26 ( H - S ) 1.83 ( H - V ) 1.41 ( L - V ) 33.5 ( L - V ) [ 0.70 , 1.49 ] ( H - V ) | 4.91 × 10 2 | ( L - V ) 1.23 × 10 1 Table A5
C C L O 1.41 ( H - S ) 2.33 ( H - S ) 1.18 ( L - V ) 31.1 ( L - V ) [ 2.60 , 2.12 ] × 10 1 ( L - V ) | 2.21 × 10 2 | ( H - S ) 1.23 × 10 1
P - P E R 0.76 ( H - V ) 1.80 ( H - V ) 1.28 ( L - V ) 28.8 ( L - V ) [ 0.72 , 1.49 ] ( H - V ) | 3.79 × 10 2 | ( L - V ) 1.21 × 10 1
C - P E R 0.90 ( H - V ) 1.70 ( L - S ) 0.80 ( L - V ) 21.6 ( L - V ) [ 2.59 , 1.41 ] × 10 1 ( L - V ) | 1.15 × 10 2 | ( H - S ) 1.21 × 10 1
O C L - 5 C C L O 3.51 ( H - S ) 5 . 18 ( L - S ) 4.26 ( L - V ) 69.5 ( L - V ) [ 2.52 , 2.02 ] × 10 1 ( L - V ) | 3.89 × 10 2 | ( H - S ) 1.64 × 10 1 Table A6
C - P E R 1.01 ( H - V ) 1.77 ( L - V ) 1.03 ( H - S ) 28.8 ( L - V ) [ 2.29 , 1.07 ] × 10 1 ( L - V ) | 1.21 × 10 2 | ( H - S ) 1.61 × 10 1
Motility pattern M O T - I I
O C L - 1 P - C L O 0.71 ( H - S ) 1.74 ( H - S ) 1.88 ( L - V ) 16.3 ( L - V ) [ 1.10 , 1.11 ] × 10 1 ( L - V ) | 5.96 × 10 3 | ( L - V ) 1.27 × 10 1 Table A7
C - C L O 0.68 ( H - S ) 1.59 ( H - S ) 1.75 ( L - V ) 16.2 ( L - V ) [ 1.12 , 0.88 ] × 10 1 ( L - V ) | 8.71 × 10 4 | ( H - S ) 1.27 × 10 1
P - P E R 0.40 ( H - S ) 1.41 ( H - S ) 1.15 ( H - F ) 7.47 ( L - V ) [ 1.00 , 1.14 ] × 10 1 ( L - V ) | 3.86 × 10 3 | ( L - V ) 1.25 × 10 1
C - P E R 0.40 ( H - S ) 1.39 ( H - F ) / ( L - S ) 1.29 ( L - V ) 7.33 ( L - V ) [ 9.80 , 8.99 ] × 10 2 ( L - V ) | 2.92 × 10 4 | ( H - S ) 1.25 × 10 1
O C L - 2 P - C L O 1.26 ( H - S ) 1.48 ( H - F ) 1.61 ( L - F ) 32.9 ( L - V ) [ 2.68 , 3.50 ] × 10 1 ( L - V ) | 1.50 × 10 2 | ( L - V ) 1.48 × 10 1 Table A8
C - C L O 1.32 ( H - S ) 1.55 ( L - V ) 1.54 ( L - V ) 36.3 ( L - V ) [ 2.91 , 2.70 ] × 10 1 ( L - V ) | 1.21 × 10 2 | ( L - V ) 1.48 × 10 1
P - P E R 0.76 ( H - S ) 1.27 ( H - F ) 1.1 ( L - V ) 19.7 ( L - V ) [ 2.15 , 3.29 ] × 10 1 ( L - V ) | 1.32 × 10 2 | ( L - V ) 1.45 × 10 1
C - P E R 0.89 ( H - S ) 1.41 ( H - V ) 1.01 ( L - V ) 20.1 ( L - V ) [ 2.78 , 2.71 ] × 10 1 ( L - V ) | 7.83 × 10 3 | ( H - S ) 1.45 × 10 1
O C L - 3 P - C L O 3.07 ( H - S ) 4.38 ( H - S ) 2.30 ( L - V )   54.8 ( L - V )   [ 2.04 , 3.21 ] ( H - V ) | 1.70 × 10 1 | ( H - S ) 2.57 × 10 1 Table A9
C - C L O 2.85 ( H - S ) 3.61 ( H - S ) 3.97 ( L - F ) 57.7 ( L - V ) [ 4.24 , 2.02 ] × 10 1 ( L - V ) | 3.41 × 10 2 | ( H - S ) 2.57 × 10 1
P - P E R 1.10 ( H - S ) 1.67 ( H - F ) 0.93 ( H - V ) 27.1 ( L - V )   [ 2.10 , 3.04 ] ( H - V ) | 1.60 × 10 1 | ( L - V )   2.52 × 10 1
C - P E R 1.00 ( H - S ) 1.61 ( L - V ) 0.97 ( L - V ) 28.1 ( L - V ) [ 4.13 , 1.93 ] × 10 1 ( L - V ) | 1.1 × 10 2 | ( H - S ) 2.52 × 10 1
O C L - 4 P - C L O 3.26 ( H - S ) 3.27 ( H - S ) 1.79 ( H - S ) 41.9 ( L - V )   [ 2.24 , 3.94 ] ( H - V ) | 4.40 × 10 1 | ( L - V )   2.35 × 10 1 Table A10
C - C L O 3.14 ( H - S ) 2.00 ( H - S ) 1.92 ( L - F ) 81.6 ( L - V ) [ 6.19 , 3.47 ] × 10 1 ( L - V ) | 6.13 × 10 2 | ( L - V ) 2.35 × 10 1
P - P E R 1.42 ( H - S ) 2.00 ( L - F ) 1.19 ( L - S ) 34.7 ( L - V )   [ 2.21 , 3.81 ] ( H - V ) | 2.60 × 10 1 | ( L - V )   2.31 × 10 1
C - P E R 1.49 ( H - S ) 1.55 ( L - V ) 1.04 ( L - V ) 50.0 ( L - V ) [ 5.85 , 3.17 ] × 10 1 ( L - V ) | 2.64 × 10 2 | ( H - S ) 2.31 × 10 1
O C L - 5 C - C L O 6.35 ( H - S ) 7.00 ( H - S ) 4.02 ( H - S ) 88.1 ( L - V )   [ 8.91 , 3.08 ] × 10 1 ( H - V ) | 7.07 × 10 1 | ( L - V )   3.12 × 10 1 Table A11
C - P E R 1.26 ( H - S ) 0.95 ( L - V ) 0.73 ( L - V ) 36.3 ( L - V ) [ 8.65 , 2.48 ] × 10 1 ( L - V ) | 3.48 × 10 2 | ( L - V ) 3.05 × 10 1
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Ginzburg, I. The Lattice Boltzmann Method with Deformable Boundary for Colonic Flow Due to Segmental Circular Contractions. Fluids 2025, 10, 22. https://doi.org/10.3390/fluids10020022

AMA Style

Ginzburg I. The Lattice Boltzmann Method with Deformable Boundary for Colonic Flow Due to Segmental Circular Contractions. Fluids. 2025; 10(2):22. https://doi.org/10.3390/fluids10020022

Chicago/Turabian Style

Ginzburg, Irina. 2025. "The Lattice Boltzmann Method with Deformable Boundary for Colonic Flow Due to Segmental Circular Contractions" Fluids 10, no. 2: 22. https://doi.org/10.3390/fluids10020022

APA Style

Ginzburg, I. (2025). The Lattice Boltzmann Method with Deformable Boundary for Colonic Flow Due to Segmental Circular Contractions. Fluids, 10(2), 22. https://doi.org/10.3390/fluids10020022

Article Metrics

Back to TopTop