Scaling laws for near-wall flows of thixo-elasto-viscoplastic fluids in a millifluidic channel Open Access

A

Flow curve hysteresis loop area (Pa/s)

A^*

Corrected hysteresis loop area (Pa/s)

AR

Aspect ratio

Bi

Bingham number

D-OCT

Doppler optical coherence tomography

EVP

Elastoviscoplastic

f

Acquisition frequency (Hz)

G

Pressure gradient (Pa/m)

G′

Storage modulus (Pa)

G″

Loss modulus (Pa)

H

Duct height (m)

h

Duct half-height (m)

LDA/V

Laser Doppler anemometry/velocimetry

LPT

Lagrangian particle tracking

ℓ

Length-scale (m)

MRV

Magnetic resonance velocimetry

n

Refractive index

OCT

Optical coherence tomography

PAA

Polyacrylamide

PEO

Polyethylene oxide

PIV

Particle image velocimetry

PTV

Particle tracking velocimetry

p

Pressure (Pa)

Q

Volumetric flow rate (m³/s)

q

Bousinnesq solution parameter

Re

Reynolds number

TEVP

Thixo-elasto-viscoplastic

t

Time (s)

U

Velocity (m/s)

UVP

Ultrasound velocity profiling

u,v,w

Velocity components (m/s)

u_s

Slip velocity (m/s)

W

Duct width (m)

WSO

Weak strain overshoot

x,y,z

Cartesian coordinates (m)

YSF

Yield stress fluid

β_q

Boussinesq solution parameter (deg)

γ

Shear strain amplitude

$γ ̇$

Shear rate (s⁻¹)

Δω

Doppler frequency shift (Hz)

δ

Loss tangent argument

ε

Strain

η

Viscosity (Pa s)

$η 0$

Newtonian viscosity (Pa s)

η^*

Complex viscosity (Pa s)

θ

Beam-flow inclination (deg)

λ₀

Central wavelength (m)

ν

Kinematic viscosity (m²/s)

ρ

Density (kg/m³)

σ

Shear stress (Pa)

σ_y

Yield stress (Pa)

τ_w

Wall shear stress (Pa)

υ

Particle velocity (m/s)

Φ

Phase (rad)

ω

Frequency (Hz)

Bulk/B

Bulk flow

hys

Hysteresis

R

Reference

rheo.

Rheometry

s

Slip

up

Ramp-up flow curve

w

Wall

Flow of complex fluids, including but not limited to colloidal suspensions, gels, polymeric solutions, and biological samples, is vital in many industrial applications, such as pharmaceutics, food processing, additive manufacturing, and bio-printing.¹ Contrary to Newtonian fluids whose material properties are well understood, often enabling a priori prediction of their flow behavior given initial and boundary conditions, complex fluids can exhibit a variety of complicated properties such as nonlinear rheology, inhomogeneity, phase separation, particle migration, and transient material functions, which make their flow evaluation challenging.² Complex fluids can exhibit a wide range of properties, for example shear-thinning for extrusion of polymers,³ yield-stress in applicability of paints against gravity,⁴ particle-migration in biomedical applications requiring the segregation of certain particles in vascular flows,⁵ and so on. Complex fluids may show rheological properties such as shear thinning/thickening, viscoelasticity, viscoplasticity, thixotropy, and rheopexy,⁶ which make their material characterization and flow evaluation difficult, especially when the fluids show a combination of these properties at varying levels.⁷ When flowing, the structure of such fluids can be greatly influenced at micro, meso, or macro levels. Although the importance and application scope of complex fluids are high and yet increasing, the coupling between rheological properties and flow behavior is far from understood.

While viscoelastic fluids show both elastic and viscous contributions under deformation,⁶ viscoplastic fluids show irreversible plastic deformation.⁸ The viscoplastic limit is given by the yield stress of the material, which is the stress beyond which the material undergoes plastic deformation. A pure viscoelastic behavior is shown by many polymeric solutions such as polyethylene oxide (PEO)⁹ and polyacrylamide (PAA),¹⁰ whereas many fluids exhibit viscoelasticity in combination with viscoplasticity, for instance aqueous solutions of Carbopol¹¹ and cellulose derivatives.¹² Fluids with a viscoplastic material response, also known as yield stress fluids (YSF), are becoming increasingly important for industrial applications such as bio-printing and food processing. Several consumer products such as paints, creams, pastes, and food products are designed to exhibit yield stress to improve their ease of handling, application, and texture. In addition to viscoplasticity, thixotropy is a complex rheological material behavior. Thixotropic fluids show a transient material response, where the microstructure and consequently rheological properties of the material under shear increasingly breaks down with time and is gradually recovered as the shear is removed.¹³ Many widely used fluids, including food products such as yogurt and ketchup, consumer products such as paints and gels, and industrial fluids such as grease and lubricating oil, exhibit thixotropy.

While viscoelasticity can be rheologically characterized by performing standard oscillatory tests on a rheometer,⁶ many different techniques have been attempted to characterize and measure the yield stress and thixotropy of fluids.^14,15 The thixotropic response of materials can be probed in a rheometer by performing hysteresis loop,¹¹ startup,¹⁶ and recovery tests,¹⁷ while the yield stress can be measured in oscillatory shear strain sweep tests,^18,19 startup tests,²⁰ creep tests,²¹ and by fitting viscoplastic constitutive relations such as the Herschel–Bulkley model to the flow curve.²² Based on the fluid under consideration, test configuration and initial conditions, different measurement protocols can result in varying material parameters for a given sample. In addition to the inherent difficulty associated with the measurement of these properties, thixotropy manifests itself together with yield stress in most cases, which makes studying these properties individually even more challenging.²³ Owing to the material complexity, ever increasing use and industrial relevance of rheological complex fluids, many attempts have been made to study the effects of viscoelasticity (interconnected elastic and viscous nature of the fluid), yield stress (thresholding shear stress for fluidization of the material, σ_y) and thixotropy on the flow behavior in different test configurations. This includes flow in pipes and ducts,^{20–22,24,25} capillary rheometry,^26,27 microfluidic channels,^28–30 agitators, and mixers.^31,32

In contrast to single-phase, Newtonian (e.g., air and water) fluid flows, there are many additional constraints in experimental measurements of the flow fields when a second phase is introduced to the system, or the internal structure of the fluid under study blocks certain means of accessibility in measurement (e.g., optical access).¹³ The second phase, either in the form of droplets or solid particles added to transparent fluid, deflects light upon crossing the interfaces and the medium itself and makes optical measurement techniques incapable of capturing flow information with accuracy.³³ Moreover, for complex fluids the presence of internal networks and structures between the microscopic elements of the fluid,³⁴ often leads to translucency and refractive index variations, if not full opaqueness that renders most modern flow measurement techniques inapplicable. Although there have been advancements in optical measurement techniques such as the refractive index matching schemes in particle image velocimetry (PIV) and particle tracking velocimetry (PTV) measurements to solve the presence of the secondary phase/macroscopic particles to some extent,³⁵ all those methods are only operational for a narrow range of fluid/flow conditions.

