Stochastic Density Functional Theory on Lane Formation in Electric-Field-Driven Ionic Mixtures: Flow-Kernel-Based Formulation.

Simulation and experimental studies have demonstrated non-equilibrium ordering in driven colloidal suspensions: with increasing driving force, a uniform colloidal mixture transforms into a locally demixed state characterized by the lane formation or the emergence of strongly anisotropic stripe-like domains. Theoretically, we have found that a linear stability analysis of density dynamics can explain the non-equilibrium ordering by adding a non-trivial advection term. This advection arises from fluctuating flows due to non-Coulombic interactions associated with oppositely driven migrations. Recent studies based on the dynamical density functional theory (DFT) without multiplicative noise have introduced the flow kernel for providing a general description of the fluctuating velocity. Here, we assess and extend the above deterministic DFT by treating electric-field-driven binary ionic mixtures as the primitive model. First, we develop the stochastic DFT with multiplicative noise for the laning phenomena. The stochastic DFT considering the fluctuating flows allows us to determine correlation functions in a steady state. In particular, asymptotic analysis on the stationary charge-charge correlation function reveals that the above dispersion relation for linear stability analysis is equivalent to the pole equation for determining the oscillatory wavelength of charge-charge correlations. Next, the appearance of stripe-like domains is demonstrated not only by using the pole equation but also by performing the 2D inverse Fourier transform of the charge-charge correlation function without the premise of anisotropic homogeneity in the electric field direction.

1. Introduction

Many industrial processes involve the transport of colloidal particles under external fields. For example, particles are driven by stirring or shearing, whereas other typical driving forces arise from gravity and an external electric field. The response of particles to the driving forces is a subject of considerable practical interest. In particular, electric-field-driven particles, with which we are concerned in this paper, play significant roles in biological ion channels, micro/nanofluidic devices for environmental and biomedical applications, and electrolyte-immersed porous electrodes for electrochemical applications [1,2,3]. We consider binary mixtures of symmetric ions under external electric fields as the electric-field-driven particles, which include not only oppositely charged colloidal mixtures but also electrolytes and room-temperature ionic liquids. Recently, the binary ionic mixtures under external fields are increasingly attracting much attention due to their diverse applications not only in chemistry and biology [1] but also in renewable energy devices such as batteries, supercapacitors, and separation media [3].

The driven binary mixtures, in which two populations of particles are driven in opposite directions, undergo an out-of-equilibrium transition [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18]. For example, colloidal particles driven by a strong external field self-organize into strongly anisotropic stripe-like or layered structures (i.e., lanes), thereby representing a prototype of nonequilibrium phase transition [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19]. The underlying mechanism of lane formation has been ascribed to the competition between thermodynamic tendency to mix binary colloids and kinetic preference to segregate colloids of the same kind for reducing collisions due to oppositely driven migrations [4,5,6,7,8,9,10]. Experimental and simulation studies have demonstrated that driven colloids provide a testbed for pattern formation occurring in many nonequilibrium systems; for the driven binary mixtures cover a variety of driven systems ranging from colloidal suspensions and electrolytes to active matter composed of autonomously moving agents such as pedestrians [4,5].

Experimental and simulation studies have also observed the band formation of like-charged colloidal particles, other than laning [11,12,13,14,15,16]. The bands are aligned in a direction non-parallel to the applied field direction and lead to a jammed state where the particles block each other’s motion. We have thus obtained dynamic phase diagrams of steady states of laned, jammed, and mixed structures formed by oppositely charged particles under a DC or AC electric field [5,6,7,8,9,10,11,12,13,14,15,16]. Especially under oscillatory electric field, laning generally occurs for a high enough field strength and a low oscillation frequency, whereas jammed and other non-laned structures emerge depending on the magnitude of the driving field and its oscillatory frequency [11,16]. The key observables to detect the emergence of such various structures in steady states are correlation functions; however, there are a few studies of electric-field-driven systems based on the correlation function analysis [20,21,22], and no attempts to address lane formation have been made. To perform the correlation function analysis, the stochastic density functional theory (DFT) [22,23,24,25,26,27,28,29,30,31,32,33] for lane formation needs to be developed.

