popov

Evaporative deposition patterns: Spatial dimensions of the deposit Yuri O. Popov * Department of Physics, University of Chicago, 5640 S. Ellis Avenue, Chicago, Illinois 60637, USA sReceived 4 August 2004; published 28 March 2005 d A model accounting for the finite spatial dimensions of the deposit patterns in evaporating sessile drops of a colloidal solution on a plane substrate is proposed. The model is based on the assumption that the soluteparticles occupy finite volume and hence these dimensions are of steric origin. Within this model, the geo-metrical characteristics of the deposition patterns are found as functions of the initial concentration of thesolute, the initial geometry of the drop, and the time elapsed from the beginning of the drying process. Themodel is solved analytically for small initial concentrations of the solute and numerically for arbitrary initialconcentrations of the solute. The agreement between our theoretical results and the experimental data isdemonstrated, and it is shown that the observed dependence of the deposit dimensions on the experimentalparameters can indeed be attributed to the finite dimensions of the solute particles. These results are universaland do not depend on any free or fitting parameters; they are important for understanding evaporative depo-sition and may be useful for creating controlled deposition patterns. DOI: 10.1103/PhysRevE.71.036313 PACS number ssd: 47.55.Dz, 68.03.Fg, 81.15. 2z I. INTRODUCTION The problem of the so-called coffee-drop deposit has re- cently aroused great interest. The residue left when coffeedries on the countertop is usually darkest and hence mostconcentrated along the perimeter of the stain. Ringlike stains,with the solute segregated to the edge of a drying drop, arenot particular to coffee. Mineral rings left on washed glass-ware, banded deposits of salt on the sidewalk during winter,and enhanced edges in water color paintings are all examplesof the variety of physical systems displaying similar behav-ior and understood by coffee-drop deposit terminology. Understanding the process of drying of such solutions is important for many scientific and industrial applications,where the ability to control the distribution of the soluteduring the drying process is at stake. For instance, in thepaint industry, the pigment should be evenly dispersed afterdrying, and segregation effects are highly undesirable. Also,in protein crystallography, attempts are made to assembletwo-dimensional crystals by using evaporation-driven con-vection f1–3g, and hence solute concentration gradients should be avoided. On the other hand, in the production ofnanowires f4gor in surface patterning f5gperimeter- concentrated deposits may be of advantage. Recent impor-tant applications of this phenomenon related to DNAstretch-ing in a flow have emerged as well f6g. For instance, a high- throughput automatic DNA mapping was suggested f7g, where fluid flow induced by evaporation is used both forstretching DNA molecules and depositing them onto a sub-strate. Droplet drying is also important in the attempts tocreate arrays of DNA spots for gene expression analysis. Ringlike deposit patterns have been studied experimen- tally by a number of groups. Difficulties of obtaining a uni-form deposit f8g, deformation of sessile drops due to a sol- gel transition of the solute at the contact line f9,10g, stick- slip motion of the contact line of colloidal liquids f11,12 g, multiple ring formation f13g, and the effect of ring formation on the evaporation of the sessile drops f14gwere all reported. The evaporation of the sessile drops sregardless of the solute presence dhas also been investigated extensively. Constancy of the evaporation flux was demonstrated f15,16 g, and the change of the geometrical characteristics scontact angle, drop height, contact-line radius dduring drying was measured in detail f17–20 g. The most recent and complete experimental effort to date on coffee-drop deposits was conducted by Deegan et al. f21–24 g. Most experimental data referred to in this work originate from observations and measurements of this group.They reported extensive results on ring formation and dem-onstrated that these could be quantitatively accounted for.The main ideas of the theory of solute transfer in such physi-cal systems have also been developed in their work f21g.I t was observed that the contact line of a drop of liquid remainspinned during most of the drying process. While the highestevaporation occurs at the edges, the bulk of the solvent is concentrated closer to the center of the drop. In order toreplenish the liquid removed by evaporation at the edge, aflow from the inner to the outer regions must exist inside thedrop. This flow is capable of transferring all of the solute tothe contact line and thus accounts for the strong contact-lineconcentration of the residue left after complete drying. Thistheory is very robust since it is independent of the nature ofthe solute and only requires pinning of the edge during dry-ingswhich can occur in a number of possible ways: surface roughness, chemical heterogeneities, etc. d. Among other things, we will reproduce some of its results in this work. Mathematically, the most complicated task is related to determining the evaporation rate from the surface of thedrop. An analogy between the diffusive concentration fieldsand the electrostatic potential fields was suggested f25,26 g, so that an equivalent electrostatic problem can be solved in-stead of the evaporation problem. Important analytical solu- *Present address: Department of Physics, University of Michigan, 500 E. University Ave., Ann Arbor, MI 48109. Email address:[email protected] REVIEW E 71, 036313 s2005 d 1539-3755/2005/71 s3d/036313 s17d/$23.00 ©2005 The American Physical Society 036313-1 tions to this equivalent problem in various geometries were first derived by Lebedev f25g, and a few useful consequences from these analytical results were later reported in Ref. f27g. In this work, we continue development of the theory of solute transfer and deposit growth. Most previous works ad-dress the issue of the deposit mass accumulation at the dropboundary; however, they treat the solute particles as if theydo not occupy any volume, and hence all the solute can beaccommodated at the one-dimensional singularity of the con-tact line. In reality, the solute deposit accumulated at theperimeter has some thickness, and the shape of the residue ina round drop resembles a ring rather than the infinitely thincircumference of a circle. The earlier efforts were aimed atdescribing how the massof the contact-line deposit grows with time and how it depends on such geometrical character-istics of the drop as its radius sfor circular drops f21,22 gdor its opening angle and the distance from the vertex sfor pointed drops f28,29 gd. Little attempt has been made to de- scribe the geometrical characteristics of the contact-line de-posit itself, for instance, the width and the height of thedeposit ring. At the same time, there is solid experimentaldata f23,24 gon various geometrical characteristics of the ring and their dependence on time, the initial solute concen-tration, and the drop geometry. Here we develop a simplemodel that addresses this lack of understanding of the geo-metrical properties of the contact-line deposit and accountsfor the finite size of the deposit ring. We attribute the finitevolume of the deposit simply to the finite size of the soluteparticles, i.e., we assume the particles do occupy some vol- ume and hence cannot be packed more densely than a certainconcentration. The model is solved in the simplest case ofthe circular geometry both analytically and numerically, andthe results of the two methods are compared with the experi-mental data of Refs. f23,24 gsand with each other d. It turns out that this model is sufficient to explain most of the col-lected data. It should be noted that the model is as universaland robust sin its range of validity das the zero-volume theory of Deegan et al. f21,22 gsince it is based on essen- tially the same physical principles. The notion that the profile of the deposit could be found by the simple assumption that the solute becomes immobi-lized when the volume fraction reaches a threshold wasoriginally suggested by Dupont f30g. Efforts to create a model were conducted by Deegan f23,24 gwho formulated some physical assumptions, wrote them down mathemati-cally, and obtained some early-time exponents. Here wepresent the entire problem, including its full formulation andits analytical and numerical solutions. In the next section, weformulate the model, describe the system, and address someissues of the geometry and the evaporation rate. Then, wederive the governing equations from the conservation ofmass and later solve them analytically for small initial con-centrations of the solute and numerically for arbitrary initialconcentrations of the solute. A discussion section concludesthis work. II. MODEL, ASSUMPTIONS, AND GEOMETRY A. System We consider a sessile droplet of solution on a horizontal surface ssubstrate d. The nature of the solute is not essentialfor the mechanism. The typical diameter of the solute par- ticles in Deegan et al. f21–24 gwas of the order of 0.1–1 µm; we will assume a similar order of magnitude throughout thiswork. For smaller particles diffusion becomes important; forlarger particles sedimentation may play an important role. The droplet is bounded by the contact line in the plane of the substrate. This smacroscopic dcontact line is defined as the common one-dimensional boundary of all three phasessliquid, air, and solid substrate d. We will restrict our attention to the case of round drops, which is both of most practicalimportance and the easiest to treat mathematically. We assume that the droplet is sufficiently small so that the surface tension is dominant, and the gravitational effects canbe neglected. Mathematically, the balance of the gravita-tional force and the surface tension is controlled by the ratioof the smaximal dhydrostatic pressure rghmaxto the Laplace pressure 2 shmax/Ri2, where ris the fluid density, gis the gravitational constant, sis the surface tension at the liquid- air interface, Riis the drop radius in the plane of the sub- strate, and hmaxis the maximal height of the drop. For the typical experimental conditions this ratio rgRi2/2sis quite small sabout 0.25 d, and thus gravity is indeed unimportant and the surface shape is governed mostly by the surface ten-sion. Our treatment will produce the main-order term in theexpansion in this parameter, and since the parameter value isnot an order of magnitude smaller than 1, it may be neces-sary to construct correctional terms for better quantitativeagreement. For the present purposes, even the main termturns out to be sufficient to obtain agreement with the experi-mental results. Experimentally, the contact line remains pinnedduring most of the drying process. Therefore, we do not assume thatthe contact angle ubetween the liquid-air interface and the plane of the substrate is constant in time. A strongly pinnedcontact line can sustain a wide range of smacroscopic dcon- tact angles. The pinning mechanism can be described as self-pinning, i.e., pinning by the deposit brought to the contactline by the hydrodynamic flows caused by evaporation. Apinned contact line entails fluid flow toward that contact line.The “elasticity” of the liquid-air interface fixed at the contactline provides the force driving this flow. We will deal with small contact angles s u!1das this is almost always the case in the experimental realizations, including the experiments of Ref. f23gstypically, umax ,0.1−0.3 d. It will also be seen necessary to assume that the contact angle is small in order to obtain any analytical resultsin a closed form. A drop with a small contact angle is nec-essarilythin, i.e., its maximal height is much smaller than its radius and the slope of the free surface is small su=hu!1d. Thus, we consider small contact angles, or, equivalently, thin drops. We also consider slowflows, i.e., flows with low Rey- nolds numbers, which amounts to the neglect of the inertialterms in the Navier-Stokes equation. The free surface is described by the local mean curvature, which is spatially uniform at any given moment of time, butchanges with time as the droplet dries. Ideally, the surfaceshape should be considered dynamically together with theflow field inside the drop. However, as was shown earlierf29,31 g, for flow velocities much lower than the characteris-YURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-2 tic velocity v*=s/3hswhere sis the surface tension and h is the dynamic viscosity d, which is about 24 m/s for water under normal conditions, one can consider the surface shapeindependently of the flow and use the equilibrium result at any given moment of time for finding the flow at that time.Equivalently, the ratio of the viscous forces to the capillary forces is the capillary number Ca= hv˜/sswhere v˜is some characteristic value of the flow velocity, which is of the orderof 1−10 mm/sd, and this number is of the order of 10−8 −10−7under typical experimental conditions. Thus, the cap- illary forces are by far the dominant ones. B. Geometry and surface shape Cylindrical coordinates sr,f,zdwill be used throughout this work, as they are most natural for the geometry of inter- est.The origin is chosen in the center of the circular footprintof the drop on the substrate. The coordinate zis always nor- mal to the substrate, and the substrate is described by z=0, withzbeing positive on the droplet side of the space. The coordinates sr, fdare the polar radius and the azimuthal angle, respectively, so that the contact line is described by r=Ri, whereRiis the radius of the drop footprint. Due to the axial symmetry of the problem and our choice of the coordi-nates, no quantity depends on the azimuthal angle f. Our model pictures the drop as a two-component system sthe components being “the fluid” and “the solute” d, which has two “phases”: “the liquid phase” in the middle of thedrop and “the deposit phase” near the contact line. Bothcomponents are present in both phases, and the differencebetween the phases lies only in the concentration of the sol-ute in each phase. In the deposit phase, the volume fractionof the solute pis high and fixedin both space and time.Thus, pis just a constant number; one can think of it as comparable to the close-packing fraction or unity. sThe case of p=1 may seem to be special as there is no fluid in the deposit phase;however, for small initial concentrations of the solute thiscase will be seen to lead to exactly the same main-orderresults. dIn the liquid phase, the volume fraction of the solute xvaries in space and changes with time, and it is relatively small compared to p. The initial volume fraction xi=xs0dis constant throughout the drop; at later moments the solute gets redistributed due to the flows, and the concentrationbecomes different in different parts of the liquid phase. Thevolume fraction of the fluid is then s1−pdin the deposit phase and s1− xdin the liquid phase. Note that we do not require xi!pso far, although we do assume xi,p. It should also be emphasized that we do not presume there is any real“phase difference” between the so-called phases: one phaseis just defined as having the maximal reachable solute frac-tionpsthe solute cannot move in this phase dwhile the other phase is characterized by a lower solute fraction xsin thisphase the solute can move and hence its concentration can change in time and space d. Apart from this difference, the phases are essentially identical. The idea that the solute losesits mobility when its concentration exceeds some thresholdwas suggested by Dupont f30g. Since the drop is thin, we employ the vertically averaged flow velocity vsrd=1 htsrdE 0htsrd ussr,zddz, s1d whereussr,zdis the in-plane component of the local three- dimensional velocity usr,zd, andhtsrdis the thickness of the drop at distance rfrom the center. By making this approxi- mation, we implicitly assume that there is no vertical segre-gation of the solute, and thus we turn our model into aneffectively two-dimensional one. This is done mostly forsimplicity and is not expected to affect our main conclusionsssee Sec. V d. Within this model, it is natural to assume that the boundary between the phases is vertical. Thus, the par-ticles get stacked uniformly at all heights when they arebrought to the phase boundary by the hydrodynamic flowvsrd. This boundary can be pictured as a vertical wall at some radius Rstdfrom the center of the drop, and this wall propagates from the contact line flocated at R i=Rs0dgtoward the center of the drop. Figure 1 illustrates the mutual location of the two phases, and Fig. 2 schematically shows the timeevolution of the drying process and growth of the depositphase. The geometrical parameters of the model are shown in Fig. 3. The radius of the drop is R i, the radius of the phase boundary is Rstd, andRs0d=Ri. The height of the phase boundary is Hstd, and the initial condition is Hs0d=0 .I nt h e liquid phase, we conveniently split the total height of the free surfacehtsr,tdinto the sum of Hstdandhsr,td. FIG. 1. Mutual location of the two ‘‘phases’’in the drying drop: Lis ‘‘the liquid phase,’’ and Dis ‘‘the deposit phase.’’ FIG. 2. Time evolution of the deposit phase growth: side view sleftdand top view sright d. Only the deposit phase is shown. Thick- ness of the ring is exaggerated compared to the typical experimentalresults. FIG. 3. Geometry of the problem. Vertical scale is exaggerated in order to see the details; typically H!Riandh!Ri.EVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-3 SinceHis independent of r, the function hsrdsatisfies the Young-Laplace equation sthe statement of the mechanical equilibrium of the liquid-air interface d 2K=−Dp s, s2d where sis the surface tension, Dpis the pressure difference across the liquid-air interface, and Kis the mean curvature of the surface, uniquely related to the surface shape hby differ- ential geometry. For typical drying conditions, Dpandhvary with time slowly. As was shown earlier f29,31 g, it is suffi- cient to find the equilibrium surface shape first, and thendetermine the velocity field for this fixed functional form ofhwith time being just an adiabatic parameter, instead of solving for all the dynamical quantities simultaneously.Thus, the right-hand side of Eq. s2ddoes not depend on the local coordinates of a point within the drop salthough it does depend on time d, and the equation itself expresses the global condition of spatial constancy of the mean curvaturethroughout the drop. It defines the equilibrium surface shape at any given moment of time. The solution to Eq. s2dwith boundary condition hsRd=0 is just a spherical cap, and hence the shape of the upper part of the drop sabove the dashed line in Fig. 3 dis just a spheri- cal cap: hsr,td= ˛R2std sin2ustd−r2−Rstdcotustd. s3d Here ustdis the angle between the liquid-air interface and the substrate at the phase boundary, and functions Rstdandustd are related via the right-hand side of Eq. s2d: Rstd=2s Dpstdsinustd. s4d In the limit of small contact angles, u!1, the preceding expression adopts an even simpler form: hsr,td=R2std−r2 2Rstdustd+Osu3d. s5d Note that we do not assume that ustdandhsr,tdare neces- sarily positive at all times: both can be negative at later dry- ing stages, and the shape of the liquid-air interface may beconcave. Both convex and concave solutions for hsr,tdare consistent with Eq. s2d; the right-hand side of this equation can have either sign. By definition, both ustdandhsr,tdare positive when the surface is convex sand hence they are posi- tive at the beginning of the drying process dand negative when the surface is concave. The initial value of ustdcoin- cides with the initial contact angle ui=us0d. Clearly, there are three unknown functions of time in this geometry: ustd,Rstd, andHstd. However, these quantities are not independent of each other. Since we assume that the solute particles fill up the entire space between the substrateand the liquid-air interface when being brought to the phaseboundary, the three geometrical functions are related by theconstraintdH dt=−udR dt. s6d Physically, the angles between the liquid-air interface and the substrate are identical on both sides of the phase boundarys u=udH/dRud, and hence hsrdand its first derivative are con- tinuous past this boundary. Thus, there are actually only two independent functions of time, ustdandRstd. Condition s6d was first introduced by Deegan f23,24 g. The geometrical definitions above allow one to determine the volumes of each of the two phases. The volume of theliquid phase is simply V L=E 0Rstd fhsr,td+Hstdg2prdr=2pSR3u 8+R2H 2D+Osu3d. s7d Taking into account relation s6d, an infinitesimal variation of this volume can be expressed via the infinitesimal variationsof uandR: dVL=pR4 4dSu RD+2pHRdR. s8d The first term is responsible for the motion of the liquid-air interface, and the second term corresponds to the inwardshift of the phase boundary. It is also straightforward to ob-tain an expression for the differential of the volume of thedeposit phase, which has only the term related to the inwardshift of the phase boundary: dV D=−2 pHRdR. s9d We will use the last two expressions in the following section. We will also adopt the notation that the subscripts LandD refer to the liquid and deposit phases, respectively. C. Evaporation rate In order to determine the flow caused by evaporation, one needs to know the flux profile of the liquid leaving eachpoint of the surface. This quantity is independent of the pro-cesses going on inside the drop and must be determined priorto considering any such processes. The functional form of the evaporation rate Jsrdsdefined as the evaporative mass loss per unit surface area per unit timeddepends on the rate-limiting step, which can, in prin- ciple, be either the transfer rate across the liquid-vapor inter-face or the diffusive relaxation of the saturated vapor layerimmediately above the drop. We assume that the rate-limiting step is the diffusion of the saturated vapor. Indeed,the transfer rate across the liquid-vapor interface is charac-terized by a time scale of the order of 10 −10s, while the diffusion process has characteristic times of the order of Ri2/DswhereDis the diffusion constant for vapor in air and Riis the characteristic size of the drop d, which is of the order of seconds for water drops under typical drying conditions.The diffusion-limited evaporation rapidly attains a steadystate. Indeed, the ratio of the time required for the vapor-phase water concentration to adjust to the changes in thedroplet shape sR 2/Ddto the droplet evaporation time tfis ofYURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-4 the order of sns−n‘d/r<10−5, wherensis the density of the saturated vapor just above the liquid-air interface, n‘is the ambient vapor density, and ris the fluid density f27g, i.e., the vapor concentration adjusts rapidly compared to the evapo-ration time. As the rate-limiting process is the diffusion, the vapor densitynabove the liquid-vapor interface obeys the diffu- sion equation. Since the process is quasisteady, this diffusionequation reduces to the Laplace equation „ 2n=0. s10d This equation is to be solved together with the following boundary conditions: sadalong the surface of the drop the air is saturated with vapor and hence nat the interface is the constant density of the saturated vapor ns,sbdfar away from the drop the density approaches the constant ambient vapordensityn ‘, and scdthe vapor cannot penetrate the substrate and hence ]zn=0 at the substrate outside of the drop. Having found the vapor density, one can obtain the evaporation rateJ=−D„n, whereDis the diffusion constant. This boundary problem is mathematically equivalent to that of a charged conductor of the same geometry at constantpotential if we identify nwith the electrostatic potential and Jwith the electric field. Moreover, since there is no compo- nent ofJnormal to the substrate, we can further simplify the boundary problem by considering a conductor of the shapeof our drop plus its reflection in the plane of the substrate inthe full space instead of viewing only the semi-infinite spacebounded by the substrate sFig. 4 d. This reduces the number of boundary conditions to only two: sadn=n son the surface of the conductor, and sbdn=n‘at infinity. The shape of the conductor sthe drop and its reflection in the substrate dis now symmetric with respect to the plane of the substrate and re-sembles a symmetrical double-convex lens comprised of twospherical caps. This equivalent electrostatic problem of find-ing the electric field around the conductor at constant poten-tial in the infinite space is much simpler than the originalproblem in the semi-infinite space. The reflection techniquefor finding the evaporation field on the basis of the analogybetween the diffusion and the electrostatics was originallyused by Deegan et al. f21,22 g. Even in the circular geometry the equivalent problem is still quite complicated despite the visible simplicity. We con-sider an object whose symmetry does not match the symme-try of any simple orthogonal coordinate system of the three-dimensional space. In order to solve the Laplace equation,one is forced to introduce a special coordinate system sthe so-called toroidal coordinates dwith heavy use of special functions. The full solution to this problem is provided in theAppendix. The evaporation rate depends only on the overall shape of the drop, and evaporation occurs in the same fashion fromboth phases. We assume that the evaporation is not influ-enced by any motion of the solute inside the drop, and thenecessary amount of fluid can always be supplied to the re-gions of highest evaporation near the contact line. Physically,high evaporation near the edge is what brings the solute tothe contact line, and we assume that presence of the depositdoes not obstruct the motion of the fluid sFig. 5 d. Since the drop is thin and the contact angle is small, we will use theexpression Jsrd=2 pDsns−n‘d ˛Ri2−r2s11d for the evaporation rate sderived in the Appendix for no- solute drops in the limit u!1d, which has a reciprocal square-root divergence near the contact line as intuitivelyexpected from the electrostatics. The real situation may bedifferent from that assumed above when pis large or com- parable to 1, and the edge of the area where the evaporationoccurs may be located near the boundary of the phases in-stead of the contact line. However, for small initial concen-trations of the solute, the main-order result will be insensi-tive to the exact location of the singularity of the evaporationrate: whether it is located at the contact line or near theboundary of the phases. We will further comment on thiscase of the “dry deposit” when we obtain the full system ofequations. III. PRINCIPAL EQUATIONS A. Global conservation of fluid The essence of the entire theory can best be summarized in one sentence: “It is all about the conservation of mass.”Indeed, as we will see by the end of this section, all threegoverning equations obtained here represent the conservationof mass sor volume din one form or another. We start from the global conservation of fluid in the drop. Since the amount of solute within the drop does not changeduring the drying process, the change of the entire drop vol- FIG. 4. Illustration of the analogy between the evaporation rate Jfor a liquid drop and the electric field Efor a conductor. Consid- eration of the drop sor conductor dand its reflection in the plane of the substrate significantly simplifies the boundary problem. FIG. 5. Presence of particles in the deposit does not obstruct fluid evaporation at the edge of the drop. All the necessary fluid issupplied, and it is this motion of the fluid that brings the particles tothe deposit phase. Also shown schematically is the fact that theboundary between the phases is vertical and the particles getstacked at full height between the substrate and the free surface ofthe drop.EVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-5 ume is equal to the change of the amount of fluid. This fluid gets evaporated from the surface, and the total change of thefluid volume equals the amount evaporated from the surface: dV tot=dVF surf. s12d By convention, the superscripts FandSrefer to the fluid and the solute components, respectively swhile the subscripts L andDcontinue to denote phases d. The total change of the drop volume is the sum of the volume changes of eachphase: dV tot=dVL+dVD=pR4 4dSu RD, s13d wheredVLanddVDwere found in the preceding section fEqs. s8dands9dg. The volume of fluid evaporated from the surface can be determined from the known evaporation rate: UdVF dtU surf=−E 0RiJsrd r˛1+s]rhd22prdr=−4Dsns−n‘dRi r, s14d where ris the fluid density.We neglected the gradient of hsrd with respect to unity swhich is always legitimate for thin drops dand used Jsrdof Eq. s11d. Thus, Eqs. s12d,s13d, and s14dyield the first main differential equation of this section: R4d dtSu RD=−16Dsns−n‘dRi pr. s15d This equation represents the global conservation of fluid in the drop and relates the time dependencies of ustdandRstd. B. Local conservation of mass The next equation represents the localconservation of mass. There are two components in the liquid phase, andhence we write a separate equation for each of them. Since afree particle of the appropriate size reaches the speed of theflow in about 50 ns in water under normal conditions f32g, the solute particles are simply carried along by the flow, andthe velocities of each component are identical at each pointwithin the liquid phase fand equal to the depth-averaged fluid velocityvdefined in Eq. s1dg. The local conservation of fluid can be written in the form =·fs1− xdsh+Hdvg+J r˛1+s=hd2+]tfs1−xdsh+Hdg=0, s16d where xis the volume fraction of solute at a given point within the liquid phase, and each of the quantities sh +Hd,J,x, andvis a function of distance rand time t.fWe drop the s=hd2part of the second term everywhere in this work since it is always small compared to unity for small contact angles. gThis equation represents the fact that the rate of change of the fluid amount in a volume element scolumn d above an infinitesimal area on the substrate sthird term dis equal to the negative of the sum of the net flux of fluid out ofthe column sfirst term dand the amount of fluid evaporatedfrom the surface element on top of that column ssecond termd; Fig. 6 illustrates the idea. A similar equation can also be written for the local conservation of solute, but without the evaporation term: =·f xsh+Hdvg+]tfxsh+Hdg=0. s17d Adding the two equations and employing the linearity of the differential operations, one obtains =·fsh+Hdvg+J r+]tsh+Hd=0. s18d This relation could have been obtained if we considered only one component with volume fraction 1 in the liquid phase,and this equivalence should be of no surprise: when the sol-ute moves in exactly the same fashion as the fluid does, anydifferentiation between the two is completely lost sfrom the point of view of the conservation of volume d. Note that if evaporation were too intense, this equivalence would nothold, as there might be an insufficient amount of fluid com-ing into a volume element, and the solution could get com-pletely dry sonly the solute component would be left d.W e implicitly assume this is not the case for our liquid phasewhere the solute fraction is relatively small and the evapora-tion is not too strong. In circular geometry, due to the symmetry, the flow is radial and independent of f. Thus, Eq. s18dcan be resolved with respect to the radial component of the velocity: vrsr,td=−1 rsh+HdE 0rSJ r+]th+]tHDrdr. s19d Straightforward integration with hsr,tdof Eq. s5d,dH/dtof Eq.s6d, andJsrdof Eq. s11d, and employment of Eq. s15dfor dsu/Rd/dtyield vrsr,td=2Dsns−n‘d prRi r˛1−sr/Rid2−f1−sr/Rd2g2 sRu/2df1−sr/Rd2g+H. s20d This expression for the flow velocity at each point rwithin the liquid phase in terms of the time-dependent geometrical FIG. 6. Conservation of mass: the liquid-vapor interface is low- ered exactly by the amount of fluid evaporated from the surfaceplus the difference between the outflow and the influx of fluid fromthe adjacent regions.YURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-6 characteristics of the drop ustd,Rstd, andHstdis a direct consequence of the local conservation of mass. With the velocity in hand, we can compute the time it takes an element of fluid initially located at distance rifrom the center to reach the contact line. First, only the particlesinitially located near the contact line reach that contact line.As time goes by, the particles initially located further awayfrom the contact line and in the inner parts of the drop reachthe contact line. Finally, the particles initially located in theinnermost parts of the drop si.e., near its center dreach the contact line as well. The more time elapses, the more par-ticles reach the contact line and the larger the area is wherethey were spread around initially. One can view this processas inward propagation of the inner boundary of the set of theinitial locations of the particles that have reached the contactline by time t. As is easy to understand, the velocity of this front is equal to the negative of the vector of the fluid veloc-ity at each point sthe fluid and the particles move together toward the contact line while this front moves away from it;hence the minus sign d.We label by r istdthe initial location of the solute particles that reach the phase boundary sand be- come part of the deposit ring dat timet. Since the solute particles from the outer areas of the drop reach the depositphase sooner than the particles from the inner areas, thisfunction is monotonically decreasing, and its derivative issimply related to vrfound in the preceding paragraph fEq. s20dg: dri dt=−vrsri,td. s21d Thus, the second principal equation of this section is dri dt=−2Dsns−n‘d prRi ri˛1−sri/Rid2−f1−sri/Rd2g2 sRu/2df1−sri/Rd2g+H. s22d This equation relates ristdto the time dependencies of the geometrical parameters of the drop fustd,Rstd, andHstdg. C. Global conservation of solute The volume of solute in the deposit phase VDSat timetis equal to the volume of solute outside the circle of radius ristd at time 0 fsince all the solute between ristdandRibecomes part of the deposit by time tg.The latter volume can be found by integrating hsr,0dover the area swept by the fluid on its way from rito the contact line and multiplying the result by the initial volume fraction of solute xi: VDS=xiE riRi hsr,0d2prdr=VSF1−Sri RiD2G2 , s23d whereVS=pxiRi3ui/4 is the total volume of solute in the drop. On the other hand, the volume of solute in the depositphase is just the constant fraction pof the volume of the entire deposit phase: V DS=pVD. s24d Equating the right-hand sides of these two equations, taking the time derivatives of both sides, and making use of thealready determined dVDof Eq. s9d, we obtain the third prin- cipal equation of this section: xiRi3uiF1−Sri RiD2Gd dtF1−Sri RiD2G=−4pHRdR dt.s25d This equation represents the global conservation of solute in the drop. Thus, we have four unknown functions of time, ustd,Rstd,Hstd, andristd, and four independent differential equations for these functions, Eqs. s6d,s15d,s22d, and s25d. In reality, we need only three of these functions, ustd,Rstd, andHstd; however, there is no simple way to eliminate ristd from the full system and reduce the number of equations. Having solved this system of equations, we will be able tofully characterize the dimensions of the deposit phase anddescribe the evolution of the deposit ring. The following sec-tion is devoted to the details and the results of this solution. Here we will only comment on how this system changes in the case of the completely dry deposit p=1. In this case there is no evaporation from the surface of the deposit phase,and the effective edge of the evaporating area is somewherein the vicinity of the phase boundary. Assuming the samereciprocal square-root divergence of the evaporation rate atr=Rinstead of r=R ifwhich mathematically means substitu- tion ofRin place of Riin Eq. s11dgand conducting a deri- vation along the lines of this section, one can obtain a verysimilar system of four differential equations.These equationswould be different from Eqs. s6d,s15d,s22d, and s25din only two minor details. First, Eqs. s15dands22dwould lose all indicesiat all occasions of R isi.e., one should substitute R for allRiin both equations d. Second, pshould be set to 1 in Eq.s25d.Apart from these details, the two systems would be identical.As we will see in the following section, this differ-ence between the two systems is not important in the mainorder in a small parameter introduced below, and thus this“dry-deposit” case does not require any special treatment,contrary to the intuitive prudence. IV. RESULTS A. Analytical results in the limit of small initial concentrations of the solute So far we have not introduced any small parameters other than the initial contact angle ui!1. In particular, Eqs. s6d, s15d,s22d, and s25dwere obtained without assuming any relation between pandxiother than the nonrestrictive con- dition xi,p. In order to find the analytical solution to this system, we will have to assume that xi!p. Then, we will solve the same system of differential equations numericallyfor an arbitrary relation between xiandp. The assumption xi!pphysically means that the solute concentration in the liquid phase is small—it is much smallerthan the concentration of close packing or any other compa-rable number of the order of 1. This is the case for mostpractical realizations of the ring deposits in experiments andobservations: the solute concentration rarely exceeds 10% ofvolume, and in most cases it is far lower. If the volumefraction of the solute is small, then the solute volume is alsoEVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-7 small compared to the volume of the entire drop. Hence, the deposit phase, which consists mostly of the solute, must alsohave a small volume compared to the volume of the entiredrop. Thus, if the initial volume fraction xiis small, then the dimensions of the deposit ring must be small compared tothe corresponding dimensions of the entire drop. Let us now introduce a parameter ethat is small when xi/pis small. We do not fix its functional dependence on xi/pfor the moment: e=fSxi pD!1, s26d wherefis an arbitrary increasing function of its argument. Then we postulate that the ring width is proportional to thisparameter: Rstd=R if1−eW˜stdg, s27d whereW˜stdis an arbitrary dimensionless function and we explicitly introduced the dimensionality via Ri. Obviously, W˜s0d=0. So far we simply wrote mathematically that the ring width is small whenever the initial volume fraction of the solute is small. Next, we introduce a dimensionless vari-able for the angle ustd: ustd=uiu˜std, s28d where both ustdanduiare small, while the newly introduced function u˜stdis arbitrary fin particular, u˜s0d=1g. Due to the geometricalconstraint s6d,theheightofthering Hstdmustbe linear in the small parameters eanduiand directly propor- tional to the only dimensional scale Ri: Hstd=euiRiH˜std, s29d whereH˜stdis yet another dimensionless function of time fH˜s0d=0g, related to the functions W˜stdandu˜stdby an ex- pression similar to Eq. s6d. The last dimensionless variable is introduced in place of the fourth unknown function ristd: V˜std=1−Sristd RiD2 , s30d with the initial condition V˜s0d=0. Thus, we introduced four new dimensionless variables in place of the four original ones and explicitly separated their dependence on small pa-rameters eandui. Finally, we define the dimensionless time tas t=t tf, s31d wheretfis a combination of system parameters with the dimensionality of time: tf=prRi2ui 16Dsns−n‘d. s32d In the limit xi/p!0 this combination represents the time at which all the solute reaches the deposit phase; for finite xi/p it does not have such a simple interpretation.Substitution of all the definitions of the preceding para- graph into the original system of equations s6d,s15d,s22d, ands25dand retention of only the leading and the first cor- rectional terms in eyield the following simplified system of equations: dH˜ dt=u˜dW˜ dt, s33d du˜ dt+eu˜dW˜ dt−3eW˜du˜ dt=−1, s34d dV˜ dt=˛V˜−V˜2f1−4 eW˜sV˜−1−1dg 2u˜V˜f1−eW˜s2V˜−1−1dg+4eH˜, s35d xi pV˜dV˜ dt=4e2H˜dW˜ dts1−eW˜d. s36d As is apparent from the last equation, the parameter e2must be proportional to xi/p. Since the separation of the ring width into eandW˜in Eq. s27dis absolutely arbitrary, the parameter eis defined up to a constant multiplicative factor. Therefore, we setthis factor in such a way that e2isequalto xi/p: e=˛xi p. s37d This fixes the function ffrom the original definition s26d. The differential equations in the system s33d–s36dare still coupled. However, in the main szeroth dorder in e, the equa- tions clearly decouple: the second equation can be solved with respect to u˜stdindependently of all the others, then the third equation can be solved with respect to V˜stdindepen- dently of the first and the fourth, and finally the first and the fourth equations can be solved together as well. Thus, onecan obtain the following main-order solution to the system ofequations above with the appropriate initial conditions: u˜std=1− t, s38d V˜std=f1−s1−td3/4g2/3, s39d H˜std=˛1 3FBS7 3,4 3D−Bs1−td3/4S7 3,4 3DG, s40d W˜std=E 0t1 8H˜st8df1−s1−t8d3/4g1/3 s1−t8d1/4dt8. s41d HereBsa,bd=e01xa−1s1−xdb−1dxis the complete beta func- tion,Bzsa,bd=e0zxa−1s1−xdb−1dxis the incomplete beta func- tion sa.0,b.0, and 0 łzł1d, and the integral in the last equation cannot be expressed in terms of the standard el-ementary or special functions. In a similar fashion, systemsof equations of the higher orders in ecan be written fonly the first-order corrections are kept in the system s33d–s36dg, andYURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-8 the higher-order terms can also be constructed up to an arbi- trary order. A system of equations similar to our system s33d–s36d was presented by Deegan in f23,24 g. However, some terms of the first order in concentration were missing and no ana-lytical solution to the system of equations was obtained inthose works. Here we derive the equations in a systematicway and provide their analytical solution. How do our results s38d–s41dtranslate into the original variables? The first two of them fEqs. s38dands39dgrepro- duce earlier results. In terms of the dimensional variables Eq.s38drepresents the linear decrease of the angle between the liquid-air interface and the substrate at the phase boundarywith time: ustd=uiS1−t tfD s42d fplotted by the solid line in Fig. 8 sadbelow g. This is a direct analog of Eq. sA11dfor the contact angle in the no-solute case, as is clear from the definition of tffEq.s32dg. So, the angle uin the case of the finite-volume solute depends on time in exactly the same fashion as the contact angle in theno-solute case does. This expression also provides an inter-pretation of t f: it is the time at which the free surface of the liquid phase becomes flat. Before tfthis surface is convex, aftertfit becomes concave and bows inward suntil it touches the substrate d. Thus,tfis generally notthe total drying time. In the limit xi/p!0 the height of the deposit ring is going to zero and the two times are the same. For finite values of thisparameter the total drying time is longer than the time atwhich the liquid-air interface becomes flat. Equation s42dhas been verified in the experiments f23,24 gwhere the mass of the drop was measured as a function of time sFig. 7 d. Since the mass of a thin drop is directly proportional to u, these results confirm the linearity of ustdduring most of the drying process. The second equation s39dhas a direct analog in the case of the zero-size solute particles. In the original variables itcan be rewritten asS1−t tfD3/4 +F1−Sristd RiD2G3/2 =1, s43d which is identical to Eq. s3.24dof Ref. f31gobtained for the zero-volume solute. Clearly, ri=Riwhent=0, andri=0 when t=tf. According to Eqs. s23dands30d, the fraction of solute in the deposit phase VDS/VSis VDS VS=V˜2=F1−S1−t tfD3/4G4/3 s44d fplotted by the solid line in Fig. 8 sbdg. This fraction is 0 at t=0 and becomes 1 at t=tf.Thus,tfcan also be interpreted as the time at which all the solute particles become part of thedeposit phase. So far, the results of this finite-volume modelcoincide with the results of the zero-volume case consideredearlier f21,22 g. However, the third and the fourth equations fEqs. s40dand s41dgrepresent additional results. In the dimensional vari- ables they yield the height of the phase boundary Hand the width of the deposit ring W;R i−R, respectively, Hstd=˛xi puiRiH˜St tfD, s45d Wstd=˛xi pRiW˜St tfD, s46d where the functions H˜stdandW˜stdare plotted in Figs. 8 scd and 8 sddsthe solid curves d. These results provide the sought dependence of the geometrical characteristics of the depositring on all the physical parameters of interest: on the initialgeometry of the drop sR ianduid, on the initial solute con- centration sxid, and on the time elapsed since the beginning of the drying process std. If the time is considered as a pa- rameter, they can also be used to obtain the geometrical pro- file of the deposit si.e., the dependence of the height on the width d, which we plot by the solid line in Fig. 9. Note that the vertical scale of this plot is highly expanded compared tothe horizontal scale since there is an extra factor of ui!1i n the expression for the height; in the actual scale the height ismuch smaller than it appears in Fig. 9. It is straightforward to obtain the asymptotics of H˜stdand W˜stdfor early and late drying stages.At early times, both the height and the width scale with the drying time as a power law with exponent 2/3: H<˛xi puiRis3td2/3 27/3f1+Ostdg s t!1d, s47d W<˛xi pRis3td2/3 27/3f1+Ostdg s t!1d. s48d Thus, at early times H<uiW, which can also be deduced directly from Eq. s6dwithout obtaining the complete solution above. fThe early-time exponent 2/3 was first obtained by Deegan f23gwithout deriving the full time dependence s40d ands41d.gAt the end of the drying process, the height and the width approach finite values swhich, apart from the di- FIG. 7. Mass of a drying drop as a function of time. Experimen- tal results, after Refs. f23,24 g. The line running through the data is a linear fit. sCourtesy Robert Deegan. dEVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-9 mensional scales, are universal, i.e., constants dand do so as power laws of stf−tdwith two different exponents: H<˛xi puiRiFH˜s1d−s1−td7/4 14H˜s1d+Os1−td5/2Gs49d s1−t!1d, W<˛xi pRiFW˜s1d−s1−td3/4 6H˜s1d+Os1−td3/2Gs1−t!1d, s50d whereH˜s1dandW˜s1dare simply numbers: H˜s1d=˛1 3BS7 3,4 3D<0.297, s51d W˜s1d=E 011 8H˜stdf1−s1−td3/4g1/3 s1−td1/4dt<0.609. s52d Clearly,dH/dW=uis1−tdand hence vanishes when t!1. This fact can also be observed in the flattening of the ana- lytical graph in Fig. 9 at late times. The dependence of the height and the width on the radius of the drop Ri, while intuitively obvious ssinceRiis the only scale in this problem with the dimensionality of the length d, has been verified in experiments f23,24 g.Alinear fit has beenobtained for the dependence of the ring width on the radius, in exact agreement with our findings. Comparison to the experimental data for the dependence on the initial concentration of the solute is slightly less trivial. Our results predict that both the height Hand the widthWscale with the initial concentration as xi1/2sat least, in the leading order for small concentrations d.The same scal- ing prediction was also made by Deegan f23,24 g. However, his experimental results show a different exponent of xi: the values 0.78±0.10 and 0.86±0.10 were obtained for two dif- FIG. 8. Results for the dependence of the geometrical characteristics of the drop on time t. In each plot, the solid curve is the analytical result in the limit xi/p!0, while the other curves are the numerical results. Different numerical curves correspond to different initial concentrations of the solute; values of parameter xi/pare shown at each curve. sadAngle ubetween the liquid-air interface and the substrate at the phase boundary. sbdVolume fraction of the solute in the deposit phase VDS/VS.scdHeight of the phase boundary Hsin units of uiRi˛xi/pd.sddWidth of the deposit ring Wsin units of Ri˛xi/pd. FIG. 9. Deposit ring profile: dependence of the height of the phase boundary Hon the width of the deposit ring W. The solid curve is the analytical result in the limit xi/p!0; the other curves are the numerical results. Different numerical curves correspond todifferent initial concentrations of the solute; values of the parameter xi/pare shown at each curve. The vertical scale is different from the horizontal scale by a factor of ui!1.YURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-10 ferent particles sizes sFig. 10 d. Why the difference? The an- swer lies in the fact that the width measured in the experi-ments f23,24 gis not the full width of the ring at the end of the drying process, but rather the width of the ring at depin- ning. Depinning is the process of detachment of the liquid phase from the deposit ring sFig. 11 d. This detachment was observed experimentally in colloidal suspensions but has notbeen explained in full theoretically yet. 1An important obser- vation, however, is that the depinning time si.e., the time at which the detachment occurs and the ring stops growing d depends on the initial concentration of the solute. This de-pendence was also measured by Deegan sFig. 12 d, and the resulting exponent was determined to be 0.26±0.08. Thus,the width of the ring at depinning W dscales with the initial concentration of the solute xias Wd~xi1/2W˜Std tfD~xi0.5W˜sxi0.26±0.08 d, s53d wheretdis the depinning time std/tf~xi0.26±0.08d.As is appar- ent from Fig. 12, the typical values of the depinning time are of the order of s0.4−0.8 dtf. In this time range, the function W˜stdis virtually linear fthe analytical curve in Fig. 8 sddg. Therefore, the dependence of Wdonxihas an overall expo- nent of the order of 0.76±0.08. It is now clear that bothexperimental results 0.78±0.10 and 0.86±0.10 fall withinthe range of the experimental uncertainty of this approximatepredicted value, and the theoretical dependence of the ringwidth on the initial concentration agrees with the experimen-tal results quite well. Note that Deegan f23,24 gdid not report direct measure- ments of the height of the deposit ring. The height was 1While the full explanation is yet to be developed, the naive rea- son for the depinning seems relatively straightforward. The pinning force depends only on the materials involved and is relatively in-sensitive to the value of the contact angle. At the same time, thedepinning force is simply the surface tension, which is directed along the liquid-air interface and which increases as the contactangle decreases ssince only the horizontal component of this force is important d. Thus, the relatively constant pinning force cannot compensate for the increasing depinning force, and after the contactangle decreases past some threshold, the depinning force wins andcauses the detachment. FIG. 10. Ring width normalized by the drop radius vs initial concentration of the solute for two different particle sizes. Experi-mental results, after Refs. f23,24 g. The two data sets are offset by a factor of 5 to avoid mixing of the data points related to the differentparticle sizes. The lines running through the data are linear fits inthe double-logarithmic scale, which upon conversion to the linearscale yield power laws with exponents 0.78±0.10 and 0.86±0.10.sCourtesy Robert Deegan. d FIG. 11. A photographic sequence demonstrating a depinning event. Experimental results, after Refs. f23,24 g. The view is from above, and the solid white band in the lower part of the frame is thering; the rest of the drop is above the ring. The time between thefirst and the last frames is approximately 6 s; the major axis of thehole is approximately 150 µm.sCourtesy Robert Deegan. d FIG. 12. Depinning time normalized by the extrapolated drying time vs initial concentration of the solute. Experimental results,after Refs. f23,24 g. The line running through the data is a linear fit in the double-logarithmic scale, which upon conversion to the linearscale yields a power law with exponent 0.26±0.08. sCourtesy Rob- ert Deegan. dEVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-11 calculated from the data in hand, and thus a direct compari- son to the experimental data is not available for the height. The square-root dependence of the height and the width on the concentration is in good agreement with generalphysical expectations. Indeed, the volume of the deposit ringis roughly proportional to the product of the height and thewidth. On the other hand, the height is of the same order ofmagnitude as the width since the ratio of the two is of theorder of uiswhich is a constant d. Thus, both the height and the width scale approximately as the square root of the ringvolume. Finally, the volume of the deposit ring is propor-tional to the initial volume fraction of the solute: the moresolute is present initially, the larger is the volume of thedeposit ring at the end. Therefore, both the height and thewidth must scale as the square root of the initial volumefraction. It is rewarding that our complex calculation leads tothe same results as this simple physical argument. Thus, the complete analytical solution to our model is available in the limit xi/p!0, and this solution compares favorably with the experimental data. Since the main-ordersolution in xi/pis perfectly adequate, the difference between the original system of equations and the one for the “com-pletely dry” case is not important. Indeed, the main-orderresults are identical in both cases, because one case is differ-ent from the other only by the presence of Rinstead of R iin a few places in the main equations, and this difference is ofthe correctional order in xi/p. B. Numerical results for arbitrary initial concentrations of the solute Apart from approaching the original system of equations s6d,s15d,s22d, and s25danalytically, we also solve it numeri- cally. During this numerical procedure we do not presumethat xi/pis small, nor do we expand any quantities or equa- tions in eor any other small parameters. Our main purpose is to reproduce the results of the first part of this section and todetermine the range of validity of our analytical asymptotics. The typical values of xi/pin most experimental realiza- tions are of the order of 0.001-0.01 and thus, only the con-centrations below approximately 0.1 are of practical interest.sNote that xi/p=0.1 corresponds to a quite substantial value of the small parameter e=˛0.1<0.32. dThus, we will con- centrate on this range of xi/pwhen describing the results despite the fact that the numerical procedure can be sand has been dconducted for any ratio xi/p. The general trend is il- lustrated well by the results in this range of concentrations.In the case of xicomparable to pour model is not expected to produce any sensible results, as the entire separation of thedrop into the two phases sthe liquid phase and the deposit phase dis based on the assumption that the mobility of the solute is qualitatively different in the two regions.When xiis comparable to pthe two phases are physically indistinguish- able, while the model still assumes they are different. We present our numerical results for the same quantities sand in the same order das in our analytical results s42dand s44d–s46d. Since for arbitrary xi/pthe timetfisnotexactly the total drying time, there is a question of where sat what timedto terminate the numerical curves. By convention, weterminate all the curves in all the graphs at the value of t/tf when all the solute reaches the deposit phase. In our model, it turns out that the time the last solute particles reach thedeposit ring and the time the center point of the liquid-airinterface touches the substrate are about the same. For all theinitial concentrations, the two times were numerically foundto be within 0.1% of each other, and the curves are termi-nated at exactly this moment. Of course, in reality a smallfraction of solute should stay in the liquid phase as long asthe liquid phase exists, and so the moment the last soluteparticles reach the deposit phase should be afterthe moment the center point touches the substrate; however, the amountof solute remaining in the liquid phase at touchdown is in-significant, and practically all the deposit has alreadyformed. Numerical results for angle uas a function of time are shown in Fig. 8 sad.All the curves behave almost linearly sas expected d; however, the slope increases with concentration: formation of the ring in the drops with more solute finishesfaster son the relative scale of t fd. The end of each curve demonstrates the value of the angle utat the moment the liquid-air interface touches the substrate. The analytical ex-pression for this angle is ut=−2H/Rfor a thin drop. The absolute value of this angle increases with concentration,which is quite natural since for small concentrations theheight of the ring grows as a square root of the concentrationwhile the radius of the liquid phase does not change substan-tially. Clearly, the numerical results converge to the analyti-cal curve when xi/p!0. Growth of the volume fraction of solute in the deposit phaseVDS/VSwith time is shown in Fig. 8 sbdfor various solute concentrations. This graph reconfirms the observationof the preceding paragraph that the solute transfer happensfaster sin units of t fdfor denser colloidal suspensions. All curves are terminated when volume fraction VDS/VSbecomes equal to 1. sThe apparent termination of the curve for xi/p =0.1 earlier than that is an artifact of the plotting software. d As the corresponding analytical results do, the numericalplots of Figs. 8 sadand 8 sbdshould presumably hold true independently of the geometrical details of the solute accu-mulation in the ring swhich cannot be expected from the following plots for the ring height and width d. The next two graphs represent the numerical results for the height fFig. 8 scdgand the width fFig. 8 sddgof the deposit ring as functions of time. The ring profile, i.e., the depen-dence of the height on the width, is also shown in Fig. 9.Asthe graphs depict, the ring becomes wider and lower sin the reduced variables dfor higher initial concentrations of the sol- ute. Since the volume of the ring is roughly proportional tothe product of the height and the width, the decrease inheight must be of the same magnitude as the increase inwidth. This can be qualitatively observed in the graphs. As a final piece of the numerical results, we create a double-logarithmic plot for the dependence of the height andthe width on the initial concentration of the solute sFig. 13 d. The predicted square-root dependence on the initial concen-tration is seen to hold true for volume fractions up to ap-proximately 10 −1/2pfor the height and up to approximately 10−3/2pfor the width. The deviations for higher volume frac- tions are due to the increasing role of the correctional termsYURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-12 inecompared to the main-order terms represented by the solid lines. In this graph, as in all the results of this section,it is clear that our main-order analytical results provide anadequate description of all the functional dependencies in therange of the initial concentrations of experimental impor-tance s0.001-0.01 d. In general, our numerical results complement and rein- force our analytical results, providing a cross-check of bothmethods. V. DISCUSSION Both the analytical results of Eqs. s42dands44d–s46dand the numerical graphs of Figs. 8 and 9 may be reproducedexperimentally giving validation to the proposed model.Measurements of the profiles in Fig. 9 should be particularlyeasy to conduct ssince there is no time dependence involved d and may confirm or refute the predicted robustness and theuniversality of the deposition profiles. While the main principles of the proposed model were laid down by Deegan f23,24 g, its analytical solution for small concentrations and its numerical solution for arbitrary con-centrations are obtained here. The availability of the exactanalytical solution demonstrated that the theoretical scalingof the deposit width at depinning with solute concentrationindeed agrees with both measured values of the exponents.The earlier estimate of Refs. f23,24 gbased solely on the early-time exponent had an overlap with only one of theconcentration exponents. Our numerical results also quantifythe range of solute concentrations where the predictedsquare-root dependence of the width holds true. In general,Deegan’s results were not sufficient to obtain the proper scal-ing of the width with time anywhere beyond the early dryingstages; the results of this work provide that time scaling at alldrying stages. All the presented results suggest that the de-posit ring profile and its growth can be fully accounted for onthe basis of the finite volume of the solute particles only andthat the governing functional dependences are universal. One may notice that the curves in Fig. 9 end at some positive snonzero dheight. This indicates that the solute is exhausted before the profile curves have a chance to return tothe substrate, and the final shape of the deposit ring must have a vertical wall at its inner side. We believe this is anartifact of our model, which is inherently two dimensionalwhen flows inside the drop are concerned. Thus, the verticaldistribution of the solute was assumed homogeneous sthe phase boundary is vertical and the particles get stacked uni-formly at all heights d, and we used the depth-averaged ve- locity s1dthroughout this work. This is equivalent to assum- ing that vertical mixing is complete. This assumption is quiteimportant, and the results are expected to get modified if thetrue three-dimensional velocity profile is used instead of thedepth-averaged velocity. We expect that if a three-dimensional model were built and the dependence on zwere taken into account for all the quantities then the discontinu-ous wall of the phase boundary would get smoothed and theheight would continuously return to zero. The question re-mains whether such a model would be solvable analytically. Our model relies on the assumption that solute mobility is different in the so-called liquid and deposit phases. In es-sence, we assume that the mobility is 0 in the deposit phaseand 1 in the liquid phase. This assumption, while artificial inits nature, seems relatively reasonable when applied to thissystem. Indeed, in physical situations near close packing, theloss of mobility typically occurs over a quite narrow range ofthe concentration values, and hence our assumption shouldwork satisfactorily when the difference between xiandpis orders of magnitude. The higher the initial concentration isand the closer the two values are, the less this assumptionholds true and the more artificial the difference between thetwo phases is. Thus, the validity of any model based on thisseparation of the mobility scales decreases for higher initialconcentrations of the solute. The model assumes that the free-surface slope between the liquid and the deposit phases is continuous. In fact, as-sumption s6dexpressing this continuity is one of the four basic equations of this work. This assumption seems quitenatural as well. Indeed, if the liquid is present on both sidesof the phase boundary, the change in the slope of its freesurface would cost extra energy from the extra curvature atthe phase boundary, since the liquid-air interface possesseseffective elasticity. The presence of this extra energy sor the extra pressure dat the location of the phase boundary is not justified by any physical reasons as all the processes are slowand the surface is in equilibrium. In equilibrium, the surfaceshape must have constant curvature past the phase boundarysince the entire separation into the two phases is quite artifi-cial as discussed above. The presence of the particles belowthe liquid-air interface does not influence the surface tension,and thus the liquid surface sand its slope dshould be continu- ous at the phase boundary. If the density of the particlesmatches the density of the liquid swhich was the case in the experiments d, nothing prevents the particles from filling up the entire space between the substrate and the liquid-air in-terface, thus providing the growth of the upper edge of thedeposit phase alongthe liquid-air interface. This is particu- larly true for the thin drops discussed here, where verticalmixing is intensive, where the free surface is nearly horizon-tal, and where the problem is essentially two dimensional.However, the equality of the slopes on both sides of thephase boundary does not seem inevitable, and one may think FIG. 13. Numerical results: log-log plot of the dependence of the height of the phase boundary Hand the width of the deposit ring Won the initial volume fraction of the solute xi. The main-order analytical results H~˛xi/pandW~˛xi/pare also provided for comparison.EVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-13 of situations when it does get violated. One example might be the late drying times, when the deposit growth is very fastfFig. 8 sbdgand hence the deposition may occur in some non- regular manner inconsistent with this slow-process descrip-tion. Other examples may be related to gravity sslightly un- equal densities of the particles and the fluid dor convection. This assumption can possibly be checked experimentally,and if condition s6dis found violated, an equivalent con- straint dependent on the details of the deposit-growth mecha-nism must be constructed in place of Eq. s6d. Another inherent assumption of our model is related to the evaporation rate Jsrd. Presence of the solute inside the drop was assumed not to affect the evaporation from its sur- face. This is generally true when the evaporation is not toofast and the deposit phase is not too thick and not too con-centrated. When these conditions are not obeyed, the pres-ence of a thick or concentrated layer of the solute in the wayof the liquid moving from the phase boundary to the contactline may create a strong viscous force. This viscous forcewould prevent the necessary amount of fluid from being sup-plied to the intense-evaporation region near the contact line. Generally, we assumed throughout this work that the viscousstresses are not important, and this is valid whenever v !s/3h. In the deposit phase, the velocity is large due to the proximity to the contact-line divergence of the evaporationrate, and the effective viscosity is large due to the high con-centration of the solute.Thus, this condition may get violatedand the viscosity may become important in the depositphase, slowing down the supply of the liquid and ultimatelymaking the deposit dry. Obviously, this affects the evapora-tion rate, and the functional form of the evaporation profilechanges. The simple assumption that the evaporation ratestays of the same functional form, but with the divergence atthe phase boundary satRdinstead of the contact line satR id, was shown above not to affect our main-order results. Thus,our results appear to be relatively insensitive to the exactlocation of this divergence within the snarrow ddeposit phase. sIn reality the evaporation edge would be somewhere be- tween the original contact line and the phase boundary, i.e.,the real situation is intermediate between the two consid-ered. dHowever, the deposit could modify the evaporation rateJsrdin other ways. When there is a dry deposit ring just outside the liquid phase, the entire functional form of Jmay change, and the Laplace equation for an equivalent electro-static problem must be solved anew with additional bound-ary conditions responsible for the presence of the dry soluterim and the modified evaporation at the edge.As we show inthe Appendix, this is the most complicated part of the prob-lem, and the mathematics can become prohibitively complex.Thus, finding the exact form of Jmay be a formidable task. One way around is in creating such evaporating conditionsthat the functional profile is simpler, for instance, Jis just a constant. This would be more difficult to control experimen-tally, but would be much easier to treat analytically. Theunavailability of the exact analytical form for Jseems to be the biggest open question in this class of problems f29,31 g. The equilibrium surface shape of the liquid phase is a spherical cap s5d. This is a rigorous result valid during most of the drying process. However, when hs0,tdbecomes nega- tive and exceeds Hstdin its absolute value si.e., when thecenter point touches the substrate d, the surface shape is no longer spherical. Moreover, a new element of the contact lineis introduced in the center of the drop in addition to theoriginal contact line at the perimeter, and the entire evapora-tion profile gets modified in addition to the modified surfaceshape, thus influencing all the other quantities. Our treatmentdoes not account for the small fraction of the drying processoccurring after this touchdown swhich is a change in topol- ogy of the free surface, and thus requires a separate treatmentafter it has happened d. First of all, the amount of liquid re- maining in the drop at this moment is of the order of e compared to the original volume, and hence it would notmodify our main-order analytical results. Second, as our nu-merical calculations show, at touchdown practically all thesolute is already in the deposit phase, and the remainingamount of solute in the liquid phase is insignificant. Thus,within our model, the remainder of the drying process cannotmodify the deposit ring substantially, and hence this neglectof the late-time regime seems well justified. Experimentally,the inner part of the deposit ring is different from our pre-diction swhich is a vertical wall dand appears to have a spread shelf. The presence of this tail in the deposit distribu-tion can be caused by several features absent in our model.Its inherent two-dimensionality may be one of these short-comings sas discussed above d; the account of the dynamical processes occurring after the deposit phase has already beenformed se.g., avalanches of the inner wall dmay be another missing feature. The absence of a treatment of the late-timeregime may be among these reasons influencing the finaldistribution of the deposit as well. A more detailed accountof the effects of this late-time regime might be required inthe future. ACKNOWLEDGMENTS The author acknowledges indispensable help and advice of Thomas A. Witten and valuable input from Todd F. Du-pont and Robert R. Deegan. This work was supported in partby the National Science Foundation MRSEC Program underGrant No. DMR-0213745. APPENDIX: EVAPORATION RATE The purpose of this section is to obtain the evaporation rate from the free surface of a round sessile drop on thesubstrate. Since presence of the solute is irrelevant to thispurpose sat least for its low concentrations d, one may assume that the solute is simply absent and the drop is just purewater.We first consider the generic problem with an arbitrarycontact angle u, and then find the appropriate limit of interest u!1. If the radius of the drop footprint on the substrate is Ri, then its surface shape for an arbitrary contact angle is givenby hsr,td=˛Ri2 sin2ustd−r2−Ricotustd. sA1d Our results are presented here in their closed analytical form and some of them correct earlier expressions of this two-century-old problem.YURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-14 Our task involves solution of the equivalent electrostatic problem sthe Laplace equation dfor a conductor of the shape of the drop plus its reflection in the plane of the substrateskept at constant potential, as a boundary condition d.I nt h e case of the round drop the shape of this conductor resemblesa symmetrical double-convex lens comprised of two spheri-cal caps. The orthogonal coordinates that match the symme-try of this object sso that one of the coordinate surfaces co- incides with the surface of the lens dare called the toroidal coordinates s a,b,fd, where the coordinates aand bare related to the cylindrical coordinates randzby r=Risinh a cosh a−cos b,z=Risinb cosh a−cos b, sA2d and the azimuthal angle fhas the same meaning as in the cylindrical coordinates. Solution of the Laplace equation intoroidal coordinates involves Legendre functions of frac-tional degree and was derived in a book by Lebedev f25g. The electrostatic potential or vapor density is independent ofthe azimuthal angle fand reads nsa,bd=n‘+sns−n‘d˛2scosh a−cos bd 3E 0‘cosh utcoshs2p−bdt cosh ptcoshsp−udtP−1/2+itscosh addt. sA3d Herensis the density of the saturated vapor just above the liquid-air interface sor the potential of the conductor d,n‘is the ambient vapor density sor the value of the potential at infinity d, andP−1/2+itsxdare the Legendre functions of the first kind sthey are real valued d. The surface of the lens is described by the two coordinate surfaces b1=p−uandb2 =p+u, and the bderivative is normal to the surface. The evaporation rate from the surface of the drop is thereforegiven by Js ad=D1 hb]bnsa,bdb=2p+b1 =Dcosh a−cos b Ri]bnsa,bdb=3p−u, sA4d whereDis the diffusion constant and hb=Ri/scosh a −cos bdis the metric coefficient in coordinate b.fNote that an incorrect expression for Jwith a plus sign in the metric coefficient was used in Eq. sA2dof Ref. f21g.gThus, an exact analytical expression for the absolute value of the evapora-tion rate as a function of ris available: Jsrd=Dsn s−n‘d RiF1 2sinu+˛2scosh a+cos ud3/2 3E 0‘cosh ut cosh pttanhfsp−udtgP−1/2+itscosh adtdtG, sA5d where the toroidal coordinate ais uniquely related to the polar coordinate ron the surface of the drop:r=Risinh a cosh a+cos u. sA6d The expression sA5dis valid for an arbitrary contact angle u and corrects an earlier expression of Ref. f27gfEq.s28dg where a factor of ˛2 in the second term inside the large square brackets is missing. The expression for the evaporation rate is not operable analytically in most cases, as it represents an integral of anontrivial special function swhich, in its turn, is an integral of some simpler elementary functions d. In most cases, it is necessary to have recourse to asymptotic expansions in thecontact angle uin order to obtain any meaningful analytical expressions. However, there is one exception to this generalstatement. An important quantity is the total rate of watermass loss by evaporation dM/dt, which sets the time scale for all the processes. This total rate can be expressed as anintegral of the evaporation rate sdefined as the evaporative mass loss per unit surface area per unit time dover the surface of the drop: dM dt=−E AJsrd˛1+s„hd2rdrd f =−E 0Ri Jsrd˛1+s]rhd22prdr, sA7d where the first integration is over the substrate area Aoccu- pied by the drop. This expression actually involves tripleintegration: one in the expression above as an integral ofJsrd, another in the expression for Jsrdas an integral of the Legendre function of the first kind, and the third as an inte- gral representation of the Legendre function in terms of theelementary functions. However, it is possible to simplify theabove expression significantly and reduce the number of in-tegrations from three to one. Investing some technical effortand using Eq. s2.17.1.10 dof Ref. f33g, one can obtain a substantially simpler result that does not involve any specialfunctions at all: dM dt=−pRiDsns−n‘dFsinu 1+cos u +4E 0‘1+cosh2 ut sinh2 pttanhfsp−udtgdtG.sA8d This result together with the expression for the total mass of water M=rE 0Ri hsr,td2prdr=prRi3cos3u−3cos u+2 3sin3u sA9d swhere ris the water density dprovide a direct method for finding the time dependence of ufor anarbitrary value of the contact angle. Combining the time derivative of the lastexpression with result sA8d, one can obtain a single differen- tial equation for uas a function of time t:EVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-15 du dt=−Dsns−n‘d rRi2s1+cos ud2Fsinu 1+cos u +4E 0‘1+cosh2 ut sinh2 pttanhfsp−udtgdtG.sA10d Having determined the dependence ustdfrom this equation, one can obtain the time dependence of any other quantity dependent on the contact angle, for instance, the time depen-dence of the mass from relation sA9d, or any other geometri- cal quantity considered earlier. In practice, however, the analytical calculations in closed form cannot be conducted any further for arbitrary contactangles, and we will use the limit of small contact angles inall the subsequent analytical calculations. Besides being thelimit of our interest and most practical importance, this limitisalsoperfectlyadequateevenforquitesubstantialangles,aswill be seen in a moment. Expanding the right-hand side of Eq. sA10din small u,w e immediately obtain that the contact angle decreases linearly with time in the main order of this expansion: u=uiS1−t tfD, sA11d where we introduced the total drying time tfdefined in terms of the initial contact angle ui=us0d: tf=prRi2ui 16Dsns−n‘d. sA12d In the main order, the total rate of water mass loss is constant and the water mass also decreases with time linearly: M=prRi3ui 4S1−t tfD. sA13d This linear time dependence during the vast majority of the drying process was directly confirmed in experimentsf23,24 g; see Fig. 7. The dependence of the evaporation rate sA8don radius slinearity in R idwas also confirmed experi- mentally and is known to hold true for the case of diffusion-limited evaporation f34g. In Fig. 14, we plot the exact numerical solution for Mstd based on Eqs. sA9dandsA10dfor several values of the initial contact angle uitogether with the small-angle asymptotic of Eq.sA13d. In this figure, Miis the initial mass of water in the drop defined by the prefactor in Eq. sA13d.fNote that tfis not the total drying time for each ui; instead, it is just the combination of the problem parameters defined in Eq. sA12d, which coincides with the total drying time only when ui!0.gFigure 14 demonstrates that the small-angle approxi- mation works amazingly well up to angles as large as 45°,and therefore no precision or generality is lost by working inthe limit of small contact angles for the typical experimentalvalues of ui. Lastly, we note that the large-angle corrections may be responsible for the observed nonlinearity of theexperimentally measured dependence Mstd, as is clear from the comparison of Fig. 14 stheory dand Fig. 7 sexperiment d. The expression for the evaporation rate sA5dbecomes particularly simple in the limit of small contact angles. Em-ploying one of the integral representations of the Legendrefunction in terms of the elementary functions fEq.s7.4.7 dof Ref. f25gg, it is relatively straightforward to obtain the fol- lowing result: Jsrd=Dsn s−n‘d Ri2 pcosha 2su!0d, sA14d which, upon the identification cosh a=sRi2+r2d/sRi2−r2dfor u=0, can be further reduced to Eq. s11d. Thus, for thin drops the expression for the evaporation rate reduces to an ex-tremely simple result featuring the reciprocal square-root di-vergence near the edge of the drop. The same result couldhave been obtained directly if we solved an equivalent elec-trostatic problem for an infinitely thin disk instead of thedouble-convex lens. It is particularly rewarding that after allthe laborious calculations the asymptotic of our result is inexact agreement with the predictions of a textbook ssee Ref. f35gfor the derivation of the reciprocal square-root diver- gence of the electric field near the edge of a conducting planein the three-dimensional space d. Equation s11dis the result we were looking for in our case of the thin circular drops. For the sake of completeness, it is also interesting to note the opposite limit of the expression sA5d, when the surface of the drop is a hemisphere s u=p/2d. In this limit, a similar calculation can be conducted, and the uniform evaporation rate is recovered: Jsrd=Dsns−n‘d Risu!p/2d. sA15d This result is also in perfect agreement with the expectations; the same result could have been obtained if we directlysolved the Laplace equation for a sphere sthe hemispherical drop and its reflection in the substrate d. The uniform evapo- ration rate is a result of the full spherical symmetry of such asystem. Similar exact results can also be obtained for a fewother discrete values of the contact angle se.g., for u=p/4d. FIG. 14. Numerical results: dependence of water mass Mon timet. Different curves correspond to different initial contact angles; values of parameter uiare shown at each curve. The ana- lytical result fEq.sA13dgin the limit ui!0 is also provided sthe solid curve d.YURI O. POPOV PHYSICAL REVIEW E 71, 036313 s2005 d 036313-16 f1gN. D. Denkov, O. D. Velev, P. A. Kralchevsky, I. B. Ivanov, H. Yoshimura, and K. Nagayama, Langmuir 8, 3183 s1992 d. f2gA. S. Dimitrov, C. D. Dushkin, H. Yoshimura, and K. Na- gayama, Langmuir 10, 432 s1994 d. f3gC. D. Dushkin, H. Yoshimura, and K. Nagayama, Chem. Phys. Lett.204, 455 s1993 d. f4gT. Ondarcuhu and C. Joachim, Europhys. Lett. 42, 215 s1998 d. f5gJ. Boneberg, F. Burmeister, C. Shafle, P. Leiderer, D. Reim, A. Fery, and S. Herminghaus, Langmuir 13, 7080 s1997 d. f6gR. G. Larson, T. T. Perkins, D. E. Smith, and S. Chu, Phys. Rev. E55, 1794 s1997 d. f7gJ. P. Jing, J. Reed, J. Huang, X. Hu, V. Clarke, J. Edington, D. Housman, T. S. Anantharaman, E. J. Huff, B. Mishra, B. Por-ter, A. Shenkeer, E. Wolfson, C. Hiort, R. Kantor, C. Aston,and D. C. Schwartz, Proc. Natl. Acad. Sci. U.S.A. 95, 8046 s1998 d. f8gA. B. El Bediwi, W. J. Kulnis, Y. Luo, D. Woodland, and W. N. Unertl, in Hollow and Solid Spheres and Microspheres— Science and Technology Associated with Their Fabricationand Application , edited by M. Berg, T. Bernat, D. L. Wilcox, Sr., J. K. Cochran, Jr., and D. Kellerman, Mater. Res. Soc.Symp. Proc. No. 372 sMaterials Research Society, Pittsburgh, 1995 d, p. 277. f9gF. Parisse and C. Allain, J. Phys. II 6, 1111 s1996 d. f10gF. Parisse and C. Allain, Langmuir 13, 3598 s1996 d. f11gE. Adachi, A. S. Dimitrov, and K. Nagayama, in Film Forma- tion in Waterborne Coatings , edited by T. Provder, M. A. Win- nik, and M. W. Urban sAmerican Chemical Society, Washing- ton, D.C., 1996 d, p. 419. f12gE. Adachi, A. S. Dimitrov, and K. Nagayama, Langmuir 11, 1057 s1995 d. f13gL. Shmuylovich, A. Q. Shen, and H. A. Stone, Langmuir 18, 3441 s2002 d. f14gJ. Conway, H. Korns, and M. R. Fisch, Langmuir 13, 426 s1997 d. f15gK. S. Birdi, D. T. Vu, and A. Winter, J. Phys. Chem. 93, 3702 s1989 d.f16gK. S. Birdi and D. T. Vu, J. Adhes. Sci. Technol. 7, 485 s1993 d. f17gM. E. R. Shanahan and C. Bourges, Int. J. Adhes. Adhes. 14, 201s1994 d. f18gC. Bourges and M. E. R. Shanahan, Langmuir 11, 2820 s1995 d. f19gS. M. Rowan, M. I. Newton, and G. McHale, J. Phys. Chem. 99, 13268 s1995 d. f20gS. M. Rowan, G. McHale, M. I. Newton, and M. Toorneman, J. Phys. Chem. B 101, 1265 s1997 d. f21gR. D. Deegan, O. Bakajin, T. F. Dupont, G. Huber, S. R. Na- gel, and T. A. Witten, Phys. Rev. E 62, 756 s2000 d. f22gR. D. Deegan, O. Bakajin, T. F. Dupont, G. Huber, S. R. Na- gel, and T. A. Witten, Nature sLondon d389, 827 s1997 d. f23gR. D. Deegan, Phys. Rev. E 61, 475 s2000 d. f24gR. D. Deegan, Ph.D. thesis, Department of Physics, University of Chicago, 1998. f25gN. N. Lebedev, Special Functions and Their Applications , re- vised English ed. sPrentice-Hall, Englewood Cliffs, NJ, 1965 d, Chaps. 7 and 8. f26gR. G. Picknett and R. Bexon, J. Colloid Interface Sci. 61, 336 s1977 d. f27gH. Hu and R. G. Larson, J. Phys. Chem. B 106, 1334 s2002 d. f28gY. O. Popov and T. A. Witten, Eur. Phys. J. E 6,2 1 1 s2001 d. f29gY. O. Popov and T. A. Witten, Phys. Rev. E 68, 036306 s2003 d. f30gT. Dupont sprivate communication d. f31gY. O. Popov, Ph.D. thesis, Department of Physics, University of Chicago, 2003, e-print cond-mat/0312196. f32gY. O. Popov, J. Colloid Interface Sci. 252, 320 s2002 d. f33gA. P. Prudnikov, Y. A. Brychkov, and O. I. Marichev, Integrals and Series sGordon and Breach, London, 1986 d, Vol. 3. f34gJ. T. Davies and E. K. Rideal, Interfacial Phenomenon sAca- demic Press, New York, 1963 d. f35gJ. D. Jackson, Classical Electrodynamics , 2nd ed. sWiley, New York, 1975 d, Chaps. 2 and 3.EVAPORATIVE DEPOSITION PATTERNS: … PHYSICAL REVIEW E 71, 036313 s2005 d 036313-17