The dynamics of highly flexible micro- and nano-filaments are important to a variety of biological, medical, and industrial problems. The filament configuration variation and cross-stream migration in a microchannel are affected by thermal fluctuations in addition to elastic and viscous forces. Here, hydrogel nano-filaments with small bending Young’s moduli are utilized to elucidate the transitional behavior of elastic Brownian filaments in an oscillatory microchannel flow. A numerical model based on chain elastic dumbbells similar to the Rouse–Zimm model accounting for elastic, viscous, and random Brownian forces is proposed and implemented. In addition, a theoretical model to describe the average orientation–deformation tensor evolution for an ensemble of filaments in an oscillatory flow is proposed. The results are compared with the evolution observed in the experiments.

Dynamics of flexible micro- and nano-filaments is of importance to a variety of biological and industrial systems.1,2 Deformation of such filaments determines rheological properties of polymer solutions and transport of subcellular structures.1 Understanding the dynamics of individual polymer objects provides a direct link between the macromolecular configurations (based on molecular models) and rheological properties of polymeric systems. A comprehensive understanding of the dynamics of individual polymer macromolecules is crucial for their further application.3,4

The microscopic structure and the macroscopic response depend on both the nature of the suspended objects and the kinematics of the surrounding viscous flow. Linking mechanical and microscopic properties of the suspended objects to the macroscopic response of the suspension or polymer solution is one of the fundamental scientific challenges of soft matter physics, which remains open for a large number of important situations. When the persistence length of a filament is small in comparison to the length of the filament, thermal fluctuations can significantly affect its deformation dynamics.1 For example, it has been predicted that long fibers may perform spectacular windings to form more or less stable knots, the phenomenon of fundamental importance for biological macromolecules.5,6 This phenomenon applies to all living cells, the interior of which is a high level of molecular crowding. Preliminary simulations have shown that the crowding has a significant impact on the process of creating loops and winding of long polymer chains. This is an important information for the creation of a DNA loop in the presence and absence of protein and for folding of polypeptide chains.7,8

New experimental techniques in fabricating fibers with desired dimensions and tracking their migration and deformation have led to better capture of the dynamics of such fibers under different flow conditions.1 Coaxial electrospinning technique is a convenient method of fabricating core–shell fibers with innovative flexible hydrogel nanofilaments as a core.9 Such objects have high flexibility making them candidates for macroscopic objects, which are affected by thermal fluctuations in addition to viscous and elastic forces.9,10

There have been considerable efforts to numerically simulate the interaction of flexible filaments with fluid flow. The following approaches have been considered: the local slender-body theory, the immersed boundary method, and the bead-rod models.1 A numerical study of the three-dimensional dynamics of a flexible non-Brownian fiber in shear flow revealed that they are attracted to stationary, periodic, or quasi-periodic solutions depending on the ratio of local bending stiffness to hydrodynamic forces.11 The fiber was modeled as a chain of spherical beads.11 Reference 12 studied the effect of confinement on the configuration evolution and drift of elastic fibers. Reference 13 elucidated experimentally and numerically the dynamics of elastic filaments in shear flow. Brownian fluctuations were considered in addition to the elastic bending and viscous stresses in the numerical model developed using nonlocal slender-body hydrodynamics.13 They observed different filament dynamics including tumbling to buckling and snaking motions.13 An alternative model that incorporates all these physical effects and allows for description of the evolution of individual filament dynamics or the ensemble-averaged orientation–deformation tensor is proposed in the present work.

In particular, in the present work, we elucidate the transient deformation of the hydrogel nano-filament subjected to the elastic, viscous, and random Brownian forces in a microfluidic channel flow. In the experiments, hydrogel nano-filaments were formed by the coaxial core–shell electrospinning process14 and observed under an oscillatory flow condition. An individual flexible gelled filament is modeled as a chain of elastic dumbbells accounting for the random thermal fluctuations and viscous and elastic forces. Furthermore, the ensemble-averaged orientation–deformation tensor describing the flexible gelled filament evolution is described theoretically, and a particular case of the oscillatory flow is considered in detail. Such theoretical predictions are compared to the available experimental data.

Hydrogel filaments were formed by the coaxial core–shell electrospinning process. The hydrogel formed the core part in the core–shell fibers. The core hydrogel solution comprised monomer, cross-linker, the initiators of polymerization reaction, and a dye. N-isopropylacrylamide (NIPAAm, 97%, Sigma-Aldrich, Poland) was the monomer used, and N,N′-methylene bisacrylamide (BIS-AAm, 99.5%, Sigma-Aldrich, Poland) was used as a cross-linker. Ammonium persulfate (APS, 98%, Sigma-Aldrich, Poland) and N,N,N′,N′- tetramethylethylenediamine (TEMED, 99%, Sigma-Aldrich, Poland) were used as the initiators of polymerization reaction, and Bovine Serum Albumin (BSA) conjugated with fluorescein (BSA-FITC, Sigma-Aldrich, Poland) was used as a dye for better visualization.

A shell layer was formed from poly(L-lactide-co-caprolactone) (PLCL, 70% L-lactyde and 30% caprolactone unit, Corbion Purac, The Netherlands). Hydrogel filaments were extracted from the core–shell fibers by dissolving the shell in small amount of pure N,N-dimethylformamide (DMF, POCh, Poland). After about 30 min, deionized water was added diluting the solution to 25% DMF (v/v), and the filament shells were suspended in this mixture. The detailed information about preparation of the core and shell solutions (in particular, the concentrations of each component), as well as the fabrication of core–shell fibers by electrospinning (e.g., the parameters of the co-electrospinning process), can be found in Ref. 10.