Let us then provide a brief review on the dynamical density functional theory (DFT) developed to describe the overdamped dynamics of Brownian particles [23]. The dynamical DFT can treat a background flow by adding an advection term; however, an additional contribution needs to be included based on phenomenological arguments [9,10,17,18,19] for explaining the lane formation irrespective of whether the density functional equation used is deterministic or stochastic. The phenomenological term added to the dynamical density functional equation considers fluctuating motions around an external-field-driven migration of a particle when neglecting other particles [9,10,17,18,19]. There are two ways to add the fluctuating flow in an advection term. One method developed for sheared colloidal suspensions adds a particle-induced fluctuation flow to the velocity field [17,18,19]. A flow kernel [17,18,19,34,35,36,37,38] introduced in this approach allows us to treat non-local effects due to density fluctuations. The other treatment provides fluctuating currents transverse to the electric field for explaining the lane formation in oppositely charged colloids under external electric fields [9,10]. Both contributions capture the coupling between flow and interparticle interactions.

The above modifications for the description of lane formation belong to the deterministic DFT. An alternative approach to the dynamical DFT adopts the density functional equation with multiplicative noise, the so-called Dean–Kawasaki (DK) equation [23]. Recently, the stochastic DFT has found the usefulness of the DK equation linearized around a reference density [22,23,24,25,26,27,28,29,30,31,32,33]. The linearized DK equation allows us to compute correlation functions for density and charge fluctuations around uniform states. It is found from the correlation function analysis that density-density and charge-charge correlations are long-range correlated in the steady state [20,27,28]. The asymptotic decay of the stationary correlation functions exhibits a power-law behavior with a dipolar character, which gives rise to a long-range fluctuation-induced force acting on uncharged confining plates [27,28].

This paper aims to develop the stochastic DFT for explaining the lane formation of binary ionic mixtures. The electric-field-driven ionic mixtures are treated as the primitive model [1], implying that we consider non-equilibrium phenomena of either symmetric electrolytes or ionic liquids, rather than colloidal mixtures. The remainder of this paper is organized as follows. Section 2 provides the basic formalism for electric-field-driven binary ionic mixtures based on the deterministic and stochastic DFTs. In Section 3, we describe the purposes of reformulating and extending the previous formulation for lane formation. Furthermore, the remaining sections serve the two purposes. In Section 4, we revisit the linear stability analysis for lane formation based on the deterministic DFT [9,10,17,18], thereby unifying previous formulations of the above additional contributions due to fluctuating flows (the first purpose). Section 5 provides the stochastic DFT for lane formation. We obtain the Fourier transforms of correlation functions from the stochastic equations, thereby demonstrating that the charge-charge correlation function verifies the stability of lane structure both analytically and numerically (the second purpose). Section 6 presents a summary and conclusions.

2. Basic Formalism

2.1. Primitive Model

We consider a binary ionic mixture of cations and anions which have equal size and equal but opposite charge using the primitive model [1]. In this model, the z-valent ions in symmetric mixtures are modeled by equisized charged hard spheres of diameter immersed in a structureless and uniform dielectric medium with dielectric constant at a temperature T. The charged spheres interact via pairwise potential () where , , and denote cation–cation, cation–anion, and anion–anion interaction potentials at a separation of , respectively. Figure 1 presents a schematic of the 2D primitive model in Cartesian coordinates, illustrating that the electrophoretic force is exerted on a single ion due to normalized electric field applied in the x-direction parallel to the unit vector .

Figure 1

A schematic of the 2D primitive model of binary ionic mixture with a static electric field applied in the x-direction. The z-valent cations and anions are modeled by equisized charged hard spheres of diameter immersed in a dielectric medium with dielectric constant at a temperature T.

It is noted that this paper defines all of the energetic quantities, including the above interaction potentials, in units of . Correspondingly, the unit of , an energetic measure of electric field strength, is , and the interaction potential is represented by using the Bjerrum length , the length at which the bare Coulomb interaction between two monovalent ions is exactly .

The dynamical density functional theories focus on instantaneous concentrations of cations () and anions () which vary depending on a time t as well as a position . Here we also use the density vector defined by

While represents the number density of ions and is equal to in average, corresponds to the charge density whose average vanishes.

