PHYSICAL REVIEW E

VOLUME 55, NUMBER 2

FEBRUARY 1997

Two-order-parameter model for an oil-water-surfactant system
Shigeyuki Komura*
Department of Mechanical System Engineering, Kyushu Institute of Technology, Iizuka 820, Japan

Hiroya Kodama†
Division of Science of Materials, Graduate School of Science and Technology, Kobe University, Kobe 657, Japan~Received 26 June 1996!
Dynamics of microphase separation of an oil-water-surfactant system is investigated by means of cell
dynamical system approach on the basis of the proposed two-order-parameter time dependent GinzburgLandau model. For equal volumes of oil and water, the time evolution of characteristic length scale of domains
is investigated by changing the average surfactant concentration. The more the amount of surfactant, the slower
the dynamics. The results are analyzed by using the crossover scaling assumption. @S1063-651X~97!07702-7#
PACS number~s!: 82.70.Kj, 64.75.1g, 82.20.Wt

I. INTRODUCTION

Microemulsions which are a mixture of oil, water, andsurfactant are known to exhibit various interesting microstructures depending on the temperature or the composition@1#. When the concentration of surfactant is relatively large,
they show a rich variety of regularly ordered structures suchas the cubic phase, the hexagonal phase, and the lamellar
phase. By lowering the concentration of surfactant and if thevolumes of oil and water are not very different, microemulsions form a bicontinuous structure where a multiply connected randomly oriented monolayer of surfactants separates
oil-rich and water-rich subvolumes. There are several experimental methods for studying the convoluted structure of oil
and water, such as by measuring the conductivity @2# or the
diffusivities of molecules @3#. A direct observation of therandomly intertwined structure by using freeze-fracture microscopy has also been reported @4#. Several x-ray and neutron scattering experiments showed that the structure factorof the bicontinuous phase shows a peak at nonzero wave
vector k, indicating that there is a structure on a length scale
2 p /k @5#. It is also found that for microemulsions containing
equal volumes of oil and water, the peak position shifts tolarger values of k and peak height diminishes as the surfactant concentration is increased @6#.
When one quenches a ternary system from a high temperature homogeneous phase where the system is uniformly
mixed to a low temperature phase where a certain structureappears, the average domain size increases in time until itreaches the equilibrium size. Since the surfactants act tolower the interfacial tension and the driving force for thephase separation is decreased, it is natural to expect that the
dynamics of domain growth becomes significantly slower inthe presence of surfactants. With this perspective, several
people have investigated the effect of added surfactants on
the phase separation dynamics using different models. For
instance, Kawasaki and Kawakatsu have proposed a ‘‘hybrid
model’’ where oil and water are represented by coarse-*Electronic address: komura@iizuka.isc.kyutech.ac.jp†

Electronic address: hkodama@icluna.kobe-u.ac.jp1063-651X/97/55~2!/1722~6!/$10.00

55

grained fields and surfactants are treated microscopically @7#,
whereas Laradji et al. performed molecular dynamics simulations @8#. Both of these works have shown that the system
containing surfactants exhibits a slow nonalgebraic growthof the domains, in contrast to the pure binary systems.
Qualitatively similar results have also been reported in theother paper by Laradji and his co-workers, who proposed aphenomenological two-order-parameter Ginzburg-Landau
free energy associated with standard time dependent
Ginzburg-Landau ~TDGL! equations @9#. In their model, oneof the order parameters represents the local concentrationdifference between oil and water, while the other one represents the local surfactant concentration. Recently, Pätzold
and Dawson extended the model by Laradji et al. to incorporate the hydrodynamic effects by coupling the TDGL
equations to Navier-Stokes-type equations @10#. They foundthat, in the presence of the hydrodynamic interactions, thecrossover scaling exponent becomes different from that obtained by Laradji et al. @9#.
However, the book by Gompper and Schick criticizes thetwo-order-parameter model by Laradji et al. as not being
well defined, since the free energy of configurations withlarge surfactant concentration at the oil-water interfaces is
not bounded from below @1#. So far, a quite general expression of the two-order-parameter Ginzburg-Landau model has
also been given by Gompper and Schick @1#. In fact, themodel by Laradji et al. can be considered as one of the special cases of their expression with too much simplification.
For the above reasons, a quantitative study of the phase separation dynamics based on the well defined two-orderparameter model has not yet been done.
In this paper we shall propose a minimum two-orderparameter Ginzburg-Landau model in which the above problems are removed but which still expresses the essential features of the ternary system. Our model intrinsically includes
the preferred value of the surfactant concentration when surfactant molecules aggregate. Moreover, we required in ourmodel that when surfactants aggregate at the interface, theinterfacial tension becomes very small. Using our model, we
numerically study the dynamics of phase separation. Here we
pay attention to the case where the volumes of oil and water
are equal.
1722© 1997 The American Physical Society