The prepared suspension resulted in fragments of hydrogel filaments of random length. Fabrication of such discontinuous fragments (filaments) was possible due to the use of one of two strategies. The first is the control of the ratio of flow rates of the core- and shell-forming solutions. Theoretically, the lower the flow rate of the core solution, the shorter the hydrogel filaments, which are also found at greater distances from each other (cf. Fig. 4.8a in Ref. 20). The smaller the difference between the flow rates of the two core–shell fiber components, the longer the hydrogel core and the smaller the distance between successive filaments (cf. Fig. 4.8b in Ref. 20). Another way of fabricating core–shell fibers with discontinuous core (in the form of filaments) is to use polymerization reaction initiators, such as APS (ammonium persulfate, 98%, Sigma-Aldrich, Poland) and TEMED (N,N,N′,N′-tetramethylethylenediamine, 99%, Sigma-Aldrich, Poland). Probably, their inhomogeneous distribution in the precursor solution of the hydrogel causes it, at some places, not to polymerize, and as a consequence it remains liquid that can leak out of the forming core–shell fibers. This problem was solved using a photoinitiator, i.e., a compound (Irgacure 2959, 98%, Sigma-Aldrich, Poland) initiating the polymerization reaction under the influence of UV radiation.21 This solution, however, was not used in the present work because its purpose is to examine filaments with a length of not more than 100 µm. In addition, the lack of complete homogeneity of the solution may further affect the non-uniformity of filament shape, manifested, for example, by narrowing the filaments at random places along their length. Such a non-uniform shape, however, is typical in the case of biological objects such as proteins.

An overall view of the hydrogel does not reflect its actual behavior in the form of filaments with a nanometer cross-section and a length not exceeding significantly 100 μm. The filament length is important here because it causes a specific behavior. Experimental studies have shown that filaments over 40 μm in length reveal a significant response to environmental conditions, in particular, to the flow field. Such longer filaments have a tendency to bend and exhibit other conformational changes. In contrast, filaments shorter than 40 μm were much less affected by fluid flow. They were most often rotated, and the distance between the ends of a given filament had not changed significantly. The observations of filaments much longer than 100 μm revealed only slight changes in the position of individual segments of a given filament. They were usually a compact (coiling) structure moving holistically in a certain direction.10 

The prepared suspension was injected into a microfluidic channel fabricated from polydimethylsiloxane (PDMS) (cf. Fig. 1) with the length, width, and depth of 30 mm, 200 µm, and 60 µm, respectively. The oscillatory flow in the system was generated by using a custom-made peristaltic pump connected to one end of the microchannel (cf. Fig. 1). The oscillatory flow was implemented to simulate the intracellular fluid motions in the human body, for example, pulsatile blood flow in the human arteries. It is therefore, a natural type of flow that determines transport of biological particles in the human body. The other end of the microchannel was locked by a syringe filled with the suspension of hydrogel filaments (cf. Fig. 1). The oscillatory flow field determines the behavior of elongated objects, such as proteins, and here, hydrogel filaments are used as their model. The selected velocity range made it possible to observe conformational changes experienced by hydrogel filaments (bending, re-orientation, and migration) caused by the flow.

FIG. 1.

The experimental setup with the schematic including the peristaltic pump used to induce oscillatory flow in the microchannel with a hydrogel filament in the visualization area.

FIG. 1.

The experimental setup with the schematic including the peristaltic pump used to induce oscillatory flow in the microchannel with a hydrogel filament in the visualization area.

Close modal

For the observation of the filament behavior during experiments, an inverted epifluorescence microscope (Leica AM TIRF MC) was used. Because the filaments contained fluorescent dye, their observation was possible due to the use of a mercury lamp (Leica EL6000), a 20×/0.40 NA microscopic lens, and a high-gain EM-CCD camera (C9100-2, Hamamatsu). Only filaments located in the middle part of the channels (15 mm from the channel inlet and 30 µm above the bottom wall of the microchannel) were observed and analyzed during the experiment. To record the filament displacements in each series, 500–2000 individual images were obtained. The microfluidic system was located in a closed chamber, and the temperature was stabilized at 302 °K. The maximum flow velocity amplitude and the corresponding Reynolds number (based on the width of the channel) were in the 0.06 mm/s–0.27 mm/s and 0.01–0.07 ranges, respectively. The frequency of the applied flow was in the 0.26 Hz–0.59 Hz range.

In the experiments, the analyzed filaments have lengths ranging between 10 µm and 80 μm. The numerical model is developed for an arbitrary filament length. The case discussed in relation to the numerical simulations corresponds to the dimensionless length of 1, as specified in Sec. III.

A flexible gelled filament can be modeled as a chain of elastic dumbbells similar to the Rouse–Zimm model;15,16 cf. Fig. 2.

FIG. 2.

Schematic of the flexible gelled filament similar to the Rouse–Zimm model. Ri* being the position vector of node i.

FIG. 2.

Schematic of the flexible gelled filament similar to the Rouse–Zimm model. Ri* being the position vector of node i.

Close modal

The model consists of elastic springs freely-jointed through N nodes corresponding to the concentrated action of the drag force from the surrounding liquid. Let Ri* be the position vector of node i. Then, the unit local tangent vector τi=Ri*Ri1*/Ri*Ri1*, and the elastic force spanning two neighboring nodes is

(1)

where ai is the local cross-sectional radius, E is Young’s modulus, and subscript zero here and hereinafter corresponds to the unloaded configuration. Here and hereinafter, vectors and tensors are denoted by boldfaced characters.

The mass conservation of any element of the chain means that

(2)

where λi is the stretching ratio, λi=Ri*Ri1*/Ri,0*Ri1,0*.

Then, using Eq. (2), Eq. (1) takes the following form:

(3)

The gelled filament is a material entity that cannot intersect itself. To model this property, a localized self-repulsion force acting between two approaching elements can be artificially introduced to guarantee their mutual repulsion (self-avoidance),

(4)

where scalar multiplier Q and the exponent γ being sufficiently large to provide a significant repulsion force between two approaching nodes i and j, and the scalar factor Fr is chosen as

(5)
(6)

The total length of the filament is denoted as L0.

Consider a microchannel schematically shown in Fig. 3.

FIG. 3.

Top view of the microchannel of the half-width H, with a flexible nano-filament in the flow field. Coordinates X* and Y* are directed along the width and axis of the microchannel, respectively.

FIG. 3.

Top view of the microchannel of the half-width H, with a flexible nano-filament in the flow field. Coordinates X* and Y* are directed along the width and axis of the microchannel, respectively.

Close modal

The flow velocity profile is implied to be unaffected by the presence of the filaments and is denoted as

(7)

with the centerline velocity being U0 and V(x, t) being the velocity profile satisfying the no-slip conditions at the channel walls; j is the unit vector of the longitudinal direction; x* and t* are the dimensional lateral coordinate and time, respectively, and x and t are their dimensionless counterparts.

Each node experiences the Stokesian drag force,

(8)

where ζ is proportional to the viscosity of the surrounding liquid μ and the cross-sectional radius ai, i.e., ζ = 6πμai.

Inertia of the elastic dumbbell chain is negligibly small because of its small size. Because of its smallness, the elastic dumbbell is assumed to experience the Brownian motion, which corresponds to the following random force:

(9)

where D is the diffusion coefficient, τ is the characteristic time of a Brownian jump, and u0i is a randomly oriented unit vector.

Then, the force balance for any node in a chain [essentially, the Langevin equation because of the presence of the stochastic force (9)] takes the following form:

(10)

using Eqs. (3), (4), (8), and (9).

The dimensionless factors Fd and Fu are incorporated to make Eq. (10) applicable for all the nodes, including the end nodes i = 1 and i = N,

(11)

Render Eq. (10) dimensionless using the following scales: L0 for all lengths and L0/U0 for time, and introduce the following dimensionless groups:

(12)
(13)
(14)

Then, Eq. (10) takes the following dimensionless form:

(15)

Here and hereinafter, parameters without asterisks are dimensionless.

In the case of two-dimensional random force, and thus, the unit vector, u0i, a filament that was originally located in the flow plane will stay in the flow plane. Then, Eq. (15) has only two projections onto x and y axes, respectively,

(16)
(17)

where the angle θi corresponds to the direction of the randomly oriented (in plane) unit vector u0i with 0 ≤ θi ≤ 2π.

The initial conditions for Eqs. (16) and (17) are imposed at t = 0 as

(18)

The dimensionless flow velocity at a node position in the microchannel involved in Eq. (17) is described according to the experiment as

(19)

where the exponent α = 4.5 according to the experimental data, and

(20)
(21)

and ω = 2πf is the angular frequency (with f being frequency).

The velocity profile (19) was obtained as a correlation accurately describing the results of the finite-element simulations of the oscillatory flow inside a rectangular microchannel, which were also conducted in the present work. This correlation is also supported by experimental data [cf. Fig. 1(c)10]. This velocity profile is a plug-like shape already at the beginning of the oscillation cycle. Such plug-like velocity profiles in the x–y plane are very similar to the velocity profile typically found in large arteries or arterial tree networks.

The Langevin equations (16) and (17) with the initial conditions (18) are solved splitting the other physical processes from diffusion for a short dimensionless time interval Δt and then accounting separately for the diffusion at the same interval. At the without-diffusion interval, Eqs. (16) and (17) are integrated numerically using the Kutta–Merson method with GD = 0. After that, at the diffusion step, the node coordinates Xi and Yi are rearranged according to the diffusion process as following.

The probability of a diffusion jump of length Δr* during time Δt* is given by

(22)

where 0 ≤ P ≤ 1 is a random number.

In the dimensionless form, Eq. (22) can be rearranged as

(23)

accounting for the fact that τ = Δt. Note that the jump Δr is directed in a random direction, which determines the angle θi in Eqs. (16) and (17) as a random number from the interval 0 ≤ θi ≤ 2π.

After the diffusion step, during which all the nodes are shifted in random directions, there is the need to check whether any of the random steps led to a local self-intersection of the filament. If so happened, the random step for that particular node is reversed. This is done in the numerical scheme in addition to the self-repulsion force of Eq. (4), which is required because the latter is not implemented during the diffusion step.

Figure 4 shows the filament evolution in the case where the filament was initially aligned normally to the channel axis with the center of this filament located at X = 7.5. The dimensionless half-channel width is taken as 10, i.e., the width of the channel spans from X = −10 to X = +10. This filament evolves in a flow with the profile shown in Fig. 5(a) without flow oscillations. Figure 5(b) shows that the filament center initially drifts toward the lower values of X, i.e., toward the channel centerline, and then drifts back toward the wall (toward X = 10) at later times. The filament configurations shown in Figs. 4(a)–4(d) reveal that the filament initially stretches and aligns along the flow (toward the direction of the channel axis) due to the shear imposed by the flow. Additionally, the Brownian motion-induced fluctuations make the filament wiggle randomly. However, being aligned with the flow (in the channel axis direction), shear is not dominant in changing the alignment of the filament. One can consider that the Brownian motion results in filament wiggle to a certain radius, while the shear further orients the filament. A quasi-equilibrium is reached, where the alignment effect of the shear is countered by the random Brownian motion. As a result, it coils under the action of the restoring elastic force and the random force responsible for the Brownian motion.