2.2. Stochastic Dft: Compact Matrix Forms

The stochastic DFT is based on the formulation that adds multiplicative noise term to the deterministic density functional Equations (A1) previously used [22,23,24,25,26,27,28,29,30,31,32,33]. Therefore, the basic formalism of the stochastic DFT presented below inherits the formulation of the deterministic DFT given in Appendix A. In the stochastic DFT, the conservation equation reads where the deterministic current is expressed by Equation (A2) using the direct correlation function (DCF) between the l-th and m-th ions, and the stochastic density current is expressed as using uncorrelated Gaussian noise fields characterized by with the subscript “” representing the Gaussian noise averaging in space and time.

We provide a compact matrix form of the stochastic equation with respect to through three steps as follows: (i) we obtain stochastic currents of and from Equations (A2), (A3) and (4), (ii) we write down a matrix form without external field, and (iii) we have the target equation of -dynamics by adding the advection terms formulated from Equations (A6) and (A7).

First, we consider the relations (A4) and (A5) for the DCF, thereby transforming the sum of Equations (4) and (A2) with Equation (A3) to the linearized current as follows: where denotes the short-range part of the DCF (see also Equation (A4)), with ( or 2), and satisfies the same statistics as the relation (5) for .

Second, we rewrite Equation (3) into a compact matrix form, using and

The matrix determines restoring forces in the absence of external field ().

Third, we suppose that a fluctuating part of appears only in the y-direction, according to the previous treatments [9,10,17,18]. It follows from Equation (A6) that yielding

It is also noted that the Fourier transform of reads where we have used in the second line that the flow kernel is an odd function satisfying where we have used in the second line that the flow kernel is an odd function satisfying [17,18]. Combining Equations (8) to (14), we find the advected form of stochastic equation for under external electric field as follows: where the matrix is given by adding the advection matrix, with a new parameter being introduced (see Equations (A11)–(A16) for details). Equation (15) including the above advection term is expected to form the basis of the stochastic DFT that is capable of addressing the lane formation.

3. Our Aim

In what follows, we assess and extend the deterministic DFT on lane formation in binary ionic mixtures with external electric fields applied. There are three reasons for revisiting the deterministic approaches. First, it is necessary to clarify the consistency between previous formulations. We have introduced the flow kernel in Equation (A7) to describe a fluctuating velocity; however, this treatment has applied to one-component systems, sheared colloidal suspensions [17,18]. Therefore, the flow-kernel-based formulation for electric-field-driven mixtures needs to be developed. Second, the linear stability analysis [9,10,17,18] based on the dispersion relation (see Appendix B) has supposed that wavenumber is a real value, though wavenumber is a complex value in general [39,40,41,42]. It remains to be investigated to make the linear stability analysis while considering the imaginary part of wavenumber that determines the decay length of spatial modulation. The last reason is that the stochastic DFT, an extension of the deterministic DFT, makes it possible to investigate the mechanism of lane formation from density-density and charge-charge correlations.

Thus, we aim to obtain the correlation functions from the stochastic DFT for achieving the following two purposes.

(i) Relationship between the deterministic and stochastic DFTs—The first purpose is to understand the above dispersion relations in terms of the charge-charge correlation function. We will reveal the connection between those used in the deterministic DFT and the pole equation to find the oscillatory wavelength of charge-charge correlations.

(ii) On the uniformity of lanes in terms of correlation function analysis based on the stochastic DFT—The second purpose is to validate the approximation necessary to explain lane formation using the linear stability analysis [9,10,17,18]. Previous studies on lane phenomena have neglected the charge density modulation of the x-direction. In other words, has been assumed when investigating the charge density modulation in the transverse direction to the electric field. The correlation function analysis allows for the investigation of oscillatory decay behaviors such as the oscillatory wavelength () along the electric field. Figure 2 is a schematic of this, which illustrates the emergence of the charge modulation, or the oscillatory charge-charge correlation due to the lane formation. We perform both the asymptotic analysis of charge–charge correlations for point charges (i.e., the primitive model at ) and the 2D inverse Fourier transform of the stationary charge–charge correlation function for charged hard spheres (i.e., the primitive model for ). While the asymptotic analysis will prove that diverges at the stabilization condition of lane structure, the real-space representation of the charge–charge correlation function will clarify the decay of oscillatory correlations in lane structures as a result of the inverse Fourier transform.