In this light, tomographic methods, often originating from medical applications and primarily designated for material/tissue scanning, have found their ways into experimental fluid mechanics, since they can penetrate into opaque fluids or even solid materials by other means than optical access. Two examples are ultrasound velocity profiling (UVP), which is used to obtain fluid flow velocity measurements based on reflected sound waves in the media,³⁶ and magnetic resonance velocimetry (MRV) that has been used to capture cross-sectional velocity fields in large scale ducts and pipe flows.^37–40 

A third unorthodox method is optical coherence tomography (OCT), which is also primarily established in medical research. The millimeter-scale of the penetration depth in OCT-based techniques have made it suitable for ophthalmologists,^41,42 dermatologists, as well as research on cardiovascular system and vascular flows.^43–45 More detailed information on the applications of OCT in medical research and diagnosis could be found in the book by Drexler and Fujimoto.⁴⁶ 

In recent years, OCT has also entered fluid mechanics laboratories as a tool for resolving material structure and velocity in a wide variety of applications. Haavisto et al. used OCT as a complementary tool for rheological assessment of the properties of complex fluids.² This work is then expanded in Haavisto et al. in a near-wall flow velocity measurement using Doppler-OCT.⁴⁷ Gowda et al. used D-OCT and 3D numerical simulations for flow field assessments alongside 3D-OCT to obtain the morphology of a two-phase flow^48,49 in a flow focusing geometry, whereas Amini et al. used Doppler-OCT to obtain tomographic information of solid particle migration in a duct flow of viscoelastic fluids.⁵⁰ More information regarding the optical principles and applications of optical coherence tomography can be found in the works of Fercher,^51,52 Brett and Tearney,⁵³ and Brezinski.⁵⁴ 

Here, we use Doppler optical coherence tomography (D-OCT) to study the influence of viscoelasticity, yield stress and thixotropy on the flow field of three widely used complex fluids: stirred yogurts, Carbopol and Laponite, in a millifluidic. The choice of the fluid types was motivated by the difference in their microstructural properties and applications in industrial settings: yogurt is a biological sample, while Carbopol and Laponite are synthetic polymers differing in their rheological properties and material preparation/handling techniques. Different concentrations/types of the fluids were chosen, which along with different flow conditions such as flow rate and the state of the microstructure (stirred and unstirred samples) reveal how the flow field is affected by the rheological properties of the fluids. Starting from the rheological characterization of the fluids considered, we obtain the wall-normal velocity profiles from the flow field recorded by the D-OCT. An ex situ protocol was developed to estimate the shear stress distribution in the millifluidic channel for different flow conditions based on the calculated shear rates. Comparing the Reynolds (⁠ $Re= ρ U B h / η$⁠; the ratio of flow inertia to viscous damping) and Bingham (⁠ $Bi= σ y h / U B η$⁠; the ratio of yield stress to viscous stress) numbers, along with the stress field, helps us to understand how the flow develops inside a millifluidic channel for thixo-elasto-viscoplastic fluids. A schematic of the workflow in this study is presented in Fig. 1.

As mentioned, three different types of elasto-viscoplastic (EVP) fluids showing thixotropic signatures were chosen for this study, yogurt, Carbopol, and Laponite (Table I). Yogurt is a widely used food product with a complex microstructure of protein and casein, exhibiting mild to high thixotropy depending on processing techniques, bacterial culture and additives. Three types of yogurts having difference in composition were used, (i) a standard stirred yogurt, (ii) a yogurt containing flavor additives, and (iii) a yogurt produced with a special bacteria culture. Respective brand names, fluid type/composition and designation are mentioned in Table I. The yogurt samples were procured off the shelf and stored in a refrigerator at 2 °C. Samples were taken out of the refrigerator and kept at room temperature for at least 30 min before measurements.

Carbopol is a synthetic polymer obtained from polyacrylic acids that can be hydrated and neutralized to obtain a non-thixotropic hydrogel, when prepared using a standard mixing protocol. However, the concentration and preparation protocol can be optimized to prepare Carbopol dispersions that exhibit thixotropy. Here, Carbopol 980NF powder (Lubrizol, USA) was added to 100 ml distilled water for the preparation of 0.1%, 0.2%, and 0.6% wt./vol Carbopol solutions. The solutions were stirred using a high torque overhead stirrer attached with a propeller-blade (Heidolph Hei Torque Ultimate 400, Heidolph Elektro GmbH and Co. KG, Kelheim, Germany) at 1200 RPM for 30 min to ensure that the Carbopol powder was homogeneously dispersed in water. The resultant acidic solution was neutralized by adding 18% wt./vol NaOH solution 2.3 times the Carbopol powder mass,²⁵ while continuously stirring the solution at 1200 RPM. The 0.1% solution was then stirred for another 30 min, whereas the 0.2% and 0.6% Carbopol solutions were stirred at 2000 RPM for 5 h to induce thixotropic effects. While stirring, the beaker containing the solution was placed in a water bath to avoid heating up of the solution. The solutions were kept at rest for at least 24 h before measurements. While the 0.1% Carbopol solution is non-thixotropic, 0.2% Carbopol exhibits mild thixotropy and 0.6% Carbopol is highly thixotropic as compared to the other two concentrations.¹¹ 

Laponite is a synthetic hectorite clay composed of disk-shaped nanoparticles that naturally exhibits thixotropy when prepared as a solution. The 1%, 2%, and 3% wt./vol solutions were prepared by adding Laponite powder to 100 ml distilled water and mixing at 2000 RPM for 1 h with the overhead stirrer to homogenize the mixture. The solutions were allowed to rest for 48–72 h before starting any measurements. Similar to Carbopol, the lowest concentration (1%) of Laponite is non-thixotropic and the medium (2%) and high (3%) concentrations are increasingly thixotropic.

Additionally, to investigate the effects of the structural state of the samples on the flow field during the timeframe of the experiments, a separate batch of each of the nine samples was prepared and additionally stirred at 1700 RPM for 5 min at room temperature. These batches were then immediately transferred to a syringe for flow field measurements using D-OCT.

The OCT measurements rely on the presence of a contrast in the fluid media for flow visualization. It was observed that while the OCT could easily pickup the contrast in yogurt, it was much more difficult to visualize the flow of Carbopol and Laponite due to their transparency. Therefore, to improve OCT signal, transparent Carbopol and Laponite solutions were dyed with 0.5% wt./vol Rhodamine solution (Sigma Aldrich) during sample preparation. Through initial rheometric assessments, it has been tested that the rheological properties of the fluids were not affected by the addition of Rhodamine.

For the purpose of this study, batches that were additionally stirred before the measurements are designated as stirred samples (“S”) and the rest of the samples are designated as unstirred samples (“US”). In total, we have three types of fluids (yogurt, Carbopol, and Laponite) and three different variants/concentrations and two pre-measurement states (stirred and unstirred).

Rheological characterizations were performed using a 702e Space rotational rheometer (Anton Paar, Graz, Austria) with a custom-made profiled bob and cup geometry (CC27/P6) in a single motor transducer configuration. The bob diameter was 27 mm, and the inner cup diameter was 29 mm, with vertical profiles on the walls to minimize slip. The C-ETD200/XL measuring cell was used to maintain the temperature constant at 23 °C. Approximately 10 ml of sample was poured into the cup geometry for each measurement. While the stirred samples were quickly poured into the cup geometry after stirring to avoid material relaxation, for the unstirred samples, a 5 min rest time was allowed to clear the shear history of the material caused during sample handling and pouring. Rheological characterizations were performed on unstirred (US) samples, and a fresh sample was used for every test.