FIG. 4.

Evolution of a filament in flow (the initial dimensionless filament length equals to 1). The filament is initially aligned normally to the channel axis and centered at X = 7.5. Dimensionless time moments: (a) t = 50, (b) t = 150, (c) t = 300, and (d) t = 500. The following values of the dimensionless groups were used in the simulations: the time step Δt = 0.005; the dimensionless groups are: GE = 1.0, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200 nodes.

FIG. 4.

Evolution of a filament in flow (the initial dimensionless filament length equals to 1). The filament is initially aligned normally to the channel axis and centered at X = 7.5. Dimensionless time moments: (a) t = 50, (b) t = 150, (c) t = 300, and (d) t = 500. The following values of the dimensionless groups were used in the simulations: the time step Δt = 0.005; the dimensionless groups are: GE = 1.0, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200 nodes.

Close modal
FIG. 5.

(a) The flow velocity profile in the half-channel (the dimensionless half-channel width is taken as 10) in the case of Fig. 4. (b) The center of the filament drifts in time. The following values of the dimensionless groups were used in the simulations: the time step Δt = 0.005, GE = 1.0, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200 nodes.

FIG. 5.

(a) The flow velocity profile in the half-channel (the dimensionless half-channel width is taken as 10) in the case of Fig. 4. (b) The center of the filament drifts in time. The following values of the dimensionless groups were used in the simulations: the time step Δt = 0.005, GE = 1.0, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200 nodes.

Close modal

In Fig. 6, the filament was initially aligned normally to the channel axis with the center located at X = 7.5 (similarly to the initial condition in Fig. 4), whereas the flow was turned off. This results in the filament deformation only due to the restoring elastic forces and Brownian motion. The coiling due to the Brownian motion is seen in Fig. 6 comparing the configurations predicted at time moments t = 50 and t = 500.

FIG. 6.

Evolution of a filament (with the initial dimensionless length equal to 1) initially aligned normally to the channel axis and centered at X = 7.5. There is no flow, and the filament is subjected only to the elastic force and Brownian diffusion. Filament configurations at the dimensionless time moments t = 50 and t = 500 are shown. The following values of the dimensionless groups were used in the simulations: the time step Δt = 0.005, GE = 1.0, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200 nodes.

FIG. 6.

Evolution of a filament (with the initial dimensionless length equal to 1) initially aligned normally to the channel axis and centered at X = 7.5. There is no flow, and the filament is subjected only to the elastic force and Brownian diffusion. Filament configurations at the dimensionless time moments t = 50 and t = 500 are shown. The following values of the dimensionless groups were used in the simulations: the time step Δt = 0.005, GE = 1.0, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200 nodes.

Close modal

The total length of the filament at the end of each time step is also predicted in both cases considered in Figs. 4–6, and the evolution of the dimensionless filament length in time is depicted in Fig. 7(a). With only Brownian diffusion and the elastic restoring force (without flow), the filament stretches from the initial dimensionless length of 1 to ∼1.5 due to the random forces and then practically does not change, keeping the same gyration radius of the corresponding coil. In the case with flow, the filament initially significantly stretches due to the effect of the shear flow to a maximum dimensionless length of ∼8. However, this stretching is accompanied by the filament alignment with the flow, which diminishes the effect of shear. Accordingly, the filament is not experiencing significant stretching anymore and begins to coil, reaching a similar gyration radius as without flow in a sufficiently long time (t ∼ 450), as seen in Fig. 7. This corresponds to the image sequence in Fig. 4, where the flow stretches and orients the filament through panels (a)–(c) of Fig. 4, and then, coiling due to Brownian diffusion prevails [panel (d) of Fig. 4].

FIG. 7.

(a) Dimensionless filament length vs time in two cases: with and without flow corresponding to Figs. 4 and 6, respectively. (b) Dimensionless filament length vs time with the input parameter GE = 5 and (c) GE = 10, while the other parameters as in panel (a): Δt = 0.005, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200.

FIG. 7.

(a) Dimensionless filament length vs time in two cases: with and without flow corresponding to Figs. 4 and 6, respectively. (b) Dimensionless filament length vs time with the input parameter GE = 5 and (c) GE = 10, while the other parameters as in panel (a): Δt = 0.005, GD = 1.0, GLR = 0.1, α = 4.5, GR = 10−10, and N = 200.

Close modal

In the experiments, however, stretching might not be that large as in Fig. 7(a). The predictions with the increased values of GE = 5 and GE = 10, i.e., for stiffer filaments with increased Young’s moduli, are shown in Figs. 7(b) and 7(c), respectively. It is seen that the predicted conformations are qualitatively similar to that of Fig. 7(a), albeit the filament stretching is much lower: ∼1.7 for GE = 5 and ∼1.5 for GE = 10.

The evolution of an ensemble of flexible gelled filaments can be characterized by the average orientation–deformation tensor Λ. Consider an entire filament as an elastic dumbbell characterized by its end-to-end vector R. In the absence of the Brownian diffusion and elastic relaxation, the filament evolution is fully dictated by the surrounding flow and is described as

(24)

where Ft(t + Δt) is the gradient-of-deformation tensor from moment t to moment t + Δt.17 

Note that in the present section, all parameters are dimensional, and asterisks are not used for the sake of brevity.

For Δt → 0,

(25)

where I is the unit tensor and ∇v is the velocity-gradient tensor fully determined by the flow field in the channel. By definition, the second rank tensor Λ = ⟨RR⟩, where ⟨⟩ denotes the ensemble averaging.

Then, using Eq. (25), one obtains

(26)

Accordingly, the change of Λ from time t to time t + Δt is

(27)