Figure 2

A schematic of lane formation in a binary ionic mixture. The green and orange lanes represent aligned segregation bands of cations and anions, respectively. Correspondingly, the positive and negative signs seen on the lanes indicate that each lane is a mesoscopically charged object. The wavelengths, and , in x-and y-directions are related to wavenumbers as and (i.e., Equation (A17)), respectively. In this study, these wavenumbers are determined by Equations (34) and (37) when considering point charges.

4. Correlation Functions Determined by the Stochastic Dft

4.1. Stationary Condition of Correlation Functions

The stochastic formulation allows us to provide the Fourier transforms of correlation functions for and at equal times [20,22,26,27,28,29]. These correlation functions are defined using as

The compact form (8) of the stochastic equation for is solved to obtain [20,25,26,27,28,29] where Equations (5) and (9) provide

Plugging Equations (19) and (20) into the definition (18), we have where it follows from the relation (20) that

It has been shown that the stationary condition for the expression (21) reads [20,22,26,27,28,29]

The four matrix elements of , or the four kinds of correlation functions in Equation (18), can be determined by four simultaneous equations generated from the above stationary condition (21) (see Appendix C for details).

4.2. Obtained Forms of Stationary Correlation Functions

As derived in Appendix C, Equation (23) yields the density–density and charge–charge correlation functions at equal times, and , as follows: where using . In what follows, two limiting cases are considered for and : (i) we confirm that these converge to the equilibrium correlation functions of electrolytes at and , and (ii) we see the dispersion relation given by Equation (A15) in terms of the stationary correlation functions at , according to the approximation (A9).

Before proceeding, we need to connect the long-range part of the DCF and the Coulomb potential to provide the Poisson-like equation and the Debye–Hückel screening length. In general, is expressed as using the weight function . For instance, with is a well-known form of the 3D primitive model [26,39,40]. Thus, the Poisson equation is generalized to the finite-spread type [43,44]: when defining the Coulomb potential as

It follows that where with denoting the conventional Debye–Hückel screening length.

First, we consider equilibrium electrolytes. Since we have and at , Equation (24) is reduced to

We also have , , and at . Hence Equation (31) reads thereby confirming that the charge–charge correlation function satisfies not only the electroneutrality but also the Stillinger–Lovett second-moment condition [27,28,29].

5. Lane Formation in Terms of Charge–Charge Correlation Function

5.1. Asymptotic Behavior of Charge–Charge Correlations

The lane formation has been investigated for , according to previous studies on lane formation. This implies that the density modulation along the electric field (or the x–direction) at a given y–coordinate is negligible. Equation (24) at transforms to with being replaced by in Equation (25). Focusing on the pole equation (i.e., ) for , we have

It is noted that the above Equation (34) is identified with the key equation previously used for determining the mean wavelength of lanes, which has been referred to as the dispersion relation (A15) at in terms of the linear stability analysis [9,10,17,18] (see Appendix B).

The advantage over the linear stability analysis is that the pole Equation (34) provides the long-range behavior of density profile, or the decay length and oscillatory wavelength in the asymptotic decay of charge–charge correlation function [39,40,41,42]. In particular, the correlation function analysis ensures the stability of steady-state lane structure only when a solution to Equation (34) has a purely real wavenumber , which is the case with the above linear stability analysis.

Furthermore, the asymptotic decay analysis of charge–charge correlations allows us to validate the supposition of uniformity along the electric field (or ) in the anisotropic lane structure as given in Figure 2. To assess the validity of , we evaluate the inverse Fourier transform of stationary charge–charge correlation function as follows:

We evaluate this inverse Fourier transform using the approximate form of for and . The denominator given in Equation (24) is approximated by

Therefore, the pole equation for yields providing the wavelength when is a purely real value for , (see Equation (A16)), and .

The asymptotic analysis allows us to provide the long-range oscillatory behavior of as follows:

It should be noted that the pole Equation (34) for reads in the present case. We obtain from Equation (39)