Oscillatory shear strain and frequency sweep measurements were performed to assess the viscoelastic behavior of the samples. To determine the linear viscoelastic range, shear strain sweep tests were performed in the range $γ 0∈ 10 − 2 , 10 2$% at a constant angular frequency of 6 rad/s. Based on the identified linear viscoelastic range, angular frequency sweep tests were performed in the range $ω∈ 10 − 2 , 10 2$ rad/s at a constant strain amplitude of 0.1%.

To get an estimate of the yield stress, stress controlled steady shear tests were performed, where the shear stress was gradually increased, and the shear viscosity was measured (see Fig. S1 in supplementary material). The stress at which the viscosity drops as a response to yielding is chosen as the yield stress. This stress value was determined from the intersection points of two tangent lines drawn on the shear viscosity vs shear stress curve.

Thixotropy of the samples was determined by performing hysteresis loop tests in shear rate controlled steady shear test, wherein the (transient) shear stress, denoted $σ +$⁠, of the material was measured as the imposed shear rate is first ramped up and then ramped down in a pre-determined range. For yogurt and Laponite, the shear rate (⁠ $γ ̇$⁠) was first ramped up in the range $γ ̇∈ 10 0 , 10 3$ s⁻¹ and subsequently ramped down from $10 3 , 10 0$ s⁻¹, with a residence time of 120 and 30 s at each shear rate, respectively. For Carbopol samples, the shear rate was ramped up in the range $γ ̇∈ 10 − 2 , 10 0$ s⁻¹ and then ramped down from $10 0 , 10 − 2$ s⁻¹, with a residence time of 30 s at each shear rate. The area enclosed by the hysteresis loop (A_hys) is a measure of thixotropy. The residence times is the time after which the shear stress is measured for each shear rate; note that the shear rate might not have reached its asymptotic level at this point.

The shear stress distribution in the near wall region of the millifluidic channel was estimated by performing ex situ creep measurements of the unstirred and stirred samples. The shear stress evolution with time was measured at local shear rates calculated from the velocity profiles. While creep tests were directly performed on the unstirred samples, to obtain stirred samples, the fluids were stirred at 1700 RPM for 5 min before performing creep measurements.

The density of the samples was measured using an Ultrapyc 5000 gas pycnometer (Anton Paar, Graz, Austria) pressurized with nitrogen (N₂) gas at 20 °C. To minimize measurement errors, the samples were carefully transferred to the measurement chamber avoiding entrapment of bubbles. For each sample, the density was recorded five times with a variation less than 0.1%.

An A1 Ti-E confocal laser scanning microscope, CLSM (Nikon, Tokyo, Japan) was used to visualize the microstructure of the fluids. The most thixotropic samples, namely Y3, C3, and L3 were chosen as representatives of yogurt, Carbopol, and Laponite microstructure, respectively. A solution of 10 μg/ml Rhodamine 6G fluorescence dye (Merck, Darmstadt, Germany) in acetone was used to stain the samples. A few drops of the dye solution were applied to a glass slide and allowed to dry, followed by small amounts of the samples, which were then covered with a coverslip. The dye was excited at 488 nm using an Ar/Kr laser and the images were captured at a resolution of 1024 × 1024 pixels with 0.8 μm/pixel for yogurt and Laponite, and 0.4 μm/pixel for Carbopol.

The flow was investigated inside a rectangular millifluidic duct [Figs. 2(b) and 2(c)]. The cross-section of the duct was constant along its length (H × W = 1.00 × 3.25 mm²). The aspect ratio of the duct cross-section was chosen to ensure limited effects of the side walls, and a constant, homogeneous flow was obtained in the middle of the longer side of the cross section. The duct has a total length of approximately 80 mm, along which the measurements were performed at the streamwise position of 65 mm from the inlet. In correspondence to the flow rates, and thereby, the bulk velocities implemented throughout this study, this length has been checked and verified to guarantee a developed flow profile by measurements at positions starting from the inlet, with 2 mm increments in the steam-wise direction. The entrance effects are negligible at distances from the inlet at least one order of magnitude smaller than the position chosen for measurements (i.e., 65 mm). The flow is supplied by a syringe pump connected to the entrance of the duct through a 1 mm diameter tube. The outlet is ejected to a beaker at atmospheric conditions through a similar tube. All tests are performed at room temperature.

The rectangular duct geometry is realized through three layers of a sandwich structure [see Fig. 2(b)]: a transparent bottom plate made of 3 mm PMMA, a 1 mm stainless-steel plate in which the rectangular channel is machined, and again a 3 mm PMMA plate. These three layers are pressed between two 10 mm aluminum plates, resulting in a rigid, leakage-free structure. There are two orifices located on the bottom plate for the inlet and outlet tube connections. Through qualitative observations, it was verified that the 3D nature of the flow near the inlet and outlet orifices on the bottom surface are limited to a region shorter than 4 mm from the location of the holes for the range of Reynolds numbers investigated in the current research.

Measurements were performed at volumetric flow rates 0.4 < Q < 5.0 ml/min. This translates to the bulk velocity U_B = Q/A in the range 2.05–25.64 mm/s, where A is the channel cross-sectional area. The channel assembly was mounted on a stage with five degrees of freedom so that the exact scanning location and angular orientation with respect to the OCT beam can be adjusted (Fig. 2).

In optical coherence tomography (OCT), a low-coherence light source is used to scan through the medium. The beam is spilt into two paths, one going through the medium, while the other is a reference beam. The first step is to obtain an intensity signal showing the presence and extent of contrast points inside the medium along the beam.

Historically developed first,⁴⁶ Time-Domain OCT (TD-OCT) uses a moving reference mirror to obtain the depth information by varying the path length of the reference arm. However, as a more sensitive and recent approach, Spectral-Domain OCT (SD-OCT) uses a broadband frequency spectrum to resolve the depth. The latter approach is used in the current research.

The velocity field measurements presented in the current manuscript were performed with Doppler optical coherence tomography (D-OCT). Using the interferometry principles of Optical Coherence Tomography (OCT) to obtain the presence of contrast in medium under investigation, D-OCT utilizes the Doppler frequency to measure the velocity.⁵⁵ The relation between the velocity field and the Doppler shift in the scanned flow field is illustrated in Fig. 3.

Given the prominently axial direction of the flow field in the millifluidic flow system, one can calculate the actual, physical velocity along the duct axis (red vector denoted as υ in Fig. 3), using the velocity component captured by the backscattered light along the inclined OCT beam (the black vector labeled υ_OCT in Fig. 3).⁵⁶ For the measurements in this manuscript, the inclination angle was set to 7° as a compromise between precision in the streamwise positioning of the measurement plane, and that of the velocity readings; higher inclination only reflects optical noise away from the detector and increases the magnitude of the velocity component projected on the OCT beam direction (see Fig. 3).