where subscript 1 denotes the flow-related deformation only, and the fact that Brownian diffusion and elastic relaxation are not accounted for yet.

Now, consider the effect of the elastic relaxation and Brownian diffusion on the change of the orientation–deformation tensor Λ = ⟨RR⟩ in time. In this case, the effect of the flow should not be accounted for because it is already taken into account in Eq. (27). A node in the equivalent elastic dumbbell considered here experiences the Stokes drag from the surrounding liquid,

(28)

where ζ = 6πμR0, with R0 = bN1/2 being an equivalent cross-sectional radius (or the corresponding Flory result for a self-avoiding filament can be used here and hereinafter), b being the length of “the Kuhn segment” in the filament, and N being the number of such segments.

The random force acting on the node and responsible for the Brownian motion is given by

(29)

where u0 is a random unit vector.

The elastic force acting on the node is

(30)

where λ is the stretching ratio and eR is the unit vector directed along the dumbbell.

Using the mass conservation condition πa2λ=πa02 and considering small strains ε, i.e., 1/λ ≈ 1 − ε, one can rearrange Eq. (30) to the following Hookean form:

(31)

Accounting for the fact that εeR = R/R0, where R0 is the unloaded filament end-to-end vector length, Eq. (31) takes the form

(32)

The force balance equation, i.e., the Langevin equation for the dumbbell node based on Eqs. (28), (29), and (32), is

(33)

Introduce the elastic relaxation time θ and the random function φ(t) as

(34)

Then, Eq. (33) takes the following form:

(35)

Integrating Eq. (35) over a time interval Δt → 0 but still much longer than the frequency of variation of the random force, one obtains

(36)

where u is a dummy variable.

Accordingly,

(37)

where m is another dummy variable.

The ergodic hypothesis is invoked to identify the ensemble averaging with the time averaging for an individual elastic dumbbell. In addition, the random force, i.e., the diffusion Brownian process is considered to be Markovian. Accordingly,

(38)

Then, averaging Eq. (37) and using the first part of Eq. (38) yield

(39)

Note that Λ = ⟨RR⟩.

Using the second part of Eqs. (34) and (38) and accounting for the fact that, essentially, τ = Δt, one obtains

(40)

The latter means that

(41)

where subscript 2 signifies that the Brownian diffusion and elastic relaxation are accounted for, whereas the flow-related deformation is excluded.

It is easy to show that

(42)

where M = 2 for the two-dimensional diffusion (as will be in the present case considered below) or M = 3 in the case of the general three-dimensional diffusion.

Accordingly, Eqs. (41) and (42) yield

(43)

Merging the physical processes described separately16,18,19 in Eqs. (27) and (43), the following equation for the orientation–deformation tensor is obtained:

(44)

Note that in fluid at rest, i.e., at ∇v = 0, a filament in a while, being subjected only to the Brownian motion, inevitably should reach its equilibrium “coiled” configuration corresponding to dΛ/dt = 0,

(45)

(In the second equality here, the corresponding Flory result for a self-avoiding filament can be used.)

Equation (45) yields

(46)

Using the first part of Eq. (34), this expression for the diffusion coefficient takes the form

(47)

Accounting for the fact that Nb2=R02 and ζ = 6πμR0, Eq. (47) for the diffusion coefficient takes the following form:

(48)

Similarly, the first part of Eq. (34) for the elastic relaxation time takes the form

(49)

In addition, Eq. (46) shows that

(50)

and thus, Eq. (44) takes the final form

(51)

Consider a two-dimensional shear flow similar to that in the microchannel in Fig. 3 with a gradient-of-velocity tensor with the following matrix:

(52)

where β̇ is the shear rate.

Then, Eq. (51) yields the following equations for the individual components of the orientation–deformation tensor:

(53)
(54)
(55)
(56)

Note that Eqs. (53)–(56) are dimensional.

Equations (53) and (56) show that the diagonal components of the orientation–deformation tensor will be constant and keep their equilibrium values

(57)

if initially filaments were, in average, in equilibrium.

In the present oscillatory flow corresponding to Eqs. (7) and (19),

(58)

Then, Eq. (54) is solved with the following initial condition:

(59)

Using Eqs. (58) and (59), the solution of Eq. (54) takes the form

(60)

The available data from experiments and the corresponding details are listed in Table I, where f = ω/(2π), d0 = 2a0 is the filament diameter, and β̇max is the maximum shear rate in time at a distance X from the microchannel centerline. The analysis carried out under the same conditions (in terms of temperature and medium solution) as during the motion on hydrogel filaments under oscillatory flow allowed the determination of the average value of the diameter of hydrogel filaments (126 ± 36 nm).9 

TABLE I.

Experimental details.

Filament134578910111314
R0 (μm) 12.7 22.0 36.0 57.7 17.7 23.6 12.7 18 12.7 28.6 26.2 
d0 (nm) 134 128 142 151 112 144 116 119 132 109 112 
U0 (μm/s) 132 168 934 157 209 188 60 133 316 74 151 
x/H 0.74 0.82 0.52 0.17 0.31 0.13 0.34 0.20 0.80 0.23 0.18 
f (1/s) 0.12 0.63 0.52 0.12 0.25 0.25 0.27 0.27 0.64 0.26 0.27 
β̇max (1/s) −2.11 −3.81 −4.14 −0.02 −0.16 −0.01 −0.07 −0.02 −6.59 −0.02 −0.02 
Filament134578910111314
R0 (μm) 12.7 22.0 36.0 57.7 17.7 23.6 12.7 18 12.7 28.6 26.2 
d0 (nm) 134 128 142 151 112 144 116 119 132 109 112 
U0 (μm/s) 132 168 934 157 209 188 60 133 316 74 151 
x/H 0.74 0.82 0.52 0.17 0.31 0.13 0.34 0.20 0.80 0.23 0.18 
f (1/s) 0.12 0.63 0.52 0.12 0.25 0.25 0.27 0.27 0.64 0.26 0.27 
β̇max (1/s) −2.11 −3.81 −4.14 −0.02 −0.16 −0.01 −0.07 −0.02 −6.59 −0.02 −0.02 