Combining Equations (37) and (39), we also find

Namely, we have when forming the lane structure with its period of . Thus, it is verified analytically that each lane is uniform along the electric field the present approximation (A9) as far as point charges () are considered.

The expression (40) of , or the lane width, reveals the underlying physics of lane formation. Each lane has the energetic cost of Coulomb repulsions due to clustering of either cations or anions, which explains why lanes can be wider as is shorter and the screening of Coulomb interactions is stronger. Despite this energetic cost, the lane formation is favored because collisions due to oppositely driven migrations are reduced by segregation of cations or anions. The kinetic preference is enhanced by increasing the strength of external field; accordingly, Equation (40) implies that the lane width is larger with increase of .

5.2. Charge–Charge Correlations on 2D Cross Section of the 3D Primitive Model

The preceding subsection has analytically demonstrated that the dispersion relation based on the conventional linear stability analysis [9,10,17,18] is equivalent to the asymptotic decay analysis of the charge–charge correlation function. We have also verified that the dispersion relation applies to the emergence of a lane structure for point charges. Turning our attention to charged hard spheres of finite size, however, it remains to be validated whether we can neglect the decay of charge–charge correlations. At least, for the 3D primitive model in the absence of an electric field, theoretical and simulation studies have found that oscillatory decay of charge-charge correlations has been observed beyond the Kirkwood crossover condition [39,40,41,42]. In terms of the asymptotic decay analysis, the solution to the pole equation becomes complex at the Kirkwood crossover when considering the wavenumber-dependence of , and the imaginary part of the solution represents the finite decay length of charge-charge correlations.

To investigated the oscillatory decay behavior in the presence of an external field, we examine the stationary charge–charge correlation function concerning a 2D cross-section of the 3D primitive model. Figure 3 represents a schematic of the present 3D system. From Figure 3, we can see that the plane in Figure 1 corresponds to the cross-section formed by the x- and y-axes embedded in this 3D system. The advantage of considering the 3D primitive model is that we can use the analytical form of DCF: the long-range part is given by Equation (26) with whereas the short-range part reads in the modified mean spherical approximation [45]. Upon introducing the 3D volume fraction , Equations (26) and (44) transform Equation (25) to and

Figure 3

A schematic of the 3D primitive model in Cartesian coordinates illustrates a binary ionic mixture confined between two parallel plates. While the y-axis is perpendicular to these plates, the electric field is applied in the x-axis.

Therefore, under the simplification of , the expressions (24), (25), (45) and (46) for imply that the inverse Fourier transform of depend on the three control parameters: , , and .

Let the 3D wavevector be which has –component in addition to –and –components. However, we set , which leads to the consideration of charge-charge correlations averaged over the z–axis density distribution in Figure 3 [22]. Accordingly, we can perform the 2D inverse Fourier transform of in the Cartesian coordinates similar to those given in Figure 1. Figure 4 shows some of the results. Figure 4a,b on the left side are the results considering the presence of with . Meanwhile, Figure 4c,d on the right side ignore the fluctuating part with . The other parameters are common to the results on the left and right sides. Namely, for all of the results in Figure 4, and in Figure 4a,c, and the concentration is increased by 10% in Figure 4b,d: the screening effect of Coulomb interactions is enhanced to due to in Figure 4b,d.

Figure 4

The real-space representation of the charge-charge correlation function at is shown using heat maps where the length scale is in units of diameter . We obtain the real-space correlation function from performing the 2D inverse Fourier transform of given by Equations (24) and (25). The remaining parameter set of is (a) , (b) , (c) , and (d) .

We can draw three conclusions from comparing the results in Figure 4.

First, Figure 4a,b verify the lane formation of binary ionic mixtures in terms of charge-charge correlations. Especially in Figure 4b, we observe no decay of correlations in the electric-field direction over the length scale of 10 times the diameter of charged hard spheres. The oscillatory charge-charge correlations demonstrate that each lane of the 3D primitive model can be homogeneous in the electric-field direction, which agrees with the analytical investigations in Section 5.1.