A spectral domain optical coherence tomography (SD-OCT) Telesto II apparatus provided by Thorlabs (Thorlabs Inc., Newton, NJ, USA) has been used for flow measurements reported in this manuscript. The instrument has a central wavelength of 1310 nm and a nominal bandwidth of 270 nm. With the optics used, the depth resolution was 2.58 μm in water. The lateral resolution, although decoupled from the depth resolution, was chosen to the identical value 2.58 μm, resulting in square pixels.

All measurements were made with 500 consecutive frames with ∼1500 spanwise (the y-direction) positioned scans in each frame at acquisition frequencies ranging from 5.5 to 48 kHz, depending on the flow rates. The slight variations in the number of pixels in the spanwise dimension of the measured field of view are rooted in cropping of the frames (see Fig. 4). However, as the flow field is steady and no pulsation is introduced to the flow from the pump, the slight alteration in the time duration of the experiment does not have any influence, and the reason for capturing multiple frames is for the purpose of noise reduction by averaging. With 500 frames, the mean of different time series from the same flow was within ±2.5% deviations from the mean bulk velocity.

The velocity profiles are also averaged spatially over 15 profiles sampled below the upper wall, around the middle point of the duct in the spanwise direction, distributed with three-pixel increments (Fig. 3). This will ensure less locally disturbed profiles due to minor imperfections of the surface fabrications, as well as preventing the otiose over sampling of the same point.

Sample D-OCT measurements for the Newtonian case of non-diluted milk (as a baseline case to validate the measurement accuracy) are shown in Fig. 4. For this specific case, the volumetric flow rate is set to 0.4 ml/min, resulting in a bulk velocity of 2.05 mm/s. The OCT scan covers a larger field of view than the actual width of the duct [Fig. 4(a)]. The duct cross-section is determined, and the frames are cropped accordingly [Fig. 4(b)]. Based on Doppler shift information modulated on the phase difference [Fig. 4(c)], the Doppler frequencies, and the velocity field in the whole cross section, is obtained, using the OCT beam inclination angle [Fig. 4(d)].

Figure 5 illustrates the sensitivity of the Boussinesq analytical solution for the rectangular cross-section ducts to the number of terms q. At q = 0, the solution gives the Poiseuille flow solution for the infinite-span channel. Given the aspect ratio of the current duct AR = 3.25, and the flow regime, the deviation of the solution with the increase in the terms above q ≥ 1 is negligible, however, q = 10 has been chosen for validation, and point of reference for the Newtonian fluid flow.

Figure 6 illustrates the excellent agreement between the D-OCT measurements on a Newtonian fluid and the Boussinesq solution. Among the three measurement cases shown in Fig. 6, red symbols with Re = 0.50 correspond to the example presented in Fig. 4. In addition, two other Newtonian cases with higher Re numbers (Re = 3.24 and 6.28) are also plotted for the sake of completeness.

Newtonian fluids (power law fluids with n = 1) without slip result in the parabolic-like Boussinesq velocity profile. For the fluids studied in this work, n < 1, which results in deviations from the Boussinesq profile [to some extent illustrated in Fig. 1(b)]: the velocity close to the wall increases while the centerline velocity decreases. At this stage, it is also noteworthy that yield stress fluids tend to exhibit a solid-like plug flow [cf. the unstirred case shown with black triangles in Fig. 1(b)].

We start with a discussion about the rheological properties of the samples. In Sec. III A. the thixotropic material response of the samples are compared, followed by the viscoelasticity and the yield stress in Sec. III B. We then present the results of OCT measurements as velocity profiles of the unstirred and stirred samples at three different flow rates (Sec. III C).

Based on the full space of experimental cases, Sec. III D presents an ansatz based on obtained scaling laws to get an a priori estimate of the wall shear stress τ_w and slip velocity u_s. The velocity profiles and the material properties determined by rheological characterization are then used to build an estimate of shear stress distribution in the near-wall regions of the millifluidic channel (Sec. III E). Finally, Sec. III F shows microstructure of representative samples of the complex fluids considered. A summary of the sample types and their physical properties are shown in Table II.

The thixotropic behavior of the samples, characterized by hysteresis loop tests, is presented in Fig. 7. In a hysteresis loop test, the shear rate and time are coupled, i.e., the thixotropic response depends strongly on the time spent at each shear rate. Since a transient between structural breakdown and recovery is established, if a longer time per each shear rate is allowed, the material undergoes a more severe structural breakdown even at a constant shear rate, until it reaches a steady state. The microstructure might attain an equilibrium beyond the time frame of the experiment as a result of the competition between structural breakdown and recovery. For yogurt, the time required to achieve the steady state at each shear rate approaches the experimental timescale and is therefore not feasible for performing a hysteresis loop test. Creep tests were performed on the yogurt samples, and it was observed that for most shear rates the material did not attain absolute steady state even after 2000 s of shearing. This could be attributed to the interplay between long polysaccharide chains, proteins, and casein in the complex microstructure of yogurt, which prevents it from achieving a steady state even after long periods of shearing. Therefore, the hysteresis loop tests on yogurt [Fig. 7(a)] were made with a time of 120 s per shear rate, giving rise to transient flow curves in ramp up and ramp down. In the case of Carbopol and Laponite, the shear stress almost plateaued after 30 s of shearing at each shear rate in creep tests. Hence, a time of 30 s per shear rate was chosen to perform hysteresis loop tests on these two fluids.

As the shear rate is ramped up, the shear stress increases nonlinearly (except for L3 at low shear rates), which is typical of non-Newtonian shear-thinning fluids. However, when the shear rate is ramped down, the flow curve drops down below indicating the occurrence of irreversible structural breakdown in the experimental timescale. The level of thixotropy is here quantified by the area enclosed by the hysteresis loop, A_hys. This is the area of the closed region between the up- and down flow curves [for example, the region ABCA in Fig. 7(c)]. By this measure Y3 shows the highest level of thixotropy among the yogurt samples with a hysteresis loop area of A_hys = 23.52 kPa s⁻¹ followed by Y1 (A_hys = 9.26 kPa s⁻¹) and Y2 (A_hys = 12.06 kPa s⁻¹), which shows similar levels of thixotropy [Fig. 7(a)]. However, Y2 has a more non-monotonic flow curve during ramp up with a stress overshoot at lower shear rates due to the initial breakdown of the structure. The structural recovery during the ramp down in shear rate for Y1 and Y2 is higher than that of Y3, which records a shear stress 8 Pa lower than the initial stress level at the end of the experiment.