55

TWO-ORDER-PARAMETER MODEL FOR AN OIL-WATER- . . .

This paper is constructed as follows. In the next section,
after reviewing the model by Laradji et al. we present ourmodel and discuss its physical meanings. In Sec. III we describe our simulation method which is essentially based on
the cell dynamical system approach. The analysis of oursimulation results and discussion are given in the following
two sections.
II. MODELIn this section we shall first review the two-orderparameter model proposed by Laradji et al. @9#. In theirmodel, the local concentration difference between oil andwater is described by c (r) and the local surfactant concentration measured from a certain reference value is denoted asr (r). Their model isF5E

dr@ d ~ “ c ! 2 2a c 2 1u c 4 1b r 2 1g c 2 r 2 2s r ~ “ c ! 2 # ,~1!

where d, a, u, b, g, and s are positive constants. The
term b r 2 prevents the surfactants from forming clusters. The
local coupling term g c 2 r 2 guarantees that the local surfactant concentration remains small in the bulk where u c u islarge. The last nonlocal coupling term 2s r (“ c ) 2 favors thesurfactants to sit at oil-water interfaces. As mentioned in the
introduction, however, this free energy functional is notbounded from below for positive values of s @1#. AlthoughPätzold and Dawson reported that they did not find any numerical instabilities within the parameter values they used
@10#, we observed a diverging tendency of the surfactant concentration at the oil-water interfaces as we decreased thesimulation mesh size. The finite value of r is supported onlyby the nonzero simulation mesh size. This is obviously notreasonable from the point of view of solving continuous partial differential equations. Moreover, we found that the domains do not flow globally in the presence of a convective
macroscopic flow within their model. This problem also
seems to be related to the model intrinsic singular behavior
at the oil-water interfaces.
Here we propose a different two-order-parameter
Ginzburg-Landau free energy functional which has none ofthe drawbacks mentioned above. What we have required in
our model are that ~i! the profiles of c and r at oil-waterinterfaces do not depend on the average values of c and r~denoted hereafter as c̄ and r̄ , respectively! and that ~ii! thecoarse-graining dynamics of c based on the free energy becomes slow when the amplitude of r at the interface takes acertain saturated value. The proposed minimum model which
fulfills these requirements is
F5

Estates, i.e., r 50 and r 5 r s . The physical meaning of these
two states is as follows. The state r 50 corresponds to thecase in which the system is locally occupied either by oil orwater. There are no surfactants in the considered local volume. The state r 5 r s corresponds to the case in which thelocal volume is occupied only by surfactants. The quantity
r s can be considered to represent the density of condensedhydrocarbon chains of surfactants when they self-assemble.
It is likely that in any type of surfactant aggregates, the density of hydrocarbon chains does not change appreciably. Itshould be noted here that we did not include any gradientterm of r , say (“ r ) 2 . This is physically reasonable since theenergy cost due to the direct attachment between hydrocarbon chains and oil molecules or between hydrophilic head
and water molecules is small. As regards the coupling terms,
the local coupling term g c 2 r 2 in Eq. ~1! has been dropped
out, whereas only the nonlocal coupling term 2s r (“ c ) 2 isleft here. Due to the latter term, the state with r 5 r s tends to
occupy the narrow region around the oil-water interfaces. It
is also essential in microemulsions that the interfacial tension
vanish when the interface is saturated with surfactants@11,12#. For this reason, we chose the values of parameters tosatisfy the relation d5s r s , which plays an essential role forour requirement ~ii! in the model. We believe that Eq. ~2! isone of the minimum models which is sufficient to grasp theessential properties of microemulsions. Another possibilityof the model will be discussed in the last section.
For the time evolution of c (r,t) and r (r,t), we assume
the standard TDGL equations as in Ref. @9#. Since both cand r are conserved quantities, TDGL equations are]cdF
5Mc ¹ 2 1 h c ~ r,t ! ,
]t
dc