Second, comparison between Figure 4a,b suggests the underscreening behavior [39,40,41,42]. On the one hand, Figure 4a indicates that the stationary charge–charge correlation function converges to zero far from the origin of , thereby illustrating an oscillatory decay behavior. In Figure 4b, on the other hand, we observe little change in the heat map color along the electric field direction. In other words, the purely oscillatory behavior, which is the premise of the linear stability analysis previously made, is demonstrated in Figure 4b. This change from Figure 4a to Figure 4b suggests that the decay length is longer as the volume fraction , or the ion concentration, increases similarly to the underscreening behavior in binary ionic mixtures with no electric field applied above the Kirkwood crossover where the equilibrium charge–charge correlation function exhibits oscillatory decays [39,40,41,42]. We have confirmed such underscreening behavior with an electric field applied.

Third, the difference between the results in Figure 4 on the left and right sides reveals that anisotropic oscillatory correlations, which reflect the lane formation, disappear in Figure 4c,d because of the absence of the fluctuating flow, , given by Equations (13) and (14). It is also important to note that the scale of the color bar on the right side is times the scale on the left side. In other words, is almost zero in Figure 4c,d. The weak charge–charge correlations imply that electric-field-driven binary ionic mixtures are uniform in the absence of the fluctuating flow which arises from collisions due to oppositely driven migrations of cations and anions.

6. Summary and Conclusions

The charge–charge correlation function studied so far can be detected using X-ray and/or neutron scattering experiments [46]. We would therefore like to evaluate experimental conditions that are consistent with the numerical results in Figure 4. For example, we consider as an room-temperature ionic diluted with propylene carbonate. It follows that in correspondence with the supposition that in Figure 4. Also, the parameters, and (or ), used in Figure 4a read V/m and 0.31 M, respectively, for the room-temperature ionic liquid. These are plausible values according to previous simulation and experimental studies [47,48]; in particular, it is interesting to note that molecular dynamics simulations of room-temperature ionic liquids have revealed that V/m corresponds to a boundary value beyond which the ionic liquids are reorganized into nematic-like order and exhibit anisotropic dynamics [49].

Finally, we summarize the results presented so far, according to the two purposes mentioned in Section 3, the section of our aim.

(i) Relationship between the deterministic and stochastic DFTs—The wavenumber appearing in the dispersion relation (A12) can actually be a complex number. It is appropriate for understanding the underlying physics of the complex wavenumber to see correlation functions instead of the dispersion relation. Hence, we have addressed the first purpose using the stochastic DFT for lane formation, thereby allowing us to obtain density–density and charge–charge correlation functions in a steady state. We have demonstrated that the asymptotic analysis of the charge–charge correlation function is equivalent to the linear stability analysis based on the dispersion relation (see also Appendix B). Specifically, the pole equation used in the asymptotic analysis proved equivalent to the lane stability condition obtained from the dispersion relation. The analytical framework is thus available to find the presence or absence of decay length and oscillatory wavelength in the oscillatory decay of the correlation function. In other words, it became possible to examine the spatial stability of the lane formation more precisely.

(ii) On the uniformity of lanes in terms of correlation function analysis based on the stochastic DFT—We have obtained the Fourier transform of the stationary charge–charge correlation function . Nevertheless, the previous treatments [9,10,17,18] have supposed that in advance prior to the inverse Fourier transform. The pole equation obtained from the correlation function at is an equation in which only is a variable, and it is equivalent to the linear stability condition determined from the dispersion relation, as described in the first purpose. Namely, in the previous treatments [9,10,17,18] described above, the presence or absence of lane formation is examined on the premise that the lane formation is uniform in the electric field direction. It is necessary to show the uniformity in the electric field direction itself without assuming the uniformity in the electric field direction. Therefore, we evaluated the Fourier transform of from the pole equation for point charge systems (i.e., ) where we have and that . In other words, we investigated the stability of lane formation with only Coulomb interaction at , showing that the oscillatory wavelength diverges at , or the solution to the pole equation given by Equations (34) or (39). Thus, the approximation has been validated analytically. Figure 4 also demonstrates numerically that, above the Kirkwood crossover [39,40,41,42], the oscillatory decay length observed for the 3D primitive model (i.e., ) is longer with the increase of ion density ; the underscreening behavior under external field applied remains to be investigated in more detail (see also [22]).