In Table I, R0 should be the magnitude of the ensemble-averaged tail-to-head vector at equilibrium configuration of the filament, i.e., is, essentially, the gyration radius of its coil. As an estimate, the magnitude of the tail-to-head vector in the initial orientation (at t = 0) is taken as R0 for each filament. Direct measurement of the elongation modulus on a macroscopic hydrogel yielded a value of 20 kPa.10,20 On the other hand, the Young’s modulus was determined as E = 2 kPa by experimentally studying the flow effect on filaments fixed to the channel wall.10 The calculation based on the cosine-correlation method yielded a similar order of magnitude.10 Accordingly, the value of E = 2 kPa is adopted here.

The filament diameter is found using the value of E = 2 kPa and the persistence length from the following relation:10,20

(61)

where kB is Boltzmann’s constant, T is temperature in K, and LP is the persistence length. On this basis, it was calculated that the values of the hydrogel filament diameters, which were then used in the models that are the subject of this article, were in the range of 105 nm–165 nm, which coincides with the data obtained using AFM scanning.

The filament dynamics in the experiment were recorded in the oscillatory flow. The filament profile is distinctly visible at stopped-flow conditions only, whereas throughout the transient flow stage, only in a few cases, the filament was distinctly visible. The experimental trials in which the filament was located at approximately identical distance from the microchannel axis can be lumped to find the ensemble-averaged values of Λxy. Additionally, the cases that have frequency close to each other are selected for averaging. In Figs. 8(a) and 8(b), the experimental data correspond only to the discrete stopped-flow moments: i.e., to ωt = πn (n is an integer). They are compared with the theoretical predictions for the entire time of the experiment. The experimental average values of Λxy are close to zero in all cases, and a more detailed comparison with the predicted data is impossible as the experimental results are available for only the stopped-flow moments.

FIG. 8.

(a) Average Λxy found using the data for filaments 7 and 9 in Table I. Theoretical prediction given by Eq. (60) is based on the following parameter values: R0 = 17.7 μm, E = 2 kPa, d0 = 112 nm, U0 = 209 μm/s, X/H = 0.31, and f = 0.25 1/s. (b) Average Λxy found using the data for filaments 10, 13, and 14 in Table I. Theoretical prediction given by Eq. (60) is based on the following parameter values: R0 = 28.6 μm, E = 2 kPa, d0 = 109 nm, U0 = 74 μm/s, X/H = 0.23, and f = 0.26 1/s.

FIG. 8.

(a) Average Λxy found using the data for filaments 7 and 9 in Table I. Theoretical prediction given by Eq. (60) is based on the following parameter values: R0 = 17.7 μm, E = 2 kPa, d0 = 112 nm, U0 = 209 μm/s, X/H = 0.31, and f = 0.25 1/s. (b) Average Λxy found using the data for filaments 10, 13, and 14 in Table I. Theoretical prediction given by Eq. (60) is based on the following parameter values: R0 = 28.6 μm, E = 2 kPa, d0 = 109 nm, U0 = 74 μm/s, X/H = 0.23, and f = 0.26 1/s.

Close modal

Figure 9 details the comparison of the theoretically predicted ⟨Λxy⟩ with the stopped-flow data. In Eq. (60), at the stopped-flow moments (ωt = πn), the first term on the right vanishes. For the time of the order of 1 s, the third exponential term on the right in Eq. (60) practically vanishes too. For the specific case (R0 = 17.7 μm, E = 2 kPa, d0 = 224 nm, U0 = 209 μm/s, X/H = 0.31, and f = 0.25 1/s), θ = 0.3 s [cf. Eq. (49)] and ω = 1.6 s−1 at the stopped-flow moments (ωt = πn) Λxy = ±4.08 × 10−13 m2 compared to the peak value of Λxy = ±1.8 × 10−12 m2, as shown in Fig. 9. It should be emphasized that a non-zero Λxy value at the stopped-flow moments revealed by the theory stems from the Brownian diffusion contribution.

FIG. 9.

Theoretical prediction of ⟨Λxy⟩ given by Eq. (60) is based on the following parameter values: R0 = 17.7 μm, E = 2 kPa, d0 = 224 nm, U0 = 209 μm/s, X/H = 0.31, and f = 0.25 1/s.

FIG. 9.

Theoretical prediction of ⟨Λxy⟩ given by Eq. (60) is based on the following parameter values: R0 = 17.7 μm, E = 2 kPa, d0 = 224 nm, U0 = 209 μm/s, X/H = 0.31, and f = 0.25 1/s.

Close modal

Filaments 5, 10, and 13 are visible in the transient flow, and the results are depicted in Fig. 10(a). Moreover, in Figs. 10(b) and 10(c), the transient experimental values of Λxy averaged for filaments 10 and 13 (subjected to approximately the same frequency) are compared with the theoretical prediction for ⟨Λxy⟩. For filament 5, the frequency of flow is significantly different; hence, it is not used for averaging (cf. Table I).

FIG. 10.

(a) Transient values of Λxy for filaments 5, 10, and 13 in Table I vs time. (b) The values of ⟨Λxy⟩ vs time found using the data for filaments 10 and 13 in Table I. Theoretical prediction given by Eq. (60) is based on the following parameter values: R0 = 28.6 μm, E = 2 kPa, d0 = 109 nm, U0 = 74 μm/s, X/H = 0.23, and f = 0.26 1/s. (c) The theoretical and experimental data in this panel are shown in different scales: for the theory—on the right side and for the experiment—on the left.

