1. Introduction
Droplet-based microfluidics has gained widespread attention among researchers for the past two decades, owing to its capability of precise control of minute volumes (Manga Reference Manga1996; Adamson et al. Reference Adamson, Mustafi, Zhang, Zheng and Ismagilov2006; Baroud, Gallaire & Dangla Reference Baroud, Gallaire and Dangla2010; Vladisavljević et al. Reference Vladisavljević, Khalid, Neves, Kuroiwa, Nakajima, Uemura, Ichikawa and Kobayashi2013; Nozaki et al. Reference Nozaki, Yoon, Furuya, Fujita, Sekiguchi and Shoji2021). Important to mention here is that the paradigm of droplet-based microfluidics has shown significant technical and research potential specifically due to its suitability in point-of-care diagnostics (Shamloo & Hassani-Gangaraj Reference Shamloo and Hassani-Gangaraj2020), drug delivery (Yue et al. Reference Yue, You, Yang, Niu and Zhang2020), analysing and screening of bio/chemical reaction products (Zheng & Ismagilov Reference Zheng and Ismagilov2005) and many more related areas (Moon et al. Reference Moon2010; Santos et al. Reference Santos, Trujillo-Cayado, Calero, Alfaro and Muñoz2016; Madadelahi, Ghazimirsaeed & Shamloo Reference Madadelahi, Ghazimirsaeed and Shamloo2019). Researchers have explored several aspects of droplet behaviour such as droplet generation, breakup (Link et al. Reference Link, Anna, Weitz and Stone2004; Jullien et al. Reference Jullien, Tsang Mui Ching, Cohen, Menetrier and Tabeling2009; Leshansky & Pismen Reference Leshansky and Pismen2009; Leshansky et al. Reference Leshansky, Afkhami, Jullien and Tabeling2012; Hoang et al. Reference Hoang, Portela, Kleijn, Kreutzer and Van Steijn2013) and merging/coalescence (Christopher et al. Reference Christopher, Bergstein, End, Poon, Nguyen and Anna2009), in various microfluidic devices such as a T-junction, flow-focusing junction and co-flowing junction, to name a few. Note that the breakup of a mother droplet into two sister droplets in a simple passive microfluidic device such as a T-junction, has significant engineering implications. Researchers have shown that the size of the sister droplets being generated in a T-junction can be passively controlled by adjusting the length of the downstream channel (Link et al. Reference Link, Anna, Weitz and Stone2004). However, for situations in which there exist geometric limitations, the suitability of a passive breakup methodology of droplets is highly restrained. This limitation, in particular, has led to the paradigm of the active droplet breakup mechanism, whereby external field modulated forcing such as electric (Xi et al. Reference Xi, Guo, Leniart, Chong and Tan2016), magnetic (Tan & Nguyen Reference Tan and Nguyen2011), acoustic (Schmid & Franke Reference Schmid and Franke2013), temperature (Yesiloz, Boybay & Ren Reference Yesiloz, Boybay and Ren2017) and optic (Marchand et al. Reference Marchand, Das, Snoeijer and Andreotti2012) are used for manoeuvring the fluid flow field.
It may be mentioned here that, in comparison with other fields, utilization of a magnetic field for controlling the droplet breakup phenomenon has some serious advantageous features. The magnetic field does not induce any changes in the flow field such as pH, ionic concentration and surface charge. The magnetofluidic-based droplet manipulation is usually realized with the help of a smart fluid, known as a ferrofluid. A ferrofluid is a colloidal suspension of ferro/ferri magnetic particles in a non-magnetic carrier medium (Rosensweig Reference Rosensweig1984; Odenbach Reference Odenbach2002). It is worth adding here that a ferrofluid exhibits a superparamagnetic nature i.e. on application of a magnetic field, its magnetization is comparable to any ferromagnetic particles, whereas, on removal of the applied field, the ferrofluid flow field does not display any net hysteresis (Rosensweig Reference Rosensweig1984). Ferrofluids have been successfully used in many engineering applications such as separation (Hejazian, Li & Nguyen Reference Hejazian, Li and Nguyen2015), heat transfer augmentation (Shyam et al. Reference Shyam, Mehta, Mondal and Wongwises2019), mixing (Zhu & Nguyen Reference Zhu and Nguyen2012; Kitenbergs et al. Reference Kitenbergs, Tatulcenkovs, Erglis, Petrichenko, Perzynski and Cebers2015), droplet generation (Tan et al. Reference Tan, Nguyen, Yobas and Kang2010), breakup (Bijarchi & Shafii Reference Bijarchi and Shafii2020; Bijarchi et al. Reference Bijarchi, Dizani, Honarmand and Shafii2021) and many more (Rinaldi et al. Reference Rinaldi, Chaves, Elborai, He and Zahn2005; Cunha et al. Reference Cunha, Siqueira, Cunha and Oliveira2020). Although breakup of droplet in a T-junction has been widely explored, studies pertaining to ferrofluid droplet breakup under the modulation of a magnetic field are sparse. On closer scrutiny of the available literature, it is found that researchers have investigated the implication of a uniform magnetic field for the ferrofluid droplet splitting phenomenon (Wu et al. Reference Wu, Fu, Ma and Li2013, Reference Wu, Fu, Ma and Li2014; Li et al. Reference Li, Wu, Wang, Zhu, Fu and Ma2016; Ma et al. Reference Ma, Fu, Zhang, Zhu, Ma and Li2017). It has been shown both numerically and experimentally that, under the modulation of a uniform magnetic field, the splitting phenomenon is mostly symmetric. Also, from the reported analysis in this paradigm, it is apparent that the application of a uniform magnetic field brings in sufficient control on the size of the sister droplets generated from the mother droplet.
Albeit several underlying issues of the implication of a uniform magnetic field for the splitting of the ferrofluid droplet, is well explored (Wu et al. Reference Wu, Fu, Ma and Li2013, Reference Wu, Fu, Ma and Li2014; Li et al. Reference Li, Wu, Wang, Zhu, Fu and Ma2016; Ma et al. Reference Ma, Fu, Zhang, Zhu, Ma and Li2017), research endeavours with an emphasis on exploring the implication of the non-uniformity of the magnetic field distribution for the overall droplet breakup phenomenon are sparse. It may be mentioned here that a uniform magnetic field is a mathematical assumption, while in a realistic physical scenario, a certain amount of non-uniformity is bound to exist pertaining to magnetic field driven fluidic applications. A few researchers, however, have explored the role of non-uniformity of the applied magnetic field in the overall size of the sister droplets generated in the process and have discussed the asymmetric splitting phenomenon attributed to the onset of the involved forcing gradient (Aboutalebi et al. Reference Aboutalebi, Bijarchi, Shafii and Kazemzadeh Hannani2018; Roodan et al. Reference Roodan, Gómez-Pastora, Karampelas, González-Fernández, Bringas, Ortiz, Chalmers, Furlani and Swihart2020; Bijarchi et al. Reference Bijarchi, Dizani, Honarmand and Shafii2021).
It is worth mentioning here that the literature exploring the effect of inequality of the force field gradient on the underlying breakup dynamics in the purview of droplet train flow are sparse. Due to the non-uniform distribution of magnetic flux density, the sister droplets moving downstream of the left/right branch of the T-junction after the breakup phenomenon will possess different flow time scales. The imbalance of this flow time scale in the branches is expected to affect the overall droplet splitting phenomenon and may result in generation of sister droplets of unequal sizes. The non-uniform force field alterations in droplet splitting/breakup are expected to be more fascinating in the presence of droplet train flow, attributed primarily to the imbalance between the involved flow time scales. This time scale aspect of the droplet splitting phenomenon, albeit interesting from a fluid dynamics point of view and of huge practical relevance in different applications, has not been studied in the literature to date.
In the present investigation, we present a novel way of controlling the droplet breakup phenomenon in a T-junction divergence of an Lab on a Chip (LOC) device, under the modulation of a non-uniform magnetic field. Droplets are generated in a T-junction of the LOC device, following which they are split in another T-junction divergence, located further downstream. A magnet is placed adjacent to the left branch of the T-junction divergence, thereby inducing a high magnetic flux gradient in the vicinity of the left divergence (of the T-junction). The main proposition of the present work is to establish a non-uniform distribution of the magnetic field over the T-junction divergence of the LOC device, which, in turn, will bring in an asymmetric distribution of the flow time scale of the sister droplets moving through the left and right divergences, respectively. As established through this analysis, judicious control of the flow time scale by altering the applied forcing offers controllability to independently manoeuvre the sister droplet size.
In what follows, we divide this study into four sections. In the first section, we explore the droplet breakup phenomenon in the presence of a non-uniform magnetic field. In this section, we discuss the droplet splitting phenomenon under the modulation of an external force field obtained from both experimental investigations and numerical simulations. In the second section, we numerically explore the internal hydrodynamics of an isolated droplet during splitting under the modulation of a non-uniform magnetic forcing. In the subsequent section, we experimentally investigate the morphological evolution of the droplet splitting characteristics. In the final section of this article, we attempt to develop a physical understanding of the typical droplet breakup behaviour by exploiting the various involved time scales. It may be mentioned here that the prime motivation of the numerical simulations is to establish the proof of concept of the hypothesis and also to have a qualitative understanding of the flow field in the presence of spatially non-uniform magnetic forcing.
2. Materials and methods
2.1. Experimental methods
In the present study, we use a ferrofluid solution as the dispersed phase, while silicon oil is used as the continuous phase. We employ a co-precipitation method for the preparation of the ferrofluid solution. The prepared dispersed phase, i.e. the ferrofluid solution, is composed of de-ionized water as the carrier phase, while iron-oxide nanoparticles form the suspended phase. Interested readers are referred to our recent articles where the preparation of the ferrofluid solution has been discussed (Shyam et al. Reference Shyam, Asfer, Mehta, Mondal and Almutairi2020a,Reference Shyam, Yadav, Gawade, Mehta, Mondal and Asferc; Shyam, Mondal & Mehta Reference Shyam, Mondal and Mehta2020b, Reference Shyam, Mondal and Mehta2021). Figure 1 depicts the characterization of the prepared ferrofluid sample. The ferrofluid solution, having a saturation magnetization 
$({M_s})$ of around 10 500 A m−1, exhibits superparamagnetic characteristics, as can be observed from the 
$\boldsymbol{M}\unicode{x2013} \boldsymbol{H}$ curve of figure 1, where 
$\boldsymbol{M}$ is magnetization and 
$\boldsymbol{H}$ is magnetic field intensity. It may be mentioned here that the experimental 
$\boldsymbol{M}\unicode{x2013} \boldsymbol{H}$ curve is measured using a vibrating sample magnetometer (VSM). Also, as witnessed in figure 1, the 
$\boldsymbol{M}\unicode{x2013}\boldsymbol{H}$ curve could be accurately predicted following Langevin's theory. As per the Langevin's theory, the magnetization 
$\boldsymbol{M}$ of ferrofluids could be accurately modelled using the following relations:
\begin{gather}\boldsymbol{M} = {M_s}L({\propto} )\frac{\boldsymbol{H}}{{|\boldsymbol{H}|}},\end{gather}
\begin{gather}L({\propto} ) = \textrm{coth} \propto{-} \frac{1}{ \propto },\end{gather}
\begin{gather}\propto{=} \frac{{{\mu _0}{M_d}{V_p}H}}{{{K_b}T}},\end{gather}
where 
${M_s}$, 
$L({\propto} )$, 
${\mu _0}$, 
${M_d}$, 
${V_p}$, 
${K_b}$, T are the saturation magnetization, Langevin function, permeability of vacuum, domain magnetization, volume of the magnetic core, Boltzmann constant and temperature, respectively. Important to mention here that the saturation magnetization is given by 
${M_s} = \beta {M_d}$ (where 
$\beta $ is the volume fraction). A very good match between the data predicted by VSM and Langevin's modelling approach, justifying the superparamagnetic behaviour of the ferrofluid solution (see figure 1). We show, in the left and right insets of figure 1, the sizes of the suspended nanoparticles and the variation of the zeta potential of the ferrofluid solution, respectively. Note that the ferrofluid solution has a zeta potential of around 
$- 53\ \textrm{mV}$, signifying an electrostatically stable solution (Xu Reference Xu2002).
Figure 1. Plot depicts the magnetization curve of the prepared ferrofluid sample measured experimentally and also calculated using the Langevin function. The prepared ferrofluid sample exhibits no hysteresis and is superparamagnetic in nature. Left inset illustrates the size distribution of the suspended magnetic nanoparticles in the ferrofluid sample. The average particle size was found to be around 
$13\ \textrm{nm}$. Right inset shows the electrostatic potential characteristics of the ferrofluid solution.
The density and viscosity of the ferrofluid solution are calculated to be around 
$1050\ \textrm{kg}\ {\textrm{m}^{ - 3}}$ and 
$0.00106\ \textrm{Pa}\ \textrm{s}$, respectively. While the volume fraction of the iron nanoparticles in the ferrofluid solution was around 
$3\,\%$. As already mentioned, we use silicon oil (make: Sigma-Aldrich) as the continuous phase in the present study. Accordingly, the density and viscosity of silicon oil are found to be around 
$930\ \textrm{kg}\ {\textrm{m}^{ - 3}}$ and 
$0.3\ \textrm{Pa}\ \textrm{s}$. The interfacial tension between silicon oil and ferrofluid, as measured by tensiometer (Make: Kyowa), is found to be around 
$0.012\ \textrm{N}\ {\textrm{m}^{ - 1}}$. Note that the iron-oxide nanoparticles particles are coated with surfactant (lauric acid for the present case) essentially to avoid any agglomeration, which is not unlikely to occur due to the interparticle interactions (Rosensweig Reference Rosensweig1984; Odenbach Reference Odenbach2002; Shyam et al. Reference Shyam, Mehta, Mondal and Wongwises2019). The presence of surfactant in the ferrofluid solution lowers its static contact angle on a rigid substrate, as can be clearly observed from figure 2(a).
Figure 2. (a) Contact angle of a typical ferrofluid droplet in a glass and polydimethylsiloxane (PDMS) substrate. Ferrofluid droplet has a contact angle
$(\theta )$ of 
$44^{\circ}$ and 
$60^{\circ}$ on a glass and PDMS substrate, respectively. (b) Schematic representation of the working mechanism of the proposed microfluidic platform for controlled droplet splitting. Droplets, on being generated in a T-junction, break down in a T-junction divergence under the influence of a non-uniform magnetic field, leading to asymmetric splitting. The sequence of operation in the presence and absence of a magnetic field is schematically represented. (c) Graphical representation of the inverted microscope in which the experiments are carried out. (d) Schematic of the numerically simulated two-dimensional computational domain.
The LOC-based device, consisting of a T-shaped fluidic channel, is fabricated by using a soft lithographic technique (Whitesides & Stroock Reference Whitesides and Stroock2001). The schematic representation of the microfluidic device is shown in figure 2(b). The microfluidic passage has a square cross-section of around 
$100\ \mathrm{\mu }\textrm{m}$ as width. The fabricated device has three sections: a droplet generation junction (T-junction), a straight microchannel and a T-junction divergence (in which the droplet splitting takes place). The continuous phase (silicon oil) and the dispersed phase (ferrofluid) are injected from the two inlets leading to droplet formation at the T-junction (i.e. droplet generation junction). The generated droplet, henceforth referred to as the mother droplet, then flows through the straight microchannel and further breaks down into smaller droplets (henceforth, these droplets will be referred to as sister droplets) at the T-junction divergence, as can be seen from figure 2(b). A neodymium iron boron magnet (NdFeB) is placed adjacent to the left branch (of the T-junction divergence) for the induction of the non-uniform magnetic field. The presence of the uneven magnetic field gradient ensures asymmetric droplet splitting, as can be observed from figure 2(b).
We use a Gaussmeter (Make: SES instruments) for the measurement of the magnetic flux density. The experiments are conducted in an inverted microscope (Make: Leica), a schematic of which is shown in figure 2(c). We employ high-speed imaging for recording of the droplet splitting phenomenon. We capture images of resolution 
$1920 \times 1280\ \textrm{pixel}{\textrm{s}^2}$ at a frame rate of 1000 f.p.s. The recorded images are further post-processed in the Matlab® platform by using an in-house developed code. Interested readers are referred to one of our previous works wherein the involved steps of the image processing have been discussed at length (Shyam et al. Reference Shyam, Mondal and Mehta2020b).
2.2. Numerical analysis
2.2.1. System description
We also make an effort to investigate the droplet splitting event numerically, essentially to understand the several intricate physical aspects involved with the underlying phenomenon. Although the experimental flow dynamics pertains to that of a droplet train flow scenario, through the numerical simulations, the dynamics of an isolated droplet in an immiscible liquid medium is explored. This particular exercise saves computational time significantly. Moreover, the prime intention of carrying out numerical simulation is to explore the isolated droplet breakup dynamics in the presence of a non-uniform magnetic field; as such, the inferences can be extrapolated to the present experimental scenario without compromising the flow physics of interest. In figure 2(d), we show a schematic of the computational domain consisting of a symmetric T-junction divergence and an isolated droplet surrounded by the continuous liquid. The isolated dispersed phase is driven forward by the continuous liquid. The droplet experiences breakage at the T-junction divergence. Downstream of the T-junction (left branch), the non-uniform magnetic field is applied, as can be observed from figure 2(d). The physical properties of the dispersed (denoted by suffix 
$d$) and continuous (denoted by suffix 
$c$) fluids are density, 
${\rho _d}$, 
$\; {\rho _c}$; the viscosity, 
${\eta _d}$, 
$\; {\eta _c}$; the magnetic permeability, 
${\mu _d}$, 
$\; {\mu _c}$.
To explore the role of the relaxation times on the underlying splitting phenomenon, we have also as performed an order of magnitude analysis of the Neel 
$({\tau _N})$ and Brownian relaxation 
$({\tau _B})$ time scales incurred by the ferrofluid flow pertaining to the present analysis. The Brownian and Neel relaxation times are expressed as (Rosensweig Reference Rosensweig1984; Rinaldi et al. Reference Rinaldi, Chaves, Elborai, He and Zahn2005)
\begin{equation}{\tau _B} = 3\eta {V_h}/kT;\quad {\tau _N} = \; {\tau _0}\ {\textrm{e}^{(K{V_P}/kT)}},\end{equation}
where 
$\eta ,{V_h},k,\; T,\; {\tau _0},\; K,\; {V_P}$ are the dynamic viscosity, hydrodynamic particle diameter, Boltzmann constant, temperature, Larmor frequency, anisotropy constant and volume of the magnetic particle, respectively. It may be mentioned here that the Brownian relaxation 
$({\tau _B})$ time scale denotes the time required for the whole particle to rotate and align the magnetic moment 
$\boldsymbol{m}$ along the applied field 
$\boldsymbol{H}$. Whereas the Neel relaxation 
$({\tau _N})$ time scale denotes the time period for the magnetic moment 
$\boldsymbol{m}$ inside the particle to align with 
$\boldsymbol{H}$, due to flipping of the atomic spin without any physical particle rotation.
Considering the particle hydrodynamic diameter of around 
$13\ \textrm{nm}$ (see inset of figure 1), we found 
${\tau _B}\sim {10^{ - 7}}$. Similarly, pertaining to the Neel relaxation time scale, it may be mentioned here that an anisotropic constant is dependent on the size of the particles (Rinaldi et al. Reference Rinaldi, Chaves, Elborai, He and Zahn2005). In line with seminal works reported in the literature, for the present study, where water-based ferrofluid is used, we assume 
${\tau _0}\sim {10^{ - 9}}$, 
$K = 2.4 \times {10^4}\ \textrm{J}\ {\textrm{m}^{ - 3}}$ and a surfactant layer of around 
$2\ \textrm{nm}$ (Rosensweig Reference Rosensweig1984, Reference Rosensweig1987; Mao et al. Reference Mao, Elborai, He, Zahn and Koser2011; Dalvi, van der Meer & Shahi Reference Dalvi, van der Meer and Shahi2022), following which we found 
${\tau _N}\sim {10^{ - 8}}$. Thus, we find for the present analysis 
${\tau _N} \ll {\tau _B}$, which ensures the relaxation mechanism is intrinsic in nature, i.e. mostly dominated by Neel relaxation, giving the particles superparamagnetic behaviour (Rosensweig Reference Rosensweig1984). This superparamagnetic behaviour of the particles also ensures that no lag exists between the magnetization vector and the magnetic field vector, the underlining 
$\boldsymbol{M}$ and 
$\boldsymbol{H}$ vectors are collinear in nature. As such, the magnetization considered is the equilibrium magnetization and any non-equilibrium effect is beyond the scope of the present numerical simulations (Ilg, Kröger & Hess Reference Ilg, Kröger and Hess2005).
2.2.2. Phase-field formalism
We employ the diffused interface-based phase-field method for simulating the flow dynamics of the immiscible two-phase flow systems considered in this study. It may be mentioned here that the modelling framework of the phase-field method is obtained by the minimization of the total energy of the system and is thermodynamically consistent (Cahn & Hilliard Reference Cahn and Hilliard1958). The phase-field method has been successfully used by several researchers to study the interfacial dynamics of a multi-fluid system in the presence of an external force field (Mondal et al. Reference Mondal, Ghosh, Bandopadhyay, Dasgupta and Chakraborty2013, Reference Mondal, DasGupta, Bandopadhyay, Ghosh and Chakraborty2015; Gorthi, Mondal & Biswas Reference Gorthi, Mondal and Biswas2017). The thermodynamics of a two-fluid flow system can be described by the Ginzburg–Landau free energy functional 
$F(\phi )$, which is expressed as (Jacqmin Reference Jacqmin2000)
\begin{equation}F(\phi ) = \int_\forall {\left\{ {f(\phi ) + \frac{1}{2}\gamma \epsilon |\boldsymbol{\nabla }\phi |} \right\}\,\textrm{d}\forall } ,\end{equation}
where 
$\forall $, spans over the whole fluid domain, 
$\gamma $, and 
$\epsilon $ are the interfacial tension and the interface thickness, respectively. The first term in (2.5) denotes the bulk free energy of the binary fluid system, while the second term signifies the interfacial free energy due to the presence of an interface separating the two fluids. The interaction between these two parts of the energy constructs the static equilibrium configuration of the interface. It is important to mention here that 
$f(\phi )$ models the immiscible aspect of fluid components in the chosen system. The bulk free energy density 
$(f(\phi ))$ mathematically is a double-well function and can be expressed as, 
$f(\phi ) = \gamma {({\phi ^2} - 1)^2}/4\epsilon $. It may be mentioned there that 
$f(\phi )$, achieves its minima in two stable phases i.e. the dispersed phase 
$(\phi ={-} 1)$ and the continuous phase 
$(\phi = 1)$. It is also worth mentioning here that, although numerous double-well functions exist, pertaining to phase-field modelling, the symmetric double-well function, as mentioned above, is the most used one and has been considered in many seminal works as well (Jacqmin Reference Jacqmin1999, Reference Jacqmin2000; Badalassi, Ceniceros & Banerjee Reference Badalassi, Ceniceros and Banerjee2003; Yue, Zhou & Feng Reference Yue, Zhou and Feng2010).
The minimization of the free energy of the system 
$(F(\phi ))$ along with the mass conservation of the respective phases, leads to the well-known Cahn–Hilliard equation. Moreover, the minimization of the free energy of the system also leads to the addition of a volumetric force term in the Navier–Stokes equation, as will be discussed in § 2.2.4. The spatio-temporal evolution of 
$\phi $ depends on the Cahn–Hilliard equation as given by (Cahn & Hilliard Reference Cahn and Hilliard1958, Reference Cahn and Hilliard1959; Badalassi et al. Reference Badalassi, Ceniceros and Banerjee2003)
\begin{equation}\frac{{\partial \phi }}{{\partial t}} + \boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{\nabla }\phi = \boldsymbol{\nabla }\boldsymbol{\cdot }({M_\phi }\boldsymbol{\nabla }G),\end{equation}
where 
${M_\phi }$ and G denote the interfacial mobility factor and chemical potential, respectively. Note that 
${M_\phi }$ determines the relaxation time of the interface and the time scale of the Cahn–Hilliard diffusion. The value of 
${M_\phi }$ should be chosen in such a way that the interface maintains a constant thickness and does not cause any mass loss. The chemical potential 
$(G)$ is basically the variational derivative of the free energy functional with respect to the order parameter
$(\phi )$ (Badalassi et al. Reference Badalassi, Ceniceros and Banerjee2003)
\begin{equation}\begin{aligned} G & = \delta F/\delta \phi = \partial f/\partial \phi - \gamma \epsilon {\nabla ^2}\phi \\ & = \gamma [({\phi ^3} - \phi ) - {\epsilon ^2}{\nabla ^2}\phi ]/\epsilon . \end{aligned}\end{equation}
It may also be mentioned that in the phase-field framework, any generic property, 
$\xi $ may be expressed in terms of the order parameter 
$(\phi )$, as follows (Badalassi et al. Reference Badalassi, Ceniceros and Banerjee2003):
\begin{equation}\xi = \frac{{1 - \phi }}{2}{\xi _c} + \frac{{1 + \phi }}{2}{\xi _d}.\end{equation}2.2.3. Modelling of magnetic field distribution
We calculate the magnetic field acting on the flow domain by solving the Maxwell equations as given by (Griffiths Reference Griffiths2017)
\begin{gather}\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{B} = 0,\end{gather}
\begin{gather}\boldsymbol{\nabla } \times \boldsymbol{H} = 0,\end{gather}
where 
$\boldsymbol{B}$ is the magnetic flux density and 
$\boldsymbol{H}$ is the intensity of the magnetic field. Also, we shall mention here that the magnetic stress tensor can be expressed in the following form (Rosensweig Reference Rosensweig1984; Strek Reference Strek2008):
\begin{gather}{\boldsymbol{F}_m} = \boldsymbol{\nabla }\boldsymbol{\cdot }{\boldsymbol{\mathsf{T}}_{\boldsymbol{m}}},\end{gather}
\begin{gather}{\boldsymbol{\mathsf{T}}_{\boldsymbol{m}}} = \left( {\mu \boldsymbol{HH} - \frac{{{\mu_0}}}{2}{H^2}\boldsymbol{\mathsf{I}}}\right),\end{gather}
where 
${\boldsymbol{\mathsf{T}}_{\boldsymbol{m}}}$ is the magnetic stress tensor, 
$\mu $ represents the permeability of the flow domain, which can be evaluated using the phase-field order parameter 
$\phi $, 
$\boldsymbol{\mathsf{I}}$ denotes the identity operator and 
${H^2} = |\boldsymbol{H}{|^2}$. Moreover, the magnetization 
$\boldsymbol{M}$, the magnetic field intensity 
$\boldsymbol{H}$ and the magnetic flux density 
$\boldsymbol{B}$ are related using the following constitutive relations (Jansons Reference Jansons1983):
\begin{gather}\boldsymbol{M} = \chi \boldsymbol{H},\end{gather}
\begin{gather}\boldsymbol{B} = {\mu _0}(\boldsymbol{H} + \boldsymbol{M}) = {\mu _0}(1 + \boldsymbol{\chi })\boldsymbol{H},\end{gather}
where 
${\mu _0}$, 
$\chi $ are the permeability of vacuum and magnetic susceptibility, respectively. Since the magnetic field is irrotational (i.e. 
$\boldsymbol{\nabla } \times \boldsymbol{H} = 0$), it can be expressed in terms of scaler potential, 
$\; \boldsymbol{H} ={-} \boldsymbol{\nabla }\psi $. Using these expressions in the Gauss law leads to the following equation:
\begin{equation}\boldsymbol{\nabla }\boldsymbol{\cdot }(\mu \boldsymbol{\nabla }\psi ) = 0.\end{equation}
In the present numerical simulations, the magnetic field 
$\boldsymbol{H}$ is solved using the above (2.12). The solution of (2.15) further helps in solving the magnetic force as given by (2.11)
2.2.4. Coupling of phase field and magnetohydrodynamics
We solve the continuity, the Cahn--Hilliard and Navier–Stokes equations to obtain the pressure and velocity field of the two-liquid flow system, which are expressed as (Batchelor Reference Batchelor1970; Jacqmin Reference Jacqmin1999, Reference Jacqmin2000; DasGupta, Mondal & Chakraborty Reference DasGupta, Mondal and Chakraborty2014; Mondal et al. Reference Mondal, DasGupta, Bandopadhyay, Ghosh and Chakraborty2015)
\begin{gather}\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{u} = 0,\end{gather}
\begin{gather}\rho \left[ {\frac{{\partial \boldsymbol{u}}}{{\partial t}} + \boldsymbol{\nabla }\boldsymbol{\cdot }(\boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{u})} \right] ={-} \boldsymbol{\nabla }P + \boldsymbol{\nabla }\boldsymbol{\cdot }[\eta \boldsymbol{\nabla }\boldsymbol{u} + {(\boldsymbol{\nabla }\boldsymbol{u})^{\boldsymbol{\mathsf{T}}}}] + G\boldsymbol{\nabla }\phi + {\boldsymbol{F}_{\boldsymbol{m}}}.\end{gather}
The momentum transport equation [(2.17)] couples the phase-field formalism with magnetohydrodynamics. Here, 
$G\boldsymbol{\nabla }\phi $, is the phase-field-dependent interfacial tension force, 
${\boldsymbol{F}_{\boldsymbol{m}}}$ is the magnetic body force as denoted by (2.11). The boundary conditions considered in the flow domain for the simulations are as follows: fully developed velocity is applied at the inlet. Non-viscous, pressure constraint outflow condition is maintained at the outlet; in other words, the pressure at the outlet is fixed to zero. i.e. 
$P = 0$. For the magnetic field simulations, the magnetic insulation boundary condition, i.e. 
$\boldsymbol{n}\boldsymbol{\cdot }\boldsymbol{B} = 0$, is applied at the surrounding air domain.
2.2.5. Normalization of the governing equations
In the present section, the aforementioned governing equations are non-dimensionalized to attain a profound understanding. The width 
$(l)$ of the microchannel is considered as the characteristics length scale of the flow. The following set of dimensionless variables is defined for normalization of the governing equations:
\begin{align}{u^\ast } = \frac{{\boldsymbol{u}{l^2}}}{Q},\quad {t^\ast } &= \frac{{tQ}}{{{L^3}}},\quad {P^\ast } = \frac{{P{t^4}}}{{{\rho _c}{Q^2}}},\quad {\rho ^\ast } = \frac{\rho }{{{\rho _c}}},\quad {\eta ^\ast } = \frac{\eta }{{{\eta _c}}},\quad \; {\mu ^\ast }\nonumber\\ &= \frac{\mu }{{{\mu _c}}},\quad {\nabla ^\ast } = l\nabla ,\quad \; {\boldsymbol{H}^\ast } = \frac{\boldsymbol{H}}{{{\boldsymbol{H}_0}}},\end{align}
where 
$Q,\mu $ and P denote the flow rate, relative permeability and pressure, respectively. Using the above set of non-dimensionalized variables, the governing transport equations i.e. Cahn–Hilliard–Navier–Stokes equations along with the continuity reduce to (Bai et al. Reference Bai, He, Yang, Zhou and Wang2017)
\begin{gather}\frac{{\partial \phi }}{{\partial {t^\ast }}} + {\boldsymbol{u}^\ast }\boldsymbol{\cdot }{\nabla ^\ast }\phi = \frac{3}{{2\sqrt 2 }}\frac{1}{{Pe}}{\nabla ^{{\ast} 2}}[\phi ({\phi ^2} - 1) - C{n^2}{\nabla ^{{\ast} 2}}\phi ],\end{gather}
\begin{gather}{\nabla ^\ast }\boldsymbol{\cdot }{\boldsymbol{u}^\ast } = 0,\end{gather}
\begin{align}{\rho ^\ast }\left[ {\dfrac{{\partial {\boldsymbol{u}^\ast }}}{{\partial {t^\ast }}} + {\nabla^\ast }\boldsymbol{\cdot }({\boldsymbol{u}^\ast }\boldsymbol{\cdot }{\boldsymbol{u}^\ast })} \right] & ={-} {\nabla ^\ast }{P^\ast } + \dfrac{1}{{Re}}{\nabla ^\ast }\boldsymbol{\cdot }[\eta {\nabla ^\ast }\boldsymbol{u} + {({\nabla ^\ast }\boldsymbol{u})^{\boldsymbol{\mathsf{T}}}}]\nonumber\\ & \quad + \dfrac{3}{{2\sqrt 2 }}\dfrac{1}{{CnCaRe}}[\phi ({\phi ^2} \!-\! 1) \!-\! C{n^2}{\nabla ^{{\ast} 2}}\phi ]{\nabla ^\ast }\phi \!+\! \dfrac{{B{o_m}}}{{Ca}}\boldsymbol{\nabla }\boldsymbol{\cdot }T_m^\ast . \end{align}Thus the present problem of two-liquid flow systems is characterized by the following set of dimensionless parameters:
\begin{equation}Re = \frac{{{\rho _c}Q}}{{L{\eta _c}}},\quad Cn = \frac{\epsilon }{L},\quad Ca = \frac{{\eta Q}}{{{L^2}\gamma }},\quad Pe = \frac{{\epsilon Q}}{{{M_\phi }L\gamma }}.\quad B{o_m} = \frac{{{\mu _0}\chi lH_0^2}}{\gamma }.\end{equation}
The Reynolds number 
$Re$, which interprets the relative significance of the inertia force to the viscous force; the Cahn number 
$Cn$, is the ratio of interface thickness to the characteristics length; the capillary number, 
$Ca$, is the ratio of the viscous force to the interfacial force; the interface Péclet number 
$Pe$, is the ratio of the advection of the order parameter 
$(\phi )$ to its diffusion. This interface Péclet number, specific to phase-field modelling, is defined to determine the relative strength of advection over the diffusion of the fluid components at the interfacial region. The magnetic Bond number, 
$B{o_m}$, is the ratio of the magnetic force to surface tension force (
$\chi $ is the magnetic susceptibility). It is worth mentioning here that, for the present analysis, the order of several non-dimensional numbers are considered as, 
$Re\sim O({10^{ - 3}})$, 
$Ca\sim O({10^{ - 2}})$, 
$Cn\sim O({10^{ - 2}})$, 
$Pe\sim O({10^2})$ (Kunti et al. Reference Kunti, Mondal, Bhattacharya and Chakraborty2018; Mondal & Chaudhry Reference Mondal and Chaudhry2018). The non-dimensional forms of the fluid properties are given as (Badalassi et al. Reference Badalassi, Ceniceros and Banerjee2003)
\begin{gather}{\rho ^\ast } = \frac{{1 - \phi }}{2}\; + \frac{{1 + \phi }}{2}\; {\rho _r},\end{gather}
\begin{gather}{\eta ^\ast } = \frac{{1 - \phi }}{2}\; + \frac{{1 + \phi }}{2}\; {\eta _r},\end{gather}
\begin{gather}{\mu ^\ast } = \frac{{1 - \phi }}{2}\; + \frac{{1 + \phi }}{2}\; {\mu _r}.\end{gather}
Note that 
${\rho ^\ast }$, 
${\eta ^\ast }$, 
${\mu ^\ast }$ in (2.23)–(2.25) denotes the non-dimensionalized phase-averaged density, viscosity and permeability respectively. Similarly, 
${\rho _r}$, 
${\eta _r}$, 
${\mu _r}$ in (2.23)–(2.25) indicate the ratio of the properties of the continuous phase to that of the dispersed phase. For the present numerical study: 
${\rho ^\ast }\sim O({10^0})$, 
${\eta ^\ast }\sim O({10^{ - 2}})$ and 
${\mu ^\ast }\sim O({10^1})$. The unsteady governing transport equations together with the boundary conditions in the present study are solved by using the finite element framework of COMSOL Multiphysics®. In this method, the spatial terms of the governing transport equations are initially discretized to obtain an ordinary differential equation in time, which is then time marched to obtain evolution of the flow pattern. The PARDISO solver with generalized-
$\alpha $ scheme is used for time-stepping method. The Galerkin weighted residual method, blended with higher-order up-winding and stabilized with the capability of handling cross-wind diffusion, has been deployed here for the discretization of the convection–diffusion equations. The second-order interpolation function for velocity and the first-order interpolation function for pressure are used for the calculation of gradients within the integral terms under the framework of the weak formulation. The interpolation function for pressure is one order lower than that used for velocity. In order to obtain the velocity and pressure profiles for incompressible flows, the incremental pressure correction scheme for the segregated predictor–corrector method is employed. A tolerance levels of 
${10^{ - 6}}$ for mass divergence are specified for all the numerical simulations (Mondal et al. Reference Mondal, DasGupta, Bandopadhyay, Ghosh and Chakraborty2015).
2.2.6. Grid independence study and validation
Here, we show the grid independence test to ensure the correctness of the numerical results presented in this analysis. In addition, to ensure a sharp interface limit, we also carry out the Cahn number independence test. It may be mentioned here that, in our study, we have chosen a grid size near the interface equal to the Cahn number 
$(Cn)$. Thus, in the context of the present analysis, grid independence will simultaneously indicate a Cahn number independence and vice versa. Figure 3(inset) shows the temporal evolution of the non-dimensionalized width 
$({W^\ast })$ of the droplet as it breaks up in the T-junction divergence of the LOC device, obtained for different values of grid resolution. As can be seen clearly from figure 3(a), the numerical results become independent of mesh size below 
$Cn < 0.02$. Although, for 
$Cn\sim 0.01$, a better resolution is obtained, considering the involved computational cost vis-a-vis the gained accuracy, we choose 
$Cn\sim 0.02$ for the present study. Note that this particular value of the Cahn number, 
$Cn\sim O({10^{ - 2}})$ ensures the sharp interface limits as well (Yue et al. Reference Yue, Zhou and Feng2010).
Figure 3. (a) Temporal evolution of the non-dimensionalized width 
$({W^\ast })$ of the droplet at the T-junction divergence for different grid resolutions. Grid independence test justifies Cahn number independence test; 
${{W}^{{\ast }}} = {w}({t})/{l}$, where 
${w}({t})$ and 
${l}$ denote the instantaneous width of the droplet and the width of the microchannel, respectively. Left inset snapshots show a typically resolved interface of a droplet for the Cahn numbers, 
$Cn = 0.01,\;Cn = 0.02$ and 
$Cn = 0.04$. Right inset shows the comparison between the present numerical study and (2.26).
We also compare in figure 3 (inset) our numerical results quantitatively with those obtained from the analytical relation of Bretherton (Bretherton Reference Bretherton1961). For low Reynolds number flow i.e. 
$Re \ll 1$, the velocity of a droplet moving in a slender tube with a thin film separating the droplet and the wall is given by (Bretherton Reference Bretherton1961)
\begin{equation}{\boldsymbol{u}^{\boldsymbol{\ast }}} = 1 + 1.29C{a^{2/3}}.\end{equation}We show in the inset of figure 3, an excellent agreement between our numerical results with those obtained from (2.26). A closer match between the results, seen from figure 3(inset), vouches for the correctness of the numerical modelling framework developed in this study. Besides, in an effort to further establish the credibility of the proposed numerical methodology, we have performed a separate set of benchmarking exercises pertaining to droplet breakup by using our own experimental results. However, interested readers are referred to the supplementary materials section of this article for this part, wherein benchmarking of simulated results with our experimental results is presented.
3. Results and discussions
In this section, we explore the droplet breakup phenomenon and its consequences for the alterations in the various hydrodynamical parameters of the flow field. As already mentioned, the present investigation is carried out in the flow regimes whereby the viscous and interfacial forces are dominant in comparison with the inertia forces. As such, the study is divided into four parts. In the first part, we develop a physical understanding of the ferrofluid droplet splitting mechanism in the presence of a non-uniform force field. To this end, we make use of the high-speed imaging modalities and perform numerical simulations for exploring the dynamical behaviour of the droplet under the modulation of the spatially varying force field. In the second section, we numerically explore the intricate flow physics involved with the droplet splitting dynamics. Subsequently, in the third section, we experimentally explore the morphological evolution of the mother droplet in the T-junction divergence of the LOC-device under the modulation of the non-uniform magnetic forcing. Finally, in the last part of the study, we make an attempt to develop a physical reasoning behind the characteristic behaviour of the sister droplets from the perspective of the involved time scales.
3.1. Droplet breakup: qualitative dynamics
It may be mentioned here that, depending on the capillary number 
$(Ca)$ and initial slug length 
$({l_0})$, a droplet may split into two sister droplets in the T-junction divergence following permanent obstruction or partial obstruction (Jullien et al. Reference Jullien, Tsang Mui Ching, Cohen, Menetrier and Tabeling2009; Leshansky & Pismen Reference Leshansky and Pismen2009). Moreover, based on the initial slug length and Capillary number, the droplet may not split at all (Jullien et al. Reference Jullien, Tsang Mui Ching, Cohen, Menetrier and Tabeling2009; Chen & Deng Reference Chen and Deng2017). Pertaining to the permanent obstruction case, the dispersed phase totally engulfs the T-junction divergence and blocks the flow of continuous phase liquid. On the contrary, for the partial obstruction case, a tunnel develops, which allows the motion of the continuous fluid over the dispersed phase. Therefore, in the partial obstruction case, the droplet blocks the T-junction divergence of the LOC-device incompletely. Henceforth, the word ‘droplet’ will refer to the mother droplet.
It is worth mentioning here that, for the present study, the mother droplet splits into two sister droplets following the permanent obstruction of the T-junction divergence, as can be observed from the depicted images in figure 4. We show in figure 4 the spatio-temporal evolution of the droplet splitting phenomenon in the absence of an external field. Note that the droplet breaks in three typical stages: squeezing, transition and pinch-off (Ma et al. Reference Ma, Fu, Zhang, Zhu, Ma and Li2017). Initially, as the droplet (i.e. dispersed phase) enters into the T-junction divergence, it tries to block the whole junction, as can be seen from figure 4(a). As the droplet occupies the whole junction 
$({t^\ast } = {0^ + })$, the upstream continuous phase flow gives rise to the formation of a depression in the droplet, which further leads to the change in the curvature of the neck (the circled region in 
${t^\ast } = 0.5$ in figure 4a). The depression is an effect of the squeezing pressure acting on the neck of the droplet (of the dispersed phase). The squeezing pressure is the pressure that is developed in the upstream continuous phase due to the permanent blockage of the T-Junction divergence by the dispersed phase. Since no tunnel develops for the present case (i.e. permanent obstruction case), the squeezing force becomes very large in comparison with the viscous force, thereby dictating the overall dynamics in the squeezing stage. The squeezing stage is followed by the transition stage. In the transition stage also, no tunnel formation is seen, and thus, the splitting phenomenon is dictated by the balance between the interfacial tension force and the squeezing force. Although the upstream pressure force (i.e. the squeezing force) is dominant in both squeezing and transition stages, the temporal evolution of the rear interface of the droplet is found to be different for the two individual cases, as will be discussed in detail in the forthcoming sections. The transition zone is followed by the pinch-off zone, in which the droplet gets fully detached, leading to the formation of two sister droplets, as can be observed from figure 4(a). In the succeeding discussions, we will demarcate all of these existing splitting zones appropriately.
Figure 4. Temporal evolution of the droplet splitting phenomenon observed using (a) experiments, and (b) numerical simulations. The white-coloured arrows indicate the direction of the flow. Here, 
${t^\ast } = t/{t_L}$ represents the non-dimensional time period of the droplet splitting phenomenon; 
$t,\; {t_L}$ represent the instantaneous time and the droplet splitting time, respectively. The span of the droplet splitting time starts from the moment the droplet fully enters the T-junction to the moment two dispersed sister droplets have been formed. The experiments are carried out using at a magnification of × 10 at a capture speed of 1000 f.p.s.
To attain a detailed visualization to arrive at an in-depth understanding of the experimental observation, we perform numerical simulations of the droplet splitting phenomenon. The details of the numerical methodology adopted in this analysis are already mentioned in § 2.2. As already mentioned, although the present experimental study deals with a droplet train, the splitting dynamics of an isolated droplet is investigated numerically. This particular exercise will give us qualitative insights into the droplet breakup dynamics in the absence/presence of a non-uniform magnetic field. We show in figure 4(b) the droplet splitting phenomenon in a T-junction simulated numerically. We can clearly observe the presence of various droplet splitting stages such as the squeezing, transition and pinch-off from figure 4(b). It is worth adding here that the similarity of the underlying droplet splitting dynamics between experimental observation and numerical results, observed from qualitative perspectives in figure 4(a) vis-à-vis figure 4(b) justifies the credibility of our experimental methods.
We show in figure 5(a,b) the spatio-temporal evolution of the ferrofluid droplet splitting phenomenon in the presence of a non-uniform magnetic field, obtained experimentally and calculated from numerical simulations, respectively. The presence of all the intermediate stages of the splitting phenomenon, i.e. the squeezing stage, transition stage and the necking stage, which are observed in figure 4 as well, can be clearly seen for the present case (see figure 5). The behaviour of the droplet splitting in the presence of a non-uniform magnetic field becomes similar to the case where no external force acts on the droplet flow field. The only difference is the asymmetricity produced during the droplet breakup process. Note that the asymmetricity developed in the droplet splitting process is precisely due to the involved non-uniformity in the force field gradient induced due to the applied field. Readers are referred to the supplementary materials section for a detailed distribution of the magnetic field flux density in the fluid flow domain. We can clearly observe from figure 5(a) that, as the ferrofluid droplet enters the T-junction divergence, the upstream pressure forces the advancement of the dumbbell-shaped bulge in the left and right branches. Note that, in the following discussion, we refer to these left/right dumbbell-shaped bulges of the dispersed phase as the left/right moving bulge.
Figure 5. Temporal evolution of the droplet splitting phenomenon in the presence of a non-uniform magnetic field, observed using (a) experiment, (b) numerical simulations. The white-coloured arrows indicate the direction of the flow. Here, 
${t^\ast } = t/{t_L}$ represents the non-dimensional time period of the droplet splitting phenomena; 
$t,\; {t_L}$ represent the instantaneous time and the droplet splitting time, respectively. The span of the droplet splitting time starts from the moment the droplet fully enters the T-junction to the moment two dispersed sister droplets have been formed. The experiments are carried out using a higher magnification of 10X at a capture speed of 1000 frames per second.
However, due to the non-uniform distribution of the applied magnetic flux density, the advancement (of the dumbbell-shaped bulge) becomes more aligned towards the left branch (asymmetricity is produced in the breakup phenomenon) in comparison with the right branch (refer to 
${t^\ast } = 0.5$ and 
${t^\ast } = 0.75$ of figure 5a) of the T-junction divergence. In particular, because of this reason, we observe in figure 5(a) the presence of unevenly sized sister droplets being produced in the process. As such, the larger-sized sister droplet moves in the left branch, while the smaller-sized sister droplet moves in the right branch, as can be clearly observed from figure 5(a). We show the numerically simulated droplet splitting phenomenon in figure 5(b). The asymmetric behaviour of the droplet breakup can easily be observed from figure 5(b) as well. Quite notably, the numerically simulated stages of the breakup phenomenon are in accordance with the experimental observations (see figure 5a), further justifying the credibility of our experimental observation. Moreover, it may be conjectured from figures 4 and 5 that, for a droplet split following permanent obstruction, the characteristic behaviour of the breakup pattern remains the same irrespective of whether a magnetic field is applied or not. The utilization of a non-uniform magnetic field gives rise to asymmetry in the droplet breakup phenomenon, leading to the generation of two sister droplets of unequal sizes.
3.2. Droplet breakup phenomenon: numerical perspectives
As seen from the preceding discussion, the numerical results agreed well qualitatively with the experimental observations. This endeavour provides us with the flexibility to delve deep into the droplet splitting dynamics following numerical computations. We show in figure 6(a,b) the velocity vectors inside the droplet flow field at different stages of splitting for the various cases under investigation, i.e. in the absence/presence of a non-uniform magnetic field. It may be mentioned here that 
${t^\ast } = {0^ - }$ denotes the stage of the droplet before entering the T-junction divergence. It can be clearly observed that the centre of the droplet, before fully impacting the T-junction divergence, attains a maximum velocity (see 
${t^\ast } = {0^ - }$ of figure 6a). In the absence of any external force, we can clearly observe that the droplet demonstrates two distinct velocity regimes when it impacts the T-junction divergence at 
${t^\ast } = {0^ + }$, as denoted by point A and point B in figure 6(a).
Figure 6. Temporal evolution of the droplet splitting phenomenon calculated numerically in the (a) absence and (b) presence of a non-uniform magnetic field. The white-coloured arrows indicate the direction of the flow. Here, 
${t^\ast } = t/{t_L}$ represents the non-dimensional time period of the droplet splitting phenomenon; 
$t,\; {t_L}$ represent the instantaneous time and the droplet splitting time respectively. (c) Zoomed view of the droplet splitting phenomenon at 
${t^\ast } = 0.75$ for I. with and II. without a magnet cases. The red-coloured arrow shows the direction of the flow. III. The right-most figure shows the zoomed view of the red-coloured ellipse along with the black-coloured arrows indicating the velocity vectors. (d) Experimental hump formation in the ferrofluid flow domain in the presence of a static magnet. The zoomed-in view shows the chained cluster formation of the magnetic nanoparticles in the hump-like structure. The yellow and red coloured dotted arrows indicate the direction of the ferrofluids in zones I and II, respectively.
Here, point A identifies the position whereby the droplet remains in contact with the upstream flow, point B on the other hand, denotes the position of the droplet which faces the wall of the T-junction divergence. Note that point A demonstrates the maximum velocity while point B depicts the region of minimum velocity. This is because the upstream flow attempts to squeeze the droplet, thereby ensuring that point A exhibits higher velocity. Similarly, as the droplet impacts the wall of the junction, a stagnation point exists, as indicated by point B in figure 6(a). Subsequently, the droplet, after impacting the T-junction divergence, stretches itself symmetrically in the left and right branches. This stretching leads to the migration of the local maximum points, as given by points C (refer 
${t^\ast } = 0.25$ of figure 6a). On the other hand, the position of the local minimum (point B) remains the same due to the symmetric nature of the stretching phenomenon, as can be observed from figure 6(a). It can also be observed that the magnitude of the velocity inside the droplet decreases with the temporal progression of the splitting phenomenon, as witnessed in figure 6(a). Having a closer look at figure 6(a) vis-à-vis figure 6(b), it is observed that the internal hydrodynamic behaviour inside the droplet under the influence of a non-uniform magnetic field is totally altered. The position of the maximum and minimum velocities inside the droplet are represented by points D and E. It can be clearly seen from figure 6(b) that the position of the local maximum and minimum inside the droplet are altered when compared with the corresponding points pertaining to the non-magnetic case (
${t^\ast } = {0^ + }$ of figure 6b). This altercation is due to the presence of the non-uniform magnetic flux density. The presence of a high force field gradient in the left branch ensures that the maximum (of velocity) is always aligned towards the left branch. Consequently, due to low magnetic flux density, the minimum (of velocity) is aligned towards the right branch, as can be clearly observed from figure 6(b). In particular, a high force field gradient acting in the left branch induces more ferrofluid mass to flow towards it, leading to the generation/development of large-sized sister droplets moving on the left side (refer to figure 6b). Also, the location of the stagnation point is altered due to the involved asymmetric stretching of the droplet (refer to point E of figure 6b). Quite notably, the influence of a magnetic field leads to the migration of the localized maximum and minimum to the left bulge (of the droplet), as can be seen at 
${t^\ast } = 0.5$ and 
$0.75$, respectively, from figure 6(b). In figure 6(c), we show the zoomed-in view of the droplet domain at 
${t^\ast } = 0.75$ for both the cases i.e. with/without magnet cases, essentially to draw clear distinctions between them. The presence of maxima and minima in the left moving bulge (of the droplet) under the modulation of a magnetic field can be clearly observed in figure 6(c-II). Note that the fluid facing the lower wall (i.e. the wall of the left branch facing the magnet) has a low magnitude of velocity in comparison with the fluids facing the upper wall. As such, two clear distinct zones can be observed in the left moving bulge (of the droplet) under the modulation of a non-uniform magnetic field.
To explore the implication of these typical zones, in figure 6(c-III) we plot a further zoomed-in view of the left moving bulge of the dispersed phase (i.e. the droplet) under modulation of a magnetic field (at 
${t^\ast } = 0.75$). We demarcate the fluid velocity magnitude inside the left moving bulge into two zones i.e. zone I and zone II, respectively, as shown by the pink coloured ellipse in figure 6(c-III). As is clearly visible, zone I, which is adjacent to the lower wall, experiences a velocity of lower magnitude. On the contrary, zone II, which is adjacent to the upper wall, experiences velocities of higher magnitude. Note that the lower wall experiences a higher magnetic field flux density in comparison with the upper wall. Therefore, it may be argued that, when the ferrofluid droplet reaches the T-junction divergence, the magnetically susceptible fluid rushes towards the lower wall in an attempt to approach the magnet. However, on reaching the high magnetic field flux zone, the ferrofluids (inside the left bulge) motion is highly restrained. This restrained motion of the fluid towards the magnet leads to the development of a ‘hump-like structure’ near the lower wall, as can be observed from figure 6(c-III). Due to the high magnetic field gradient existing nearby the lower wall, the magnitude of the velocity is low inside this hump-like structure, as can be clearly observed from figure 6(c-III). It may also be mentioned here that, owing to the presence of the hump-like structure, the succeeding ferrofluids (inside the left moving droplet volume) are forced towards the upper wall to satisfy the mass conservation constraint, as can be observed in zone II from the right-most snapshots of figure 6(c-III). Due to the relatively low force field gradient near the upper wall, a high velocity magnitude is encountered near the upper wall, as can be clearly seen from figure 6(c-III). Following this, it may be argued that the presence of the ‘hump-like structure’ non-trivially ensures augmented convections inside the left-directed bulge of the droplet, thereby endorsing its amplified velocity in comparison with the bulge advancing in the right branch.
In the previous discussion, we have argued the presence of a ‘hump-like structure’ inside the left moving ferrofluid droplet under the modulation of a magnetic field. To check the validity of this argument, we carried out a bright field visualization of the ferrofluid flow field in a microchannel under the influence of a static magnet. We show in figure 6(d), the ferrofluid flow field in the presence of a steady magnet placed near the lower wall. The presence of a hump-like cluster, together with the formation of nanoparticle agglomeration following a chain-like structure, can be clearly observed from figure 6(d). This observation justifies our arguments for the enhanced flow convections in the left moving bulge (see movie 4 for a clearer insight). We will show in the succeeding sub-sections that this augmented velocity in the bulge (moving in the left branch) can be significantly maneuvered in controlling the size of the generated droplet (sister droplet) particularly, for the cases involving the droplet train.
As discussed in § 3.1, when a droplet splits into sister droplets following the permanent obstruction, the balance among the upstream pressure force, magnetic force and the interfacial tension force dictates the overall splitting/breakup phenomenon. Accordingly, in order to have a comprehensive understanding on the competition of these three forces, we show in figures 7 and 8, the detailed evolution of the pressure distribution and the Laplace pressure drop across the rear droplet interface during the breakup phenomenon with permanent obstruction. The discussion pertaining to this aspect follows the results obtained from the numerical simulations performed in this analysis. It may be mentioned here that, due to the upstream pressure, the curvature of the rear interface of the droplet undergoes a temporal change from a convex to a concave shape. Therefore, we represent the entire phenomenon in figures 7 and 8, by depicting three typical stages i.e. when the droplet rear interface is convex, flat and concave (see inset of individual figures for the snapshot of the temporal instant). The order parameter distribution in figures 7 and 8 helps in identifying the rear interface as it varies between 
$\phi ={-} 1$ and 
$\phi = 1$. Accordingly, the pressure difference across the interface indicates the Laplace pressure 
$({P_c} = {P_{inside}} - {P_{outside}})$.
Figure 7. Pressure and order parameter distribution along the symmetric axis of the T-junction around the interface for the droplet in the (a) entering phase, (b) at the junction, (c) after entering the junction, in the absence of any external force field; 
${P_c} = {P_{inside}} - {P_{outside}}$, where 
${P_c}$ is the Laplace pressure, 
${P_{inside}}$ is the pressure inside the droplet (dispersed phase), 
${P_{outside}}$ is the pressure inside the upstream continuous phase, 
$\phi $ is the order parameter. Here, 
$\bar{P} = {P^\ast }/{P_{inlet}}$, 
${P^\ast }$ is the instantaneous non-dimensionalized pressure and 
${P_{inlet}}$ is the non-dimensionalized pressure at the inlet of the channel. Negative value indicates that the surface tension is oriented upstream. The blue shaded area indicates the branched channel area. The inset indicates the respective spatio-temporal location of the droplet.
Figure 8. Pressure and order parameter distribution along the symmetric axis of the T-junction around the interface for the droplet in the (a) entering phase, (b) at the junction, (c) after entering the junction, in the presence of a non-uniform external magnetic field; 
${P_c} = {P_{inside}} - {P_{outside}}$, 
${P_c}$ is the Laplace pressure, 
${P_{inside}}$ is the pressure inside the droplet (dispersed phase), 
${P_{outside}}$ is the pressure inside the upstream continuous phase, 
$\phi $ is the order parameter. 
$\bar{P} = {P^\ast }/{P_{inlet}}$, 
${P^\ast }$ is the instantaneous non-dimensionalized pressure and 
${P_{inlet}}$ is the non-dimensionalized pressure at the inlet of the channel. Negative value indicates that surface tension is oriented upstream. The blue shaded area indicates the branched channel area. The inset indicates the respective spatio-temporal location of the droplet.
It may be mentioned here that, in the absence of an external force field, when the droplet enters the T-junction, its motion becomes restrictive due to the presence of the wall of the junction. As a consequence of this restricted motion, the droplet undergoes deformation symmetrically. At this juncture, the rear interface of the droplet exhibits a convex profile, as can be observed from the inset of figure 7(a). Due to this convex profile, the capillary pressure, 
${P_c}( = {P_{inside}} - {P_{outside}}) \gg 0$, demonstrates a positive value (see figure 7a). With the progression of temporal instances, the droplet leaves the main channel fully and occupies the whole junction, and at this stage, the upstream pressure forces the rear interface of the droplet to attain a flat profile (
${t^\ast } = 0$, of figure 6a). It is due to this flat profile that the capillary pressure assumes an almost negligible value i.e. 
${P_c}( = {P_{inside}} - {P_{outside}})\sim 0$, as can be observed from figure 7(b). The upstream pressure further forces the neck of the droplet to attain a concave profile, thereby ensuring that the capillary pressure attains a negative value i.e. 
${P_c}( = {P_{inside}} - {P_{outside}}) < 0$, as can be observed from figure 7(c). Note that the dynamical evolution of the pressure remains the same qualitatively even in the presence of a non-uniform magnetic field, i.e. the rear interface of the droplet evolves from the convex to the concave shape as observed from figure 8(a–c). Up to this section, we correlate the insight gained from the numerical simulations directly to the experimental scenarios. Following these inferences, which will be used to support the experimental observations, we explore the dynamics of the droplet breakup events focusing on experimental results in the succeeding sections.
3.3. Droplet breakup: evolution of the droplet morphology
3.3.1. Evolution of droplet width
We have seen from the previous discussion that, on application of a non-uniform magnetic field, the ferrofluid droplet tends to get stretched in the direction of the applied magnetic field. We have also observed from both experimental observations and numerical simulation that the uneven stretching of interface results in an asymmetric splitting of the droplet. In figure 9, we show the variation of the non-dimensional thickness 
$({W^\ast } = w/l)$ of the droplet (dispersed phase) as it breaks both in the presence and absence of a magnetic field. It may be mentioned here that the underlying event of droplet splitting is mainly governed by the intricate competition among the interfacial force, viscous force and the magnetic force. Here, we introduce the magnetic Bond number 
$(B{o_m})$, which is used to represent the relative strength between the magnetic force and the surface tension force. The time zero indicates the moment at which the dispersed phase entirely penetrates into the T-Junction divergence. As such, the rear interface of the dispersed phase is almost flat at 
$t = 0$. We demarcate, in figure 9, the various regimes of splitting encountered by the droplet i.e. the squeezing regime, the transition regime and the pinch-off regime. The demarcation of the various regimes is identified out following the characteristic slopes exhibited by the deforming ferrofluid droplet in the T-junction divergence, as can be observed from figure 9. It may be mentioned here that evolution of the neck thickness of the droplet is such that, in the initial squeezing regime, a linear gradient is observed. Important to mention here is that this linear variation of 
${W^\ast }$ is attributed to the typical role played by the two responsible forces i.e. the interfacial force and the upstream pressure force, on the underlying phenomenon in this regime. Due to the limited deformation of the rear interface of the droplet in the squeezing regime, the role played by the interfacial force in the pertinent regime becomes miniscule. As a result, the rear interface of the droplet moves primarily due to the squeezing force of the upstream flow, with limited to no role played by the surface tension force. It is because of this force balance, a linear relationship between 
${W^\ast }$ vs t, can be observed during the squeezing regime as observed from figure 9. This is followed by the transition regime, in which the interfacial tension force resists the deformation of the rear interface. This resistance leads to a delay in the overall deformation of the rear interface of the droplet. Consequently, we observe an exponential relationship between 
${W^\ast }$ and t, as witnessed in figure 9 (see regime B). The transition regime is followed by the pinch off stage, in which the thinning rate is significantly amplified with time.
Figure 9. Temporal evolution of the non-dimensionalized droplet width 
$({W^\ast })$ in absence and presence of the magnetic field. The circles indicate the case when 
$B{o_m} = 0$, the squares indicate the case when 
$B{o_m} = 6.1$. The blue, black and red lines show the squeezing, transition and pinch-off regimes. The inset shows a schematic of the width of a typical droplet breakup process under a permanent obstruction case. Here, 
${W^\ast } = w(t)/l$, where 
$w(t)$ and l denote the instantaneous width of the droplet and the width of the microchannel; 
$R = 1$, 
$Ca ={\sim} {10^{ - 3}}$, 
$Re\sim {10^{ - 3}}$. The top inset shows the variation of the thinning rate 
$1 - {W^\ast }$ when a magnetic field is applied.
The top inset of figure 9 shows the variation of the thinning rate, 
$1 - {W^\ast }$ vs 
${t^\ast }$, for the case when a magnetic field is applied. As already mentioned, in the squeezing regime, it can be observed that, 
$1 - {W^\ast }\sim {t^\ast }$, and this is essentially due to the minimal involvement of the interfacial tension force. While, for the transition regime, we observe, 
$1 - {W^\ast }\sim {t^{{\ast} 3/7}}$, where the scaling 
$3/7$, agrees well with the theoretical solutions as proposed by Leshansky and Pismen (Leshansky & Pismen Reference Leshansky and Pismen2009). Therefore, it can be argued that the characteristic behaviour of the droplet splitting, even as observed in the regime of permanent obstruction, is independent of the fact of whether a magnetic field is applied or not. However, it can be clearly visualized from figure 9 that magnetic field has a role to play in the overall lifetime of the droplet breakup phenomenon. This further implies that, by specifically tuning the force field, we could effectively control the size of the droplet. In the forthcoming sections, we comprehensively discuss effective ways in which the size of the sister droplets can be controlled.
3.3.2. Effect of magnetic flux density
We show in figure 10, the variation of non-dimensionalized width 
$({W^\ast })$ for the various magnetic flux densities (for various 
$B{o_m}$). It can be clearly observed from figure 10 that the characteristic variation of the droplet width is the same irrespective of the applied magnetic field strength. Also, figure 10 demonstrates that the required time for splitting of the ferrofluid droplet into the sister droplets is dependent on 
$B{o_m}$. As such, it could be observed from figure 10 that the droplet splitting time 
$(t)$ varies in the following sequence: 
$B{o_m}( = 0) < B{o_m}( = 7.5) < B{o_m}( = 6.1) < B{o_m}( = 15.1)$. In the next sections, we explore this particular aspect of the droplet splitting time as modulated by the magnetic flux density in greater detail.
Figure 10. Temporal evolution of the non-dimensionalized droplet width 
$({W^\ast })$ for the various magnetic field Bond numbers 
$(B{o_m})$. Here, 
${W^\ast } = w(t)/l$, where 
$w(t)$ and l denote the instantaneous width of the droplet and the width of the microchannel. The right side inset shows the variation of the non-dimensionalized length 
$(l_l^\ast )$ of the sister droplet flowing in the left branch of the microchannel, for the various magnetic field strengths. The corresponding snapshots of the generated droplet are shown in the left-bottom inset marked by 1, 2, 3, 4. Here, 
$l_l^\ast = l(t)/{l_0}$, where 
$l(t)$ and 
${l_0}$ denotes the instantaneous length of the sister droplet and characteristic length of the mother droplet; 
$R = 1$, 
$Ca ={\sim} {10^{ - 3}}$, 
$Re\sim {10^{ - 3}}$.
The top right inset of figure 10 shows the variation of non-dimensionalized length 
$(l_l^\ast )$ of the sister droplet flowing in the left branch. It is important to mention that the magnet is placed adjacent to the left branch of the T-junction divergence. Therefore, to explore the role of magnetic field on the droplet splitting and, subsequently, on the change in the length of the sister droplets, we focus our attention on 
$l_l^\ast $. Note that 
$l_l^\ast > 0.5$ ensures that the length of the sister droplet moving to the left branch 
$(l_l^\ast )$ is higher in comparison with the length of the sister droplet moving to the right branch 
$(l_r^\ast )$. In other words, it can also be said that the application of a non-uniform magnetic field ensures that an increased fluid volume moves to the left branch when compared with the volumes of fluid moving to the right branch, precisely due to the existence of the high field gradient (in the left branch). We can clearly observe from figure 10 (right-inset) that the curve exhibits a positive slope and, on reaching a particular threshold value, it exhibits a negative slope. Note that, with an initial increase in 
$B{o_m}$, the change in length of the sister droplet in the left branch 
$(l_l^\ast )$ increases and on reaching a particular critical 
$B{o_m}$, the length of the sister droplet 
$(l_l^\ast )$ decreases. Any increase in 
$B{o_m}$ beyond the critical value reduces the size of the sister droplet moving in the left branch. The left-side inset (of figure 10) shows the corresponding snapshots of the droplet splitting phenomenon for the different cases under investigation. It can be clearly observed from the left-side inset of figure 10 that the length of the sister droplet 
$(l_l^\ast )$ is maximum for 
$B{o_m} = 7.5$.
3.3.3. Effect of flow ratio 
$({R})$
In this subsection, we explore the effect of the dispersed phase flow rate on the droplet splitting phenomenon. We have mentioned that, in the present study, the flow rate of the dispersed phase is changed while that of the continuous flow is kept constant. This particular exercise is carried out essentially to vary the length of the generated droplet (equivalently the volume), keeping the continuous phase flow velocity constant. We normalize the flow rate as the flow ratio, given as, 
$R = {Q_c}/{Q_d}$. Note that lower the value of R, higher is the slug length of the dispersed phase (i.e. the droplet). We show in figure 11 the variation of the droplet breakup time for the various values of the flow ratio, R. We can clearly see that on increasing the value of the flow ratio 
$(R)$, the splitting time decreases, and any increase in the flow ratio 
$(R )$ beyond the critical value increases the splitting time. Thus, it may be argued that there exists a particular ferrofluid slug length at which the splitting time is minimum. For the present case, the droplet experiences a minimum splitting time for 
$R = 0.66$. In the inset of figure 11, we show the variation of 
${W^\ast }$ for the different values of 
$B{o_m}$ obtained at 
$R = 0.66$. Although the characteristic behaviour of the change in width of the droplet is similar, a clear distinction between the overall droplet splitting time is evident for varying magnetic field strength (as reflected by the change in 
$B{o_m}$). The minimum splitting time can be observed for 
$B{o_m} = 7.5$, as witnessed by the inset of figure 11.
Figure 11. Droplet splitting time 
$({t^\ast })$ for the various magnetic Bond numbers 
$(B{o_m})$. The red-coloured ellipse indicates the flow ratio at which the ferrofluid droplet exhibits minimum splitting time. The inset shows the variation of non-dimensionalized width 
$({W^\ast })$ of the droplet for the various magnetic field intensities. Here, 
$Ca ={\sim} {10^{ - 3}}$, 
$Re\sim {10^{ - 3}}$.
We show in figure 12, the length of the sister droplet 
$(l_l^\ast )$ moving in the left branch after the breakup phenomenon. Needless to mention here is that the influence of a non-uniform magnetic field is significant in the left branch. As can be clearly observed from figure 12, there exists a threshold magnetic field strength (precisely 
$B{o_m} = 7.5$) and flow ratio 
$(R)$ at which 
$l_l^\ast $ is maximum. As can be observed from figure 12, any increase in the magnetic field strength beyond the threshold value decreases the length of the sister droplet 
$(l_l^\ast )$ that is migrating in the left branch. The corresponding snapshots of the events pertaining to sister droplet generation for the various cases under investigation can be observed from the inset of figure 12. The consequence of this particular insight is significant, since by tunning the strength of the magnetic field and the flow ratio, we could ensure the desired volume of slug to move into the various branches of the T-junction divergence (precisely the left/right branch). In the following section, we explore the physical reason behind this typical behaviour of the sister droplet train in a confined microfluidic passage in the presence of a non-uniform magnetic field.
Figure 12. Variation of the non-dimensional length of the sister droplet moving in the left branch 
$(l_l^\ast )$ for the various magnetic field Bond numbers 
$(B{o_m})$. Blue-coloured arrow indicates the position of the threshold 
$B{o_m}$ at which 
$l_l^\ast $ is maximum. Here, 
$R = {Q_c}/{Q_d}$, where 
${Q_c}$, 
${Q_d}$ are the flow rates of the continuous phase and the dispersed phase, respectively. The inset shows the snapshots of the sister droplet for the various cases under consideration. Here, 
$Ca ={\sim} {10^{ - 3}}$, 
$Re\sim {10^{ - 3}}$.
3.4. Mechanism of splitting: a time scale perspective
The results until now have shown that there exist a threshold magnetic field strength 
$(B{o_m})$ and flow ratio 
$(R)$ for which the size of the generated sister droplet in the left branch becomes maximum. In this section, we attempt to unearth the reason behind this typical splitting behaviour of droplets in the presence of a non-uniform magnetic field. As mentioned before in the present work, splitting takes place by virtue of permanent obstruction. Therefore, the underlying motion of the droplet (being split) moving in the left branch and the right branch simultaneously will dictate the overall splitting process.
To develop an overall understanding of the splitting phenomenon, two distinguishable parameters were identified: the velocity of the sister droplets advancing towards the left 
$({u_l})$, and the right divergences 
$({u_r})$, simultaneously. It is found that, due to the non-uniformity in the distribution of the magnetic field flux density, there exists a significant difference between 
${u_l}$ and 
${u_r}$. It is worth mentioning here that, owing to this effect, the flow time scale of the sister droplets’ motions in the left branch 
$({t_l} = l/{u_l})$ and in the right branch 
$({t_r} = l/{u_r})$ of the microfluidic channel become different.
The flow time scale refers to the time required by the individual sister droplet to travel the characteristics length 
$(l)$ in the respective branches (left/right) of the T-junction divergence. It may be reiterated that, in the absence of any external forcing, the dispersed phase (the ferrofluid droplet) on impacting the T-junction divergence (i.e. the splitting junction) blocks the whole junction and stretches itself symmetrically in the left and right branches. During the stretching phase, the dumbbell-shaped bulges of the droplet (being split) move in both the left and right branches independently. This migration event eventually culminates in the breaking of the droplet into two sister droplets. As a consequence of the symmetric motion of the sister droplets (in the absence of any external forcing), the respective time scales are also equal. Whereas, under the motion of a non-uniform magnetic forcing, the time scale of motion of the respective sister droplet (i.e. the sister droplets migrating towards the left and right branch ) is different. Due to the non-uniform distribution of the magnetic force field gradient, the asymmetry between the two time scales i.e. 
${t_l}$ and 
${t_r}$, leads to the formation of sister droplets of unequal sizes. It is understandable that, due to the presence of a high force field gradient in the left branch, 
${t_l} < {t_r}$, indicating that the time required for the dispersed phase (i.e. the sister droplet) to move a characteristics length 
$(l)$ is higher in the right branch in comparison with the left branch. Moreover, we found that the flow time scale 
$({t_l})$ in the left branch is dependent on the simultaneous effect of the magnetic field flux density and the flow ratio 
$(R)$. In figure 13, we show the variation of 
${t_l}$, for the different cases under consideration. It can be clearly observed from figure 13 that the migration time scale of the sister droplet moving in the left branch is minimum for 
$B{o_m} = 7.5$ for all the cases under consideration. As a consequence, we have observed that 
$l_l^\ast $ is maximum for 
$B{o_m} = 7.5$ (refer to figures 11 and 12 for details). This observation is also suggestive of the presence of a threshold flow ratio and magnetic field intensity at which 
$l_l^\ast $ is maximum.
Figure 13. Variation of the flow time scales of the ferrofluid sister droplet migrating in the left branch for the various magnetic Bond numbers 
$(B{o_m})$ and flow ratios 
$(R)$. Here, 
${t_l}$ indicates the flow time scale in the left and right branches of the T-junction.
We explain this typical behaviour of the times scales in the left/right branch of the T-junction divergence with the help of a graphical representation, as shown in figure 14(a,b). Although, as already discussed, it is comprehensible that the magnetic field gradient promotes the asymmetric splitting of the ferrofluid droplets, it cannot be ignored that the flow of the ferrofluid droplet train in the left branch also faces resistance due to the non-uniform distribution of the magnetic field flux density (since the magnet is placed adjacent to the left branch). It is worth mentioning here that, it is because of this resistance that the droplets moving in the left branch and the right branch have different times scales, as illustrated in figure 14(a). With an initial increase in the magnetic Bond number 
$(B{o_m})$, the size of the ferrofluid droplets moving in the left branch 
$(l_f^\ast )$ increases. With an increase in 
$B{o_m}$ beyond the critical value, 
$l_f^\ast $ decreases further. This decrease in 
$l_f^\ast $ is due to the substantial increase in the restriction of the sister droplets motion in the left branch, as can be observed from figure 14(a). This argument is further justified, since we observe an increase in 
${t_l}$ beyond critical Bond number (i.e. 
$B{o_{m,Crt}} = 7.5$), as depicted in figure 13. As a consequence of this resistance, the flow rate of the continuous phase in this particular branch decreases and, primarily for this reason, we observe a decrease in the size of the sister droplet (moving in the left branch) beyond the critical Bond number 
$(B{o_{m,Crt}})$.
Figure 14. Schematic representation of the droplet breakup dynamics for the case when (a) R is constant, 
$B{o_m}$ increases, (b) 
$B{o_m}$ is constant, R increases.
Similar observation can be seen with a change in the flow ratio 
$(R)$. Note that a low flow ratio 
$(R)$ ensures a high initial slug length 
$({l_0})$ of the mother droplet and vice versa. As discussed previously, an increase in the mother droplet length will ensure a simultaneous high resistance to the flow (in the left branch), thereby further decreasing 
$l_f^\ast $ as depicted in figure 12. On the other hand, with a reduction in initial droplet length 
$({l_0})$, 
$l_f^\ast $ increases. However, a further reduction in 
${l_0}$ (i.e. increasing flow ratio, 
$R$) beyond the critical value
$({l_{0,crt}})$, limits the amount of magnetization realized by the droplet that is being split. This limited realization of the magnetic force by the ferrofluid droplets restricts the asymmetricities involved in the droplet splitting phenomenon and leads to a decrease in 
$l_f^\ast $ further. This phenomenon can be clearly observed from the graphical representation of 14(b).
4. Conclusion
In summary, we have systematically investigated the ferrofluid droplet breakup dynamics in a T-junction divergence of an LOC device in the presence of a non-uniform magnetic field. The study is especially limited to the ‘breakup with permanent obstruction’ regimes. Firstly, we have methodically explored the droplet breakup behaviour under the modulation of a non-uniform magnetic field. With the help of numerical simulations, we have investigated, the internal hydrodynamics of the droplet under the influence of a non-uniform forcing. We found that the presence of a ‘hump-like structure’ developed inside the left moving bulge (of the ferrofluid droplet) triggers the onset of augmented flow convections inside the moving volume. By performing further investigation, we have explored the roles of the critical magnetic Bond number 
$(B{o_m})$, and the flow ratio 
$(R)$ on the underlying droplet splitting phenomenon. A critical magnetic Bond number 
$(B{o_{m,crt}})$ is found for which the length of the split slug 
$(l_f^\ast )$ becomes maximum in the left branch, while any increase in magnetic field strength beyond the critical value is seen to reduce 
$l_f^\ast $. Similarly, we have also observed the presence of a threshold flow ratio 
$({R_{crt}})$ at which the length of the split droplet 
$(l_f^\ast )$ is a maximum. To understand these characteristic behaviours of the ferrofluid droplet undergoing deformation in the presence of non-uniform force field, we have investigated the role played by the various dominant time scales. We found that, due to the imbalance between the flow time scale in the left and right branches of the T-junction divergence, asymmetric droplet breakup takes place. The non-uniformity of the magnetic field (placed adjacent to the left branch) initiates a lower flow time scale in the left branch in comparison with the right branch. Also, we unveiled that, at the critical state, the flow time scale is lowest in the left branch, thereby ensuring the largest size of the sister droplet formed therein. At the critical 
$B{o_{m,crt}}$, the length of the sister droplet (in the left branch) is increased by almost 33 % when compared with the base case, i.e. in the absence of a magnetic field. The proposed technique is expected to open up new avenues towards the development of a LOC device for achieving rapid and controlled droplet splitting, typically finding relevance in point-of-care-related diagnostics.
Supplementary material
Supplementary material is available at https://doi.org/10.1017/jfm.2022.504.
Acknowledgements
The authors acknowledge the help received from Mr D. Wankawala, a research scholar in the Microfluidics and Microscale Transport Processes Laboratory, Department of Mechanical Engineering, regarding characterization of the particle size distribution in the ferrofluid solution.
Funding
S.S. and P.K.M. gratefully acknowledge the financial grant obtained from NEWGEN IEDC. Also, P.K.M. gratefully to acknowledges the financial support provided by the DSIR, Govt. of India, through project no. DSIR/PRISM/170/2020-21 and SERB (DST), India, through Project No. MTR/2020/000034. The authors also acknowledge the CIF, IIT Guwahati for the support in characterization of ferrofluid.
Declaration of interests
The authors report no conflict of interest.