The Carbopol samples show a hysteresis loop only at low shear rates (⁠ $γ ̇∈ 10 − 2 , 10 0$ s⁻¹); at higher shear rates the up- and down curves collapse. The Carbopol samples have lower hysteresis loop areas as compared to yogurt and Laponite [Fig. 7(b)]. For C1, A_hys $∼$ 0 and the sample is non-thixotropic, i.e., the flow curves in the ramp up and ramp down cycles coincide. However, C2 (A_hys = 7.70 × 10⁻⁴kPa s⁻¹) shows mild thixotropy and C3 (A_hys = 7.01 × 10⁻³kPa s⁻¹) shows high thixotropy as compared to C2. As mentioned earlier, the Carbopol solutions are not naturally thixotropic, however, by inducing strong structural breakdown during sample preparation by increasing the mixing speed and time, higher concentrations of Carbopol (C2 and C3 in this case) can be made to show thixotropic behavior. We conjecture that this could be caused due to the intense breakdown of the material microstructure during mixing that allows the polymeric chains to untangle from their initial structure and rearrange into lower favorable energy states post mixing. This could introduce a higher elastic response to shear that resists structural breakdown at low shear rates during the ramp up cycle. However, once the material is sheared at high shear rates, structural breakdown causes the shear stress to be slightly lower in the ramp down flow curve.

Our observations are in agreement with Dinkgreve et al.,¹¹ who showed that at higher concentrations (0.25% and 0.6%), Carbopol exhibits a rheological hysteresis at low shear rates (⁠ $γ ̇∈ 10 − 2 , 10 0$ s⁻¹). Dinkgreve et al.¹¹ also studied the effect of the time per shear rate on the hysteresis loop area, which reduced as the time was decreased from 10 to 0.1 s. This confirms that the time per shear rate, when lower than the time required to reach steady state, plays a crucial role on the shear stress response of the material in a hysteresis loop test, and justifies the choice of taking 30 s as the time per shear rate in this study, which is much closer to the steady state time for the concentrations of Carbopol considered. The lower amount of Carbopol in the solution in the case of C1 sample does not allow the formation of a strong gel network [i.e., $G ′(C1)≪ G ′(C2)$] that could be rearranged to induce thixotropic response.

The Laponite samples are naturally thixotropic at mid to high concentrations and do not require any additional mixing to induce thixotropic effects. While L1 is non-thixotropic (A_hys $∼$ 0), L2 shows mild thixotropy (A_hys = 2.62 kPa s⁻¹) and L3 shows high thixotropy (A_hys = 13.47 kPa s⁻¹) as compared to L2 [Fig. 7(c)]. However, Laponite shows a monotonic evolution of the shear stress with increasing and then decreasing shear rate in the ramp up and ramp down curves, respectively, whereas the yogurt and Carbopol samples in many cases exhibit a complex shear stress response with initial overshoots.

The oscillatory strain and frequency sweep test results are shown in Fig. 8. For the special case of flow of viscoelastic fluids in a millifluidic channel, where the length-scales of the macroscopic material structure approach those of the channel, the elasticity of the fluids play a crucial role in their flow behavior. Elastic behavior of the fluids can create a normal force anisotropy within the channel and give rise to transverse forces and lateral migration of the structural elements,⁵⁷ which can have an impact on the stress distribution inside the channel. The viscoelastic behavior of the samples under consideration (Table II) is presented in Fig. 8. All the samples show elastic dominated material response or gel like behavior, $G ′> G ″$⁠, with $G ′$ and $G ″$ representing the elastic and viscous contributions of an elasto-viscoplastic material, respectively. A shear-thinning behavior and a steady dynamic moduli response can be seen when the angular frequency is increased in Figs. 8(b), 8(d), and 8(f), where the complex viscosity $η *$ is the viscosity obtained in response to the oscillatory shear inputs.

The yield stress values of the samples are given in Table II and the corresponding stress-controlled shear viscosity functions can be found in the supplementary material (Fig. S1). All three yogurt samples have similar yield stress magnitudes (⁠ $σ y∼$ 5–6 Pa) with sample Y1 showing a slightly higher yield stress. The Carbopol samples C2 and C3 have yield stresses almost twice as high as the yogurt samples, while the yield stress for C3 is much higher (⁠ $σ y$ = 40.4 Pa). In the case of Laponite, a monotonic increase in yield stress can be observed with concentration, by about 15 Pa. Although a direct comparison between steady shear and oscillatory tests is not feasible due to the different nature of the imposed deformations, it is interesting to note that the higher yield stresses of C2, C3, L2, and L3 are correlated with higher elastic responses $G ′$⁠, while the lower yield stresses of C1 and L1 are correlated with lower elastic responses, which are similar in magnitude to those of the yogurts. This could be because yield stress marks the onset of plastic deformation, for which the elastic forces have to be overcome.

A comparative summary of the yield stress determined from stress-controlled shear viscosity functions, thixotropic loop area, the (linear viscoelastic) loss tangent, $tan δ= G ″/G′$⁠, the storage moduls, $G ′$⁠, and the relative increase in $G″$ during weak strain overshoot (WSO) in strain sweep tests (herewith denoted by $ΔG″$⁠) can be found in the supplementary material (Fig. S2). The loss tangent can be taken as a measure of viscoelasticity and the storage modulus is a measure of the strength of the gel-like network. We note that while linear viscoelastic parameters, loss tangent and storage modulus, effectively characterize the unyielded state of the samples, the weak strain overshoot occurs at the linear-nonlinear transition in strain sweep tests and is thus likely related to the yielding behavior.

Overall, it appears that the only rheological parameter that distinguishes the three categories of samples, independently of the concentration in the case of Laponite and Carbopol samples, is $ΔG″$ during WSO. Thus, yogurts have no WSO, wheras Carbopol and Laponite samples show medium and high WSO. The weak strain overshoot has been observed in many soft matter systems, including soft glassy materials. Many different, sometimes contradicting, explanations for the strain overshoot, each tailored to specific materials are available in the literature, suggesting that the molecular origin and/or the physical mechanism behind the strain overshoot is not fully understood.⁵⁸ Overall, it is considered as a form of jamming of the structure prior to yielding. It has also been suggested that the overshoot arises from viscous dissipation caused due to slight rearrangement of unstable clusters in the microstructure under shear, which dies off at large shear strains, where the material kinetics is dominated by large-scale structural rearrangements.⁵⁹ Interestingly, despite having no WSO, yogurts have measurable yield stress. Both Carbopol and Laponite are highly dependent on the concentration. Generally, increasing concentration means an increase in yield stress, thixotropy, storage modulus and a decrease in damping (loss tangent). Importantly, the sample types and the concentrations cover most permutations in the rheological properties listed. For example, the yogurts, although relatively close in yield stress with L1, show a thixotropic response while L1 does not, and yoghurts have no WSO, while L1 shows considerable WSO; C2 and L2 have also comparable yield stresses, storage moduli, and loss factors; however, L2 is far more thixotropic and has a significantly higher WSO, etc.