FIG. 10.

(a) Transient values of Λxy for filaments 5, 10, and 13 in Table I vs time. (b) The values of ⟨Λxy⟩ vs time found using the data for filaments 10 and 13 in Table I. Theoretical prediction given by Eq. (60) is based on the following parameter values: R0 = 28.6 μm, E = 2 kPa, d0 = 109 nm, U0 = 74 μm/s, X/H = 0.23, and f = 0.26 1/s. (c) The theoretical and experimental data in this panel are shown in different scales: for the theory—on the right side and for the experiment—on the left.

Close modal

Furthermore, the effect of such physical parameters as Young’s modulus E, the filament diameter d0, and equilibrium magnitude of tail-to-head vector is investigated. Figure 11(a) shows the dependence of ⟨Λxy⟩ on Young’s modulus, and Fig. 11(b) shows the dependence of ⟨Λxy⟩ on the filament diameter. The amplitude of ⟨Λxy⟩ decreases with an increase in E or d0. It is seen that the magnitude of ⟨Λxy⟩ is quite sensitive to the diameter of the filament d0. This is important as the experimentally derived value of d0 is dependent on such physical parameters as Young’s modulus E and persistence length Lp; cf. Equation (61). Hence, the inaccuracies in measurements of E and Lp, which are far from being straightforward,10,20 affect the theoretically predicted values. The effect of R0 can be significant: indirectly through the value of θ and directly via Eq. (60). For example, doubling the value of R0 results in an order of magnitude increase in the value of ⟨Λxy⟩ in Fig. 11(c). Note that a potential increase in the value of R0 is limited by the maximum length of the filament, which is around 40 μm in most cases. Additionally, the filament location (X/H) and the centerline velocity U0 have direct influence on the magnitude of ⟨Λxy⟩, as shown in Eq. (60). Specifically, if the value of X/H increases by a factor γ, i.e., becomes γ(X/H), then the value of ⟨Λxy⟩ changes by the factor γα−1. Similarly if the value of U0 changes by a factor δ, i.e., becomes δU0, the value of ⟨Λxy⟩ changes proportionally to δ⟨Λxy⟩. Experimentally, it is difficult to control and measure all the above-mentioned parameters simultaneously, which can cause the disagreement of the predictions with the data, as in Fig. 10(c). Additional future experiments would be required to obtain more accurate values of these parameters. In addition, measurement of Young’s modulus and diameter of the hydrogel filament could require an independent experimental method because relying on a single method could lead to significant inaccuracies.10,20

FIG. 11.

(a) Effect of the filament Young’s modulus on the theoretical prediction of Eq. (60). The following parameter values are used: R0 = 20.0 μm, d0 = 100 nm, U0 = 100 μm/s, X/H = 0.25, and f = 0.25 1/s. The inset shows the zoomed-in view of the result for the larger values of E. (b) Effect of the filament diameter on the theoretical prediction of Eq. (60). The following parameter values are used: R0 = 20.0 μm, E = 2 kPa, U0 = 100 μm/s, X/H = 0.25, and f = 0.25 1/s. The inset shows the zoomed-in view of the result for the larger diameters. (c) Effect of the value of R0 on the theoretical prediction of Eq. (60). The following parameter values are used: d0 = 100 nm, E = 2 kPa, U0 = 100 μm/s, X/H = 0.25, and f = 0.25 1/s. The inset shows the zoomed-in view of the result for smaller R0.

FIG. 11.

(a) Effect of the filament Young’s modulus on the theoretical prediction of Eq. (60). The following parameter values are used: R0 = 20.0 μm, d0 = 100 nm, U0 = 100 μm/s, X/H = 0.25, and f = 0.25 1/s. The inset shows the zoomed-in view of the result for the larger values of E. (b) Effect of the filament diameter on the theoretical prediction of Eq. (60). The following parameter values are used: R0 = 20.0 μm, E = 2 kPa, U0 = 100 μm/s, X/H = 0.25, and f = 0.25 1/s. The inset shows the zoomed-in view of the result for the larger diameters. (c) Effect of the value of R0 on the theoretical prediction of Eq. (60). The following parameter values are used: d0 = 100 nm, E = 2 kPa, U0 = 100 μm/s, X/H = 0.25, and f = 0.25 1/s. The inset shows the zoomed-in view of the result for smaller R0.

Close modal

In this work, the dynamics of a Brownian co-electrospun nano-filament is studied experimentally, numerically, and theoretically. The individual filament dynamics is simulated numerically using the novel model proposed in this work. The filament orientation, stretching due to flow, and coiling due to Brownian motion are fully revealed by the numerical predictions. A novel theoretical model describing the evolution of the ensemble-averaged orientation–deformation tensor in an oscillatory flow is developed. The results are compared with the filament evolution observed in the limited available experiments. The theoretically predicted values differ from the experimentally observed ones by two orders of magnitude, which stems from multiple uncertainties in the parameters involved. Nevertheless, the model does show predictions qualitatively similar to those of the averaged experimental data.

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

This work was supported by the National Agency for Academic Exchange (NAWA) (Grant No. PPI/APM/2018/1/00045/U/001) and the National Science Centre (Grant No. 2015/17/N/ST8/02012). F.P. acknowledges the financial support from the Polish Ministry of Science and Higher Education obtained through scholarships for outstanding young scientists. This work was partially supported by the LIDER program funded by the National Center for Research and Development (NCBiR) under NCBiR Grant No. LIDER/28/0067/L-7/15/NCBR/2016.