~3!

]r
dF
5Mr ¹ 2 1 h r ~ r,t ! .
]t
dr

~4!

Here Mc and Mr are transport coefficients, h c and h r represent thermal noise which satisfies the fluctuation-dissipation
theorem^ h c ~ r !~ r,t ! h c ~ r !~ r8 ,t 8 ! & 522k B TM c ~ r ! ¹ 2 d ~ r2r8 !
3 d ~ t2t 8 ! ,~2!

where w, d, a, u, e, r s , and s are positive constants.
First we have added the term w(¹ 2 c ) 2 which prevents themodel from becoming unbounded. Next we replaced the local potential of r with a double-minimum potential
e r 2 ( r 2 r s ) 2 which allows the coexistence of the two bulk

~5!

where k B is the Boltzmann constant and T is the temperature.
By inserting Eq. ~2! into Eqs. ~3! and ~4!, the time evolution
equations can be written explicitly as]c
5Mc ¹ 2 @ 22a c 14u c 3 12w¹ 4 c 22 ~ d2s r ! ¹ 2 c]t
12s ~ “ c ! • ~ “ r !# 1 h c ~ r,t ! ,dr@ w ~ ¹ 2 c ! 2 1d ~ “ c ! 2 2a c 2 1u c 4

1e r 2 ~ r 2 r s ! 2 2s r ~ “ c ! 2 # ,

1723

~6!

]r
5Mr ¹ 2 @ 2e r ~ r 2 r s !~ 2 r 2 r s ! 2s ~ “ c ! 2 # 1 h r ~ r,t ! .
]t~7!
In the present work we have entirely ignored the hydrodynamic interactions which might play an important role inmicroemulsions. One of the realistic systems which corresponds to our model is a binary homopolymer mixture con-

1724

SHIGEYUKI KOMURA AND HIROYA KODAMA

taining diblock copolymers @9,13#. Since we have not included the hydrodynamic interactions, our model lacks abare viscous time scale. Nevertheless, we can choose themodel intrinsic time scale as the inverse of the initial growthrate of the most unstable mode. The wave vector dependent
initial growth rate is defined as] c k~ t !
5L ~ k! c k~ t ! ,]tc k~ t ! 5

E c

dr ~ r,t ! exp~ ik•r! .

~9!

Since the typical time scale for surfactants to assemble at theoil-water interfaces is long compared to that of the initial
phase ordering process, we neglect here the effect of surfactants on L(k). Hence we haveL ~ k! 522Mc ~ wk 6 1dk 4 2ak 2 ! ,~10!

where k5 u ku . The most unstable mode k5k 0 which maximizes L(k) is
k 20 52d1 Ad 2 13aw
.
3w

~11!

Thus the typical time scale t 0 in our model is given byt 0 51/L ~ k 0 ! .

The first term in the right-hand side of Eq. ~15! is introduced
for the sake of numerical stability @14#. In the above equations, the discretized differential operators ¹̃ and ¹̃ 2 are defined as
˜ f 5 1 „f ~ n 11,n ! 2 f ~ n 21,n ! , f ~ n ,n 11 !
“
x
y
x
y
x
y
2

~8!

where we have introduced the spatial Fourier transform ofc (r,t) as~12!

III. CELL DYNAMICAL SYSTEM APPROACH

In order to solve the above time evolution equations, weused the cell dynamical system ~CDS! approach proposed by
Oono and Puri @14#. The CDS model is a space-time discretemodel to describe a phenomenon at the mesoscopic level andproved to be an efficient algorithm for numerical simulations. Here we restricted ourselves to a two-dimensional system. Accordingly, the space coordinate is specified by thelattice point n5(n x ,n y ) in an L3L square lattice. We also
imposed periodic boundary conditions. The CDS equations
corresponding to Eqs. ~6! and ~7! arec ~ n,t11 ! 5 c ~ n,t ! 1Mc ¹̃ 2 I~ n,t ! 1Cc h ~ n,t ! ,~13!