Figure 9 shows the axial velocity magnitude for Y3. Two volumetric flow rates of 2.6 [Figs. 9(a) and 9(c) and 5.0 [Figs. 9(b) and 9(d)] ml/min are shown for stirred [Figs. 9(c) and 9(d)] and unstirred [Figs. 9(a) and 9(b)] states. The plotted velocity fields (velocity magnitude at a vertical cross-section of the duct) are limited to the top 60% of the duct half-height, and half of the duct in the crossflow direction. Since unstirred Y3 possesses yield stress, large portions of the cross-section form a plug flow; a uniform velocity in the central regions of the duct where the local stress remains below the yield stress. The unyielded plug region extends closer to the wall at z = 0 with high slip velocity magnitudes (25%–65% of the bulk velocity) as close to the wall as the depth resolution of the D-OCT permits (2.58 μm).⁶⁰ 

The flow field changes drastically if the samples are stirred [Figs. 9(c) and 9(d)]. This is apparent in both the height and spanwise directions. The plug decreases and the starting point of the plug region is approximately 30%–40% of the duct half-height from the wall. To compensate for the flow rate decrease near the surface, the maximum velocity at the central regions of the duct increases drastically. The effects of the side walls also extend further toward the central region. Figure S3 in the supplementary material includes the velocity fields of all types/flow rates for yogurt samples used in this study.

From the 2D flow fields, wall normal velocity profiles were averaged in the region below the upper surface of the duct. The resulting velocity profiles are plotted in Fig. 10. Each subplot shows the velocity profiles for one of the fluids studied in the current work. The flow rates are Q = 0.4, 2.6, and 5.0 ml/min giving bulk velocities of U_B = 2.05, 13.33, and 25.64 mm/s, respectively. The measurements are performed once for the fluid taken from rest (unstirred sample), and again after additional stirring (stirred sample).

As already observed in Fig. 9, an increase in thixotropic behavior promotes structural yielding and is expected to lead to (i) reduction, or full elimination, in the slip velocity, (ii) reduction in the size of the plug region, and thereby (iii) an increase in the yielding height in the velocity profiles, and (iv) tendency toward a parabolic, Newtonian-like velocity profile, in response to the shear history. These behaviors are observed in the velocity profiles plotted in Fig. 10 for the complete sample cases studied in this work. As the complex interplay of structure breakdown due to the said shear history, and material recovery is in place, but in different timescales, in thixotropic fluids, the actual measured data might deviate from the above listed attributes in detail. The description of the results focuses on two major aspects: the influence of flow rate and stirring as an estimation of thixotropy on the near-wall velocity profiles, and the correlations between the type of velocity profile and the rheological properties of the complex fluids investigated. The following is a detailed report on the specifics of each case (Fig. 10), and our arguments and reasoning for them.

For Y1 [Fig. 10(a)], it was observed that the unstirred and stirred samples show qualitatively the same velocity profiles at each flow rate, suggesting that the breakdown of the microstructure due to additional stirring has minimal to no effect on the flow profile. This can be validated from the hysteresis loop of sample Y1 [see Fig. 7(a)], where the shear stress recorded at the end of the hysteresis loop test is only 3 Pa lower than the initial stress response, showing that the shearing history of Y1 has only a small impact on its material behavior at low shear levels. For all flow rates, a clear deviation from the Newtonian flow profile (Boussinesq curve) can be observed, with a clear indication of the formation of a plug flow profile, a characteristic of yield stress fluids. At the lowest flow rate (Q = 0.4 ml/min) a more pronounced plug shape can be seen, as compared to the higher flow rates, presumably because the material structure gets enough time to recover from the structural breakdown caused due to shear in the channel, leading to higher viscosities compared to higher flow rate cases. Moreover, with a lower flow rate, the region where the shear stress is below the yield stress is larger, resulting in an increase in the extent of the plug region. At higher flow rates the material yields and presents a velocity profile comparable to the Newtonian case, although the profile could be uniform toward the core of the channel (i.e., plug flow), out of the field of view of the OCT scans performed in this work.

A similar trend can be observed for Y2 [Fig. 10(b)], where the lowest flow rate gives the most plug-like velocity profile. The effect of higher residence time in the channel at the lowest flow rate dominates the additional structural breakdown in stirred samples (red hallow pointers), resulting in a plug flow profile similar to the unstirred case (red solid pointers) for samples Y1 and Y2. In the case of Y3 [Fig. 10(c)], which shows the highest thixotropy (see Fig. 7), velocity profiles at higher flow rates are also in the plug regime. Due to the slightly higher yield stress of sample Y3, the velocity profiles for the unstirred samples at Q = 2.6 and 5.0 ml/min also show a pronounced plug flow, originating close to the wall. For stirred samples at the higher flow rates (black and blue hollow pointers), additional time dependent breakdown of the structure due to high thixotropy, results in a less pronounced plug profile (as already seen in Fig. 9). Although wall slip can be observed for all yogurt samples, it is especially significant for the unstirred samples of Y3.

Unlike the three yogurt samples, which show similar yield stress and viscoelastic response while having different thixotropic responses, the Carbopol and Laponite samples show a large variation in rheological properties for different concentrations (see Table II and Fig. S2 in the supplementary material). However, it is difficult to decouple the effects of these properties on the flow field and study their effects separately.

The Carbopol sample with lowest concentration C1 [Fig. 10(d)] is an elasto-viscoplastic (EVP) fluid showing no signs of thixotropy (see Fig. 7). Due to (i) the concentration being sub-critical, i.e., below the concentration at which thixotropy appears, (ii) an extremely low yield stress, and (iii) higher shear stresses, the plug regions are also limited, if occurring at all, to the center-most areas of the duct cross section. Except for the lowest flow rate Q = 0.4 ml/min, the plug region forms at z/h > 0.6, i.e., above the present measurements. For the unstirred sample at Q = 0.4 ml/min, the plug starts to reach z/h < 0.6. This is seen by the red curve approaching the vertical direction. Due to low shear stresses at this low flow rate, the unyielded plug at the core extends toward the wall. As an extension to this reasoning, the other profiles could be showing similar behavior, toward the core of the channel.

For C2 [Fig. 10(e)], having a higher yield stress, thixotropic and viscoelastic material response as compared to C1, velocity profiles at all flow rates show a clear deviation from the Newtonian regime. While both unstirred and stirred samples show a more pronounced plug profile at the lowest flow rate of 0.4 ml/min, at higher flow rates a seemingly linear velocity profile is observed which could be exhibiting a smaller plug region toward the channel core. A similar trend can be observed for C3 [Fig. 10(f)], where both stirred and unstirred samples show a clear plug profile at the lowest flow rate Q = 0.4 ml/min due to the high yield stress and higher residence time. At higher flow rates, all samples converge toward the Newtonian case and show a deviation toward a plug regime starting from z/h > 0.3. Moreover, like in the case of yogurt, higher slip velocities can be observed with the increase in the thixotropy.

The case of C2 concerns a critical concentration of the material, exhibiting transient regime attributes between non-thixotropic and thixotropic Carbopol samples. In this regard, the material breakdown and recovery do not occur homogenously in the flow field. In addition, it should be noted that microgels (such as Carbopol) are shown in the literature to have heterogeneous flow profiles in channels and ducts. It should also be noted that the normalization of the plotted velocity profiles is performed based on the a priori set value of flow rate by the pump, and not based on the integration of the flow rates obtained through the velocimetry results. Therefore, discrepancies seen in Fig. 10(e) are expected for such heterogeneous flow fields.

It can be clearly seen that for the thixotropic C2 and C3 samples, the velocity profiles are wavy and exhibit variations with length-scales approximately 50 μm, in particular at higher flow rates. One reason for the variations is bubbles trapped inside Carbopol samples with high concentrations. Bubbles are difficult to eliminate from strong gels with high yield stress, and one risks breaking up the microstructure during strong de-airing processes. Another source of fluctuation in the profiles could be attributed to the constant competition between structural breakdown and recovery. The structural breakdown caused by the shear at the wall and interlayer shear is recovered very fast due to the relatively shorter time required to reach the material steady state and also by the lower thixotropic response of the samples, i.e., smaller hysteresis loop areas. This is more prominent at higher flow rates where breakdown of the structure is faster. For a highly thixotropic material, like yogurt, the structural breakdown caused due to shear is maintained during the timeframe of the experiment and is gradually and slowly opposed by the structural recovery, leading to the longer times required to reach steady state. For a non-thixotropic sample like C1, breakdown in material structure due to yielding is much less and recovery more instant, as compared to more transient thixotropic effects, which sustain structure breakdown longer. However, for samples like C2 and C3, showing a less thixotropic response, the OCT scans the interplay between the structural breakdown and recovery in each frame and when averaged, gives rise to this velocity variation, which can also be seen for Laponite samples at higher flow rates, although much less than that for Carbopol samples. The underlying layering caused by the material breakdown and recovery time-scales, while interacting with the imposed varying shear rate profiles in different elevations from the wall, is speculated to have led to the above-mentioned velocity variations. It should be noted that over the time-scales of the experiments performed and reported in the current paper, O(10–100 s), the variations appear steady. Whether they stay steady with time-scales beyond the current study requires further investigations in future works.

Similar to the case of the non-thixotropic fluid C1, the fluid L1 [Fig. 10(g)] exhibits in velocity profiles that resemble the Newtonian flow profile, up to a normalized depth of z/h = 0.3, above which a clear deviation from the Newtonian regime can be seen in form of a plug flow. For L2 [Fig. 10(h)] the unstirred sample at the lowest flow rate (0.4 ml/min) shows a plug profile, while in the case of L3 [Fig. 10(i)], all profiles show a significant deviation from the Newtonian profile toward a plug profile, with both the unstirred and stirred samples having the largest plug at the lowest flow rate. This pattern seems to be similar to that observed in yogurts, where for Y3, having the highest yield stress among the yogurt samples, velocity profiles at higher flow rates are also in the plug regime. However, for L3, having a lower thixotropic area than Y3, a profound plug profile can only be observed for the cases of least structural breakdown, i.e., at the lowest flow rate. Like in the previous cases, the slip velocity increases with higher concentration and thixotropic response.

Overall, we can identify three types of velocity profiles: (i) Newtonian-like velocity profiles that follow the analytical Boussinesq solution relatively closely (C1 and L1); (ii) velocity profiles that are non-Newtonian with some evidence of high gradient regions near the wall at the lowest flow rates investigated (C2, C3, L2, L3); and (iii) velocity profiles that show pronounced unyielded regions (i.e., plug flow) toward the center of the millifluidic duct (Y1, Y2, Y3). Considering the linear viscoelastic, yield stress, thixotropic, and weak strain overshoot behavior of the samples analyzed (Fig. 11), as expected, velocity profiles of type (i) generally correspond to low dynamic moduli and $tan δ$⁠. In contrast, type (ii) velocity profiles occur for low $tan δ$ but higher overall dynamic moduli. In this category, the Carbopol and Laponite samples have similar $tan δ$⁠, however, in weakly nonlinear conditions (seen in strain sweep tests, Fig. 8) C2 and C3 show less pronounced weak strain overshoot compared to L2 and L3. It could be conjectured that the weak strain overshoot (i.e., microstructural jamming/yielding behavior) corresponds to diverging scaling laws in the Cole–Cole plot (Fig. 11), with critical dynamic moduli above which the velocity profiles transition from type (i) to type (ii). Finally, type (iii) velocity profiles distinguish themselves as having considerably higher $tan δ$(meaning the relative viscous contribution is higher) and consequently no WSO.

The specific selection and grouping of these four non-dimensional numbers have been done in a way to allow the following ansatz for an a priori calculation of wall shear stress τ_w and slip velocity u_s. Using Π₁ and Π₂ [Fig. 12(a)], a power law dependency is shown over more than three orders of magnitude in Π₂, and more than five orders of magnitude in Π₁. The non-dimensional numbers Π₁ and Π₂ are composed of fluid level parameters (hysteresis loop area A^*, as an indicator of thixotropy, and the yield stress σ_y, showing the plastic threshold of the material), as well as flow level parameters [shear rate $γ ̇ w$⁠, as the response of the near-wall yielded material to the imposed shear stress, and the local viscosity $η γ ̇ w$]. As the near-wall local viscosity is a function of the near-wall shear rate, which is fully characterized through shear rheometry, a fit to the flow curve alongside the scaling relationship between Π₁ and Π₂ could be used to solve for the shear rate and the near-wall viscosity (as a function of shear rate). It is then straightforward to calculate the wall shear stress as $τ w= γ ̇ w·η γ ̇ w.$ This shear stress τ_w is of utmost importance as it provides information needed for further calculations of pressure loss. It is worth emphasizing that by having the referred scaling power law, there is no need for performing further high-resolution near-wall velocimetry, as the only requirement to obtain the wall shear stress is a shear flow curve of the fluid in question.

Following the above analysis, the obtained wall shear stress τ_w can be used in the scaling relationship between non-dimensional groups Π₃ and Π₄ [Fig. 12(b)] leading to the calculation of the slip velocity u_s as the only unknown among the parameters in the Π₃/Π₄ scaling. This scaling covers almost six orders of magnitude of both Π₃ and Π₄. It should be noted that Π₃ is reminiscent of inertial to viscous behavior of the near-wall flow field, as it is a representation of Reynolds number calculated based on the slip velocity [see Eq. (3)].

To get an estimate of the stress distribution in the near wall region of the millifluidic channel, we performed ex situ rheological measurements on the fluids under consideration, mimicking the operating conditions encountered during flow inside the millifluidic channel. From the velocity profiles obtained, the local shear rate profile was calculated based on Eq. (1). Transient shear stress was measured in creep tests at constant shear rates obtained using Eq. (1) for all samples. From creep curves the stress at t = t_TS is recorded, where t_TS is the travel time of the fluid from the syringe to the test section in the millifluidic channel. The shear rate was assumed to be constant during the travel from the syringe outlet to the test section. While creep tests were directly performed on the unstirred samples, to obtain stirred samples, the fluids were stirred at 1700 RPM for 5 min before performing creep measurements. Overall, creep measurements were performed for both stirred and unstirred samples of each type/concentration of fluids at three different flow rates. In the following paragraphs we will discuss the (i) shear stress distribution of all the samples at a constant flow rate, Q = 0.4 ml/min, and (ii) the influence of flow rate on the shear stress distribution for the particular case of yogurt.

Shear stress in the near-wall region (z/h = 0.05 ± 0.02) scales up with yield stress for unstirred samples [Fig. 13(a)]. Unstirred samples with similar yield stress magnitude, such as Y1, Y2, and Y3, show qualitatively the same stress response in the timeframe of the measurement. Note that the shear stress is represented in log scale in Fig. 12(a). The drop in shear stress for stirred samples, as compared to unstirred case, is a result of the structural breakdown caused by stirring and shear in the millifluidic channel. For thixotropic samples, the structure additionally breaks down with time. Therefore, a linear increase in the stress difference (normalized by $σ y·$A_hys/A_up) between unstirred and stirred samples is observed with an increase in thixotropy [Fig. 13(b)]. The lack of time induced structural breakdown results in a less severe drop in shear stress for non-thixotropic samples, C1 and L1. Data obtained from creep tests show near steady state achieved by the stress response of non-thixotropic samples, C1 and L1 [Fig. 13(c)]. However, transience in material behavior increases with thixotropy and highly thixotropic samples (Y3, C3, and L1) display greater fluctuations in the stress response for unstirred cases. Importantly, especially in the yielding region it can be observed in Fig. 13(c), and it is well known in other cases, that the stress can have non-monotonic behavior including quasi-steady oscillations between the yielded and unyielded states.²³ 

For yogurt samples the local shear stress in the near-wall region (z/h = 0.05 ± 0.02) is compared at Q = 0.4, 2.6, and 5.0 ml/min for unstirred and stirred cases [Fig. 14(a)]. For a given yogurt sample with a constant yield stress and hysteresis loop area, $Re ℓ$ increases and $Bi ℓ$ decreases with increasing flow rate. However, no such trend can be observed when different types of complex fluids with varying levels of thixo-elasto-viscoplastic properties are compared. Samples Y1 and Y2 are equivalent in their yield stress and thixotropy and therefore exhibit similar difference in shear stress between the unstirred and stirred samples (normalized by $σ y·$A_hys/A_up) at a constant flow rate [Fig. 14(b)]. As the flow rate increases, an apparent yet slight fluctuation in the stress difference is observed for Y1 and Y2, caused by the transience in their stress response. Due to its higher thixotropy, Y3 records a larger stress difference than Y1 and Y2 at a given flow rate. The stress difference increases monotonically with the flow rate as a result of much higher thixotropy of Y3. The transience in material behavior due to thixotropy can be visualized by plotting the creep test results [Fig. 14(c)]. Unstirred and stirred Y1 and Y2 show decreasing shear stress as the structure has not reached a steady state in the timeframe of the experiment. The highly thixotropic unstirred Y3 shows greater transience in the stress response, suggesting that the breakdown of largely intact initial structure is slowed down by structural recovery. However, a clear decrease in the transience is observed for stirred Y3 samples. This is also aligned with the previously mentioned explanations regarding the effect of stirring and shear history on the breakdown of the material structure, which requires a longer time duration for recovery, than permitted in the experimental time-scales.

Confocal micrographs of Y3, C3, and L3 samples, showing the structure at rest, are presented in Figs. 15(a)–15(c), respectively. Note that the images are color corrected in post processing to highlight the distinction between the media (in black) and the actual structure (in green). It can be observed that the microstructure is very different for the different fluids. While for yogurt a highly networked structure can be seen, Carbopol has tightly packed globular structural cells and Laponite has rather well dispersed structural agglomerates. The amount of structure present in the material was estimated by calculating the area fraction of the structure (in green) with respect to the total area of the micrograph. Sample C3 has the highest area fraction (0.85), followed by Y3 (0.54) and L3 (0.41).

Factors such as lack of optical access and particle/dye migration within the flow make flow field measurements of complex fluids through conventional means (PIV/PTV) challenging and fraught with errors. In this study, a tomographic velocimetry apparatus using D-OCT was utilized to probe the near wall flow field of TEVP fluids in a rectangular millifluidic duct. The 2D velocity profiles of unstirred and stirred samples of Yogurt, Carbopol, and Laponite show great influence of rheological material functions on the flow field. For most cases;

• The size of the plug region increases with yield stress, especially at low flow rates where low shear levels and higher residence times promote structural recovery, limiting the yielded region close to the wall. An increase in structural breakdown leads to a reduction in slip velocity and the size of the plug region. Linear (⁠ $G ′, G″$⁠) and weakly nonlinear (weak strain overshoot) viscoelastic properties appear to also play a role in the development of the velocity profiles, although their contribution is difficult to assess when considering fully developed flow in the near-wall region.

• In cases, where the structural breakdown dominates, such as the stirred cases and at high flow rates, the tendency toward a parabolic Newtonian like velocity profile increases with thixotropy, due to shear as well as time dependent structural breakdown.

  Four non-dimensional numbers were developed empirically using rheological parameters such as yield stress, thixotropy (hysteresis loop area) and viscosity as a function of shear rate, and the flow parameters wall shear stress and slip velocity.

• Two scaling power laws were obtained using data from nine samples at three different flow rates, spanning over multiple decades, and can be used to predict the wall shear stress and slip velocity of complex fluids from steady shear rheological test data, without the need to perform high-resolution near-wall velocimetry.

  Moreover, an ex situ rheometric protocol, based on the flow conditions inside the millifluidic channel, was used to estimate the near wall shear stress in different flow cases of the complex fluids considered.

• Unstirred fluid samples record higher shear stresses than stirred samples due to a largely intact structural network. At a given flow rate, shear stress scales up with the yield stress and the stress drop in stirred samples increases with thixotropy due to additional time dependent structural breakdown. For highly thixotropic samples, the drop in stress due to stirring increases with flow rate.

• The near-wall region Reynolds and Bingham numbers, modified with the Herschel Bulkley formulation, increase, and decrease with increasing flow rate, respectively. A length-scale for the near-wall regions has been proposed to incorporate the fluid properties, as well as the dynamics of the imposed flow condition, and is used in the definition of the mentioned non-dimensional numbers.

See the supplementary material for Fig. S1: yield stress determination via steady shear stress tests for all nine fluids; Fig. S2: summary of key rheological parameters for all nine fluids; and Fig. S3: two-dimensional cross sections of the streamwise velocity for the yogurts (stirred and unstirred) at three flow rates.

This project has received funding from the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie Grant Agreement No. 955605 YIELDGAP and the European Research Council through Starting Grant MUCUS (Grant No. ERC-StG-2019-852529). Authors would also like to acknowledge Thorlabs Inc. for granting the permission to use technical drawings of the Telesto II SD-OCT apparatus (Part No. Series: OCTG13 and TEL211C1).

The authors have no conflicts to disclose.

Kasra Amini and Ases Akas Mishra contributed equally to this paper.

Kasra Amini: Conceptualization (equal); Formal analysis (lead); Investigation (lead); Writing – original draft (lead); Writing – review & editing (equal). Ases Akas Mishra: Conceptualization (equal); Formal analysis (equal); Investigation (lead); Writing – original draft (lead); Writing – review & editing (equal). Amit Kumar Sivakumar: Investigation (equal). Dragana Arlov: Conceptualization (equal); Funding acquisition (equal); Writing – review & editing (equal). Fredrik Innings: Conceptualization (equal); Funding acquisition (equal). Roland Kádár: Conceptualization (equal); Funding acquisition (equal); Project administration (equal); Supervision (equal); Writing – review & editing (equal). Outi Tammisola: Conceptualization (equal); Funding acquisition (lead); Project administration (lead); Supervision (equal); Writing – review & editing (equal). Fredrik Lundell: Conceptualization (equal); Funding acquisition (equal); Methodology (equal); Supervision (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.