1.
O.
Du Roure
,
A.
Lindner
,
E. N.
Nazockdast
, and
M. J.
Shelley
, “
Dynamics of flexible fibers in viscous flows and fluids
,”
Annu. Rev. Fluid Mech.
51
,
539
572
(
2019
).
2.
M. J.
Shelley
, “
The dynamics of microtubule/motor-protein assemblies in biology and physics
,”
Annu. Rev. Fluid Mech.
48
,
487
506
(
2016
).
3.
H. M.
López
,
J.-P.
Hulin
,
H.
Auradou
, and
M. V.
D’Angelo
, “
Deformation of a flexible fiber in a viscous flow past an obstacle
,”
Phys. Fluids
27
(
1
),
013102
(
2015
).
4.
D. J.
Mai
,
C.
Brockman
, and
C. M.
Schroeder
, “
Microfluidic systems for single DNA dynamics
,”
Soft Matter
8
(
41
),
10560
10572
(
2012
).
5.
K.
Sadlej
,
E.
Wajnryb
,
M. L.
Ekiel-Jeżewska
,
D.
Lamparska
, and
T. A.
Kowalewski
, “
Dynamics of nanofibres conveyed by low Reynolds number flow in a microchannel
,”
Int. J. Heat Fluid Flow
31
(
6
),
996
1004
(
2010
).
6.
S.
Kuei
,
A. M.
Słowicka
,
M. L.
Ekiel-Jeżewska
,
E.
Wajnryb
, and
H. A.
Stone
, “
Dynamics and topology of a flexible chain: Knots in steady shear flow
,”
New J. Phys.
17
(
5
),
053009
(
2015
).
7.
J.
Shin
,
A. G.
Cherstvy
, and
R.
Metzler
, “
Polymer looping is controlled by macromolecular crowding, spatial confinement, and chain stiffness
,”
ACS Macro Lett.
4
(
2
),
202
206
(
2015
).
8.
J.
Shin
,
A. G.
Cherstvy
, and
R.
Metzler
, “
Kinetics of polymer looping with macromolecular crowding: Effects of volume fraction and crowder size
,”
Soft Matter
11
(
3
),
472
488
(
2015
).
9.
P.
Nakielski
,
S.
Pawłowska
,
F.
Pierini
,
W.
Liwińska
,
P.
Hejduk
,
K.
Zembrzycki
,
E.
Zabost
, and
T. A.
Kowalewski
, “
Hydrogel nanofilaments via core-shell electrospinning
,”
PLoS One
,
10
(
6
),
e0133458
(
2015
).
10.
S.
Pawłowska
,
P.
Nakielski
,
F.
Pierini
,
I. K.
Piechocka
,
K.
Zembrzycki
, and
T. A.
Kowalewski
, “
Lateral migration of electrospun hydrogel nanofilaments in an oscillatory flow
,”
PLoS One
12
,
e0187815
(
2017
).
11.
A. M.
Słowicka
,
H. A.
Stone
, and
M. L.
Ekiel-Jeżewska
, “
Flexible fibers in shear flow approach attracting periodic solutions
,”
Phys. Rev. E
101
(
2
),
023104
(
2020
).
12.
J.
Cappello
,
M.
Bechert
,
C.
Duprat
,
O.
Du Roure
,
F.
Gallaire
, and
A.
Lindner
, “
Transport of flexible fibers in confined microchannels
,”
Phys. Rev. Fluids
4
(
3
),
034202
(
2019
).
13.
Y.
Liu
,
B.
Chakrabarti
,
D.
Saintillan
,
A.
Lindner
, and
O.
Du Roure
, “
Morphological transitions of elastic filaments in shear flow
,”
Proc. Natl. Acad. Sci. U. S. A.
115
(
38
),
9438
9443
(
2018
).
14.
A. L.
Yarin
,
B.
Pourdeyhimi
, and
S.
Ramakrishna
,
Fundamentals and Applications of Micro- and Nanofibers
(
Cambridge University Press
,
Cambridge
,
2014
).
15.
M.
Doi
and
S. F.
Edwards
,
The Theory of Polymer Dynamics
(
Clarendon Press
,
Oxford
,
1986
).
16.
A. L.
Yarin
,
Free Liquid Jets and Films: Hydrodynamics and Rheology
(
Longman Scientific & Technical
,
Harlow
;
John Wiley & Sons
,
New York
,
1993
).
17.
G.
Astarita
and
G.
Marrucci
,
Principles of Non-Newtonian Fluid Mechanics
(
McGraw-Hill
,
New York
,
1974
).
18.
A. L.
Yarin
, “
Strong flows of polymeric liquids: 1. Rheological behavior
,”
J. Non-Newtonian Fluid Mechanics
37
(
2+3
),
113
138
(
1990
).
19.
A. L.
Yarin
, “
Strong flows of polymeric liquids. Part 2. Mechanical degradation of macromolecules
,”
J. Non-Newtonian Fluid Mech.
38
(
2+3
),
127
136
(
1991
).
20.
S.
Pawlowska
, “
Diffusion and migration of nanoparticles and nanofilaments suspended in oscillating flow
,” IPPT Report No. 2,
2018
.
21.
S.
Pawłowska
,
C.
Rinoldi
,
P.
Nakielski
,
Y.
Ziai
,
O.
Urbanek
,
X.
Li
,
T. A.
Kowalewski
,
B.
Ding
, and
F.
Pierini
, “
Ultraviolet light-assisted electrospinning of core–shell fully cross-linked P(NIPAAm-co-NIPMAAm) hydrogel-based nanofibers for thermally induced drug delivery self-regulation
,”
Adv. Mater. Interfaces
7
,
2000247-1
2000247-13
(
2020
).