r ~ n,t11 ! 5 r ~ n,t ! 1Mr ¹̃ 2 J~ n,t ! 1Cr h 8 ~ n,t ! ,~14!

where I(n,t) and J(n,t) are the discrete thermodynamic
forces given by55

2 f ~ n x ,n y 21 ! …

~17!

and
¹̃ 2 f 5

12( f ~ nearest-neighbor cells!

1

1
4( f ~ next-nearest-neighbor cells! 23 f ,~18!

respectively. The noise terms in Eqs. ~13! and ~14! are givenbyh ~ 8 ! ~ n,t ! 5 h ~x8 ! ~ n x 11,n y ,t ! 2 h ~x8 ! ~ n x ,n y ,t !
1 h ~y8 ! ~ n x ,n y 11,t ! 2 h ~y8 ! ~ n x ,n y ,t ! , ~19!
where h (x8 ) and h (y8 ) are random numbers uniformly distributed in the interval @ 21,1# and C c ( r ) are the noise amplitudes taken as independent parameters in CDS @15#. The
noise perturbations are important in order to accelerate theevolution processes and prevent domains from freezing @10#.
In our simulations, we fixed the parameters to L5128,
A51.3,
W50.2,
D50.5,
S50.5,
E50.25,r s 51, Mc 5Mr 50.05, Cc 5Cr 50.02, and c̄ 50, whereasr̄ has been changed as r̄ 50.1,0.2,0.3,0.4, and 0.5. The initial distributions of c and r are specified by random uniform
distributions
intherangesand@ 20.01,0.01#@ r̄ 20.01, r̄ 10.01# , respectively, which correspond to acompletely disordered uniform state. With our choice of parameters, the most unstable mode Eq. ~11! is k 0 '0.51 andthe typical time scale Eq. ~12! is t 0 '5.03102 .
IV. RESULTS

The typical time evolutions of c and r are depicted inFig. 1 for ~a! r̄ 50.1 and ~b! r̄ 50.4. For a quantitative discussion, we have calculated the average domain size of c in
the following way. First a discrete Fourier transform ofc (n,t) is defined byc k~ t ! 5 ( c ~ n,t ! exp~ ik•n! ,n~20!

I~ n,t ! 52Atanhc 1 c 1W ~ ¹̃ 2 ! 2 c 2 ~ D2S r ! ¹̃ 2 c
˜ c !•~ “˜ r !,
1S ~ “1 ˜ 2c! .
J~ n,t ! 5E r ~ r 2 r s !~ 2 r 2 r s ! 2 S ~ “
2

~15!
~16!

with k52 p n/L and nP $ 0,1, . . . ,L21 % 2 . The structure factor is given by S(k,t)5 ^ c k(t) c 2k(t) & , where the average isover the ensemble of systems. In this paper, we calculated
the time dependent ~inverse! characteristic length scale defined by @16#

TWO-ORDER-PARAMETER MODEL FOR AN OIL-WATER- . . .

55

1725

FIG. 1. Time evolutions of cand r for ~a! r̄ 50.1 and ~b!
r̄ 50.4. The left and right columns show the spatial distributionof c and r fields, respectively.
Dark area denotes the region ofhigher values of c and r . Noticethat the typical time scale ist 0 '5.03102 .

^k~ t !&5( u ku 21 S ~ k,t !
kÞ0
(

kÞ0

u ku 22 S ~ k,t !

.

~21!

This expression provides better estimation of the characteristic length scale than that obtained by using the sphericallyaveraged structure factor @16#. In Fig. 2 ^ k(t) & is plotted as afunction of time step for different values of r̄ . Each line
corresponds to a single run. By turning off the coupling between c and r , we have also included the result of ordinary
spinodal decomposition which exhibits the well-known
t 21/3 evolution. The estimated most unstable mode k 0 is consistent with our simulation result. We recover the following
features which have also been found in the previous works.
~i! In the presence of surfactants, evolution of the pattern
becomes exceedingly slow, which appears to be almost loga-

1726

SHIGEYUKI KOMURA AND HIROYA KODAMA55

FIG. 4. Scaling plot of f (x) as a function of x. The straight line
has a slope of 1/3.
FIG. 2. Inverse characteristic length scale ^ k(t) & as a function oftime step for different values of r̄ . The straight line has a slope of
21/3.

rithmic in time for large time steps and large r̄ . ~ii! At almost
any time step, the average domain size decreases as r̄ isincreased.
In addition to these features, especially for small r̄ , we
can observe the coexistence of two types of interfaces at theearly stage of phase separation, i.e., the coexistence of saturated interface and unsaturated interface. At the saturated
interfaces, r takes a saturated value ~slightly above 1!,
whereas the value of r at the unsaturated interfaces is almost
the same as that in the bulk phase. This can be more clearlyseen in Fig. 3 where we plotted the cross section profiles of
the two fields along the straight lines drawn in Fig. 1~a!
(t550 000). As the phase separation proceeds, the system
will eventually be governed by the saturated interfaces andwill reach the equilibrium configuration. It is interesting tonote that the dynamics of phase separation in microemulsions can be understood from the point of view of motion of
one-dimensional interface ~or contour! separating saturated
and unsaturated interfaces.

The time evolution of the characteristic length scale forvarious values of r̄ has been analyzed by using the crossover
scaling assumption. Since the hydrodynamic interaction is
ignored in our simulation, we consider a similar form proposed in Ref. @9#;^ k ~ t ! & 5t 21/3 f „~ r̄ 2 r 0 ! 3 t…,~22!

where f (x) is a scaling function with f (x); const forx→0 and f (x);x 1/3 for x→`. Here, we denote the averagevalue of r in the bulk phase as r 0 , which can be considered
to correspond to the origin of the surfactant density. We
checked that this value does not depend on r̄ . Hence forsmall x, the characteristic length scale exhibits a t 21/3 behavior as in the ordinary spinodal decomposition, whereas forlarge x, it becomes proportional to the total surfactant concentration. The latter fact has been confirmed both by experiments @6# and by Monte Carlo simulations @17# and can beinterpreted in the following way; if all the interfaces are saturated, the total amount of surfactant should be proportional
to the total area ~‘‘length’’ in our simulation! of the interface
which is also proportional to the inverse characteristic
length. In Fig. 4 we have plotted the scaling function f (x) asa function of x for various values of r̄ by fixing r 0 520.1.
Although the data collapse for larger values of r̄ is satisfactory, the data for small r̄ ~such as r̄ 50.1) deviate from theuniversal behavior. We also comment that the present scaling
behavior can be observed even in the presence of noise,
which is in contrast to what has been reported previously@10#.
V. SUMMARY AND DISCUSSIONS

FIG. 3. Profiles of c and r along the straight lines drawn in Fig.
1~a! (t550 000). The solid line represents c and dashed line represents r .

In this paper, using the CDS approach, we have investigated the effect of surfactants on the dynamics of phase separation between oil and water on the basis of the proposed
minimum two-order-parameter TDGL model. We restricted
ourselves to the case where the volumes of water and oil are
equal whereas the average surfactant concentration has been

55

TWO-ORDER-PARAMETER MODEL FOR AN OIL-WATER- . . .

1727changed systematically. The time evolution of the typicallength scale of the domains is characterized by an extremely
slow dynamics. Our results can be interpreted according to
the dynamical scaling assumption.
We comment on the other possibilities of the two-orderparameter model. As regards the double-minimum local potential of r in Eq. ~2!, one may consider replacing it with asingle-well potential r 2 which appears in the original modelby Laradji et al. @see Eq. ~1!#. @The term (¹ 2 c ) 2 is always
necessary for the stability of the model.# We also examined
this case, but it turned out that the system exhibits macrophase separation rather than microphase separation within
the parameter we investigated. This means that the surfactants do not adsorb enough at the interface to suppress thephase separation dynamics as far as r (“ c ) 2 is the only included coupling term. This situation can be changed, for example, by including a local coupling term such as c 2 r according to the symmetry consideration. ~Notice that thec 2 r 2 term in the model by Laradji et al. gives a higher ordercontribution than c 2 r .! In this paper we tried to reduce thenumber of different types of coupling terms to as few as
possible. We consider that the nonlocal coupling term usedin this paper works more essentially for the slow dynamics.
Finally we discuss our choice of parameters. In principle,
d and s in Eq. ~2! are independent parameters. Although weassumed the relation d's r s according to the physical reasons in Sec. II, the following situations may be more general.
If d@s r s , the surfactant cannot cancel the interfacial tensionof the oil-water interface. In this case, the system will undergo a macrophase separation. If d!s r s , the interfacial
tension becomes essentially negative and the system willproduce more interfaces.
In fact, our primary purpose in considering the GinzburgLandau free energy for microemulsions is to investigate its
rheological behavior in ways such as by measuring the response to the applied shear flow. According to our preliminary study, we found a pronounced effect of the added surfactants on the macroscopic mechanical properties. The
details of our results will be published in the future @18#.

@1# G. Gompper and M. Schick, Self-Assembling Amphiphilic Systems ~Academic Press, London, 1994!.
@2# M. Laguës and C. Sauterey, J. Phys. Chem. 84, 3503 ~1980!; J.
Peyrelasse and C. Boned, Phys. Rev. A 41, 938 ~1990!.
@3# P. Guéring and B. Lindman, Langmuir 1, 464 ~1985!; U. Olsson, K. Shinoda, and B. Lindman, J. Phys. Chem. 90, 4083
~1986!; J-F. Bodet, J. R. Bellare, H. T. Davis, L. E. Scriven,and W. G. Miller, ibid. 92, 1898 ~1988!.
@4# W. Jahn and R. Strey, J. Phys. Chem. 92, 2294 ~1988!; P. K.
Vinson, W. G. Miller, L. E. Scriven, and H. T. Davis, ibid. 95,
2546 ~1991!.
@5# See, for instance, L. Auvray, in Micelles, Membranes, Microemulsions, and Monolayers, edited by W. M. Gelbart, D.
Roux, and A. Ben-Shaul ~Springer, New York, 1994!, p. 347.
@6# M. Kotlarchyk, S.-H. Chen, J. S. Huang, and M. W. Kim,
Phys. Rev. Lett. 53, 941 ~1984!.
@7# K. Kawasaki and K. Kawakatsu, Physica A 164, 549 ~1990!;
K. Kawakatsu and K. Kawasaki, ibid. 167, 690 ~1990!.
@8# M. Laradji, O. G. Mouristen, S. Toxvaerd, and M. J. Zucker-mann, Phys. Rev. E 50, 1243 ~1994!.
@9# M. Laradji, H. Guo, M. Grant, and M. J. Zuckermann, J. Phys.
A 24, L629 ~1991!; J. Phys. Condens. Matter 4, 6715 ~1992!.
@10# G. Pätzold and K. Dawson, Phys. Rev. E 52, 6908 ~1995!.
@11# J. H. Schulman and J. B. Montagne, Ann. N.Y. Acad. Sci. 92,
366 ~1961!.
@12# P. G. De Gennes and C. Taupin, J. Phys. Chem. 86, 2294
~1982!.
@13# T. Hashimoto and T. Izumitani, Macromolecules 26, 3631
~1993!; T. Izumitani and T. Hashimoto, ibid. 27, 1744 ~1994!.
@14# Y. Oono and S. Puri, Phys. Rev. Lett. 58, 836 ~1987!; Phys.
Rev. A 38, 434 ~1988!.
@15# S. Puri and Y. Oono, J. Phys. A 21, L755 ~1988!.
@16# A. Shinozaki and Y. Oono, Phys. Rev. E 48, 2622 ~1993!.
@17# M. Laradji, H. Guo, M. Grant, and M. J. Zuckermann, Phys.
Rev. A 44, 8184 ~1991!.
@18# H. Kodama and S. Komura, J. Phys. II France ~to be published!.

ACKNOWLEDGMENTS

We would like to thank Professor Shigehiro Komura, Dr.
G. Gompper, and Dr. J. L. Harden for their helpful discussions. We also thank Max-Planck-Institut für Kolloid undGrenzflächenforschungs for its hospitality while part of thiswork was in progress. This work is supported by the Ministry of Education, Science and Culture, Japan ~Grant-in-Aid
for Scientific Research No. 08226233, No. 08740324, and
No. 0657!.