• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Non-Stokes drag coefficient in single-particle electrophoresis:New insights on a classical problem

    2019-08-16 01:17:46MaiJiaLiao廖麥嘉MingTzoWei魏名佐ShiXinXu徐士鑫DanielOuYang歐陽新喬andPingSheng沈平
    Chinese Physics B 2019年8期

    Mai-Jia Liao(廖麥嘉), Ming-Tzo Wei(魏名佐), Shi-Xin Xu(徐士鑫),H Daniel Ou-Yang(歐陽新喬), and Ping Sheng(沈平)

    1Department of Physics,Hong Kong University of Science and Technology,Clear Water Bay,Kowloon,Hong Kong,China

    2Department of Physics,Lehigh University,Bethlehem,PA,USA

    3Department of Bioengineering,Lehigh University,Bethlehem,PA,USA

    4Centre for Quantitative Analysis and Modelling(CQAM),The Fields Institute,Toronto,Ontario M5T 3J1,Canada

    Keywords: electrophoresis,drag coefficient,vortices belt,charged colloidal particle

    1. Introduction

    When immersed in an electrolyte solution, a charged particle would be enveloped in an ionic cloud of screening counter-ions, denoted the Debye layer. Application of an external electric field to a suspension of such charged particles can result in the steady motion of the solid particulates.The physical picture underlying this phenomenon, known as electrophoresis,[1-5]dates back to Smoluchowski[6]in which the crucial element is the electroosmotic fluid flow in the Debye layer. Through clever mathematical manipulations,Smoluchowski has shown rigorously that electrophoretic mobility of the charged particle,μE,is directly proportional to the zeta potential ζ (which is directly related to the surface charge density)on the surface of the solid particle,i.e.,v∞=μEE∞,where E∞is the applied electric field,v∞is the electrophoretic velocity, and μE= -(ε/η)ζ, with η, ε being the solution viscosity and dielectric constant, respectively. The Smoluchowski relation is accurate in the limit of a ? λD,where λDis the Debye length and a the particle radius. A derivation of the Smoluchowski relation is given in Section S2 in Supplementary material.

    There have been extensive theoretical[6-14]and experimental[15-20]studies of a charged particle under the simultaneous effect of an electric field E∞and a non-electric force Fext. In particular, the electrophoretic drag coefficient γEis of interest. This is because not only the Smoluchowski electrophoretic flow field(see Section S2 in Supplementary material) differs significantly from the Stokes flow field, but also the lack of an accurate flow field solution inside the Debye layer prevents an accurate account of the actual hydrodynamic drag force on the solid surface of the charged particle.A rough evaluation of drag force on the particle surface may be established through a simple scaling argument.[14,21]In the thin-Debye layer electrophoresis, the velocity gradient in the liquid is screened beyond the λDscale, thus the drag force is(4πa2)η/v∞. As a/1 in most cases, it follows that(=6πηa, the Stokes drag coefficient), with the ratioreaching as high as 100 to 200 in some cases.However,our accurate measurements of γE,described below,does not come anywhere close to such order of magnitude deviation from the Stokes drag. This large discrepency necessitates an understanding of the physical picture underlying the electrophoretic drag.

    Based on the linearization of coupled electrohydrodynamic equations, plus the superposition of external nonelectric force Fextand the electrophoresis problem,[1,21]it is found that a charged particle behaves similarly in an electric field as in a hydrodynamic flow.[22]From force balance one obtains:where v denotes the solid particle velocity under the combined electric and non-electric (mechanical) forces. Equation (1)suggests that the mechanical force needed to stall an electrophoretic motion is a Stokes-like drag force.

    In anticipation of subsequent developments,we make two remarks in relation to Eq.(1). The first is that equation(1)can be expressed in terms of two equations by the velocity decomposition v=v∞+Δv, where the first equation, v∞=μEE∞,is the definition of the electrophoretic velocity. The second equation is then Fext-γsΔv =0. Such a decomposition is useful in demonstrating that the electrophoretic drag coefficient γEis not necessarily identical to γS. This can be easily seen by multiplying the velocity equation on both sides by γE:FE=γEv∞=γEμEE∞=QeffE∞. Here γEis simply defined to be the coefficient of proportionality between velocity and hydrodynamic viscous force, in this case both resulting from an applied electric field,and Qeff=γEμEis introduced to distinguish it from the surface charge, since Qeffis known to be much smaller than the surface charge[14,16]and represents,in the context of force balance,the solid particle’s coupling to the applied electric field.

    The second remark is related to our experimental approach of using an optical trap to hold a single charged particle in a harmonic potential and applying an AC electric field to induce periodic oscillations of the particle. The accurately measured quantities are then the amplitude of particle’s periodic motion and its phase difference with the applied AC electric field. Owing to the AC nature of the applied electric force,there are inevitably the in-phase and out-of-phase components of the force relative to the particle velocity. For Eq. (1), i.e.,Fext-γsΔv =0, the applied Fextis necessarily the in-phase component of the force, it describes that the force needed to alter, or stall, an intrinsic electrophoretic motion is Stokes by nature as argued by Long et al.[21]and thus,the equation itself has nothing to do with electrophoresis. Our simulation results(see below)supports the correctness of Long’s argument, but contradicts the non-Stokes prediction of the in-phase external force by Lizana et al.[14]We will show below in the subsequent section that in our experiment the optical trapping force is always out-of-phase relative to the phase of the velocity,whereas the applied electric force has both an in-phase component and an out-of-phase component.The in-phase component of the electric force drives the particle velocity and the out-ofphase component automatically counter-balances the out-ofphase optical trapping force. These facts account for our approach’s ability to measure the electrophoretic drag force and its drag coefficient.

    The drag coefficient is always related to the hydrodynamic drag force exerted on a surface. As Smoluchowski has successfully linked μEto the surface zeta potential, and measured drag force is given byγEv∞=γEμEE∞=QeffE∞,achieving force balance between the electrical force and drag force at the solid surface(by equating Qeffto QS)would seem to be the most convenient choice. However, the fact that Qeffis known to be much smaller than QS,[14,16]which is also confirmed by our measurements as seen below,indicates that there must be another surface,away from the solid surface,[23-26]on which the drag force and electrical force attain force balance.Two questions naturally arise: (i) How does such a surface emerge consistently from the relevant mathematical equations governing the electrophoresis, and (ii) Can one measure an electrophoretic drag coefficient γEthat is consistent with the theory prediction on such a surface?

    To address these two questions, we set out to measure γEand Qeffdirectly by experiments, and to obtain mathematically accurate simulations of the electrophoretic flow field.From experimental measurements the mobility was obtained as μE= Qeff/γE. It is shown that while γE/γSvaries only from 1 to 3 over a range of salt concentrations,the magnitude of Qeffis smaller than that of QSby many orders of magnitude. The resulting μE, however, agrees well with the literature values on similar micro-spheres, as well as that obtained by DC measurements in our experiments.[26]To provide a microscopic explanation of our experimental data, we have carried out simulation of the electrophoresis effect by numerically solving the Poisson-Nernst-Planck (PNP) equations coupled with the Navier-Stokes(NS)equation, implemented with the appropriate boundary conditions. Using only the experimentally measured μEas input, the most significant outcome is that the measured γEand Qeffare obtained at a uniquely distinctive interface between the inner and outer flow fields that is at a distance a few λD’s away from the solid surface. While the electrophoretic flow in the outer flow region follows the scaling of the Smoluchowski solution,i.e.,r-3,our numerical solution reveals an inner flow field that is carried along by the solid particle, with a highly nonlinear electro-hydrodynamic flow behavior. The interface between the inner and outer flow fields acts as the reference surface, or slip surface/plane, at which both the magnitude and trend(e.g.,with respect to salt concentration)of the drag coefficient and effective charge can be accounted for,and from which the deviations can be further explored.

    In what follows,we first present the experimental results,followed by simulations and discussion of the physical picture that emerges. We end by underscoring the closure between theory and experiment.

    2. Materials and methods

    An optical trap was effected by an IR (wavelength =1064 nm) laser coupled into an oil-immersion objective lens(100X, NA = 1.3, Olympus). A second IR laser beam(wavelength=980 nm),aligned and focused by the same objective lens to be par-focal with the trapping laser focus, was used for particle tracking. A schematic diagram of the measurement apparatus is shown in Fig.S1 in Supplementary ma-terial. Here, we use the optical trap both as a tool to control/monitor the position of the particle as well as a calibrated force sensor[27-29]with a stiffness constant ktrap= 17.9±0.1 pN/μm. Experimental details are given in Section S3 in Supplementary material.

    Polystyrene(PS,Thermo Scientific,catalog#5153A)particles with a mean radius a=0.75μm were dispersed at volume fractions below 1×10-4in solutions of varying concentrations of potassium chloride electrolyte with a Debye length λDranging from 96 nm to 4.3 nm(KCl concentration 0.01 mM to 5 mM, κa(≡a/λD) ranging from 7.8 to 174). To reduce the effect of dissolution of CO2in DI water, all electrolytes were prepared with fresh DI water filtered with resin. The colloidal solutions were inserted in a glass chamber with vertical height ~200 μm. To avoid the effects of the Fax′en drag on the particle[30-32]and the flow due to any surface electroosmosis flow near the glass surface,[33]we kept the colloid sphere 8μm-15μm above the lower glass plate.

    The electrophoretic chambers were fabricated by coating a glass microscope slide substrate with Sigmacote (Sigma-Aldrich),a solution of chlorinated organopolysiloxane in heptane that readily forms a covalent,microscopically thin film on glass, with a purpose to suppress surface charge and to minimize surface electroosmotic flow. Subsequently gold wires for applying the electric field were placed on the glass substrate with a separation of ~1 cm. An aqueous solution of dilute particles was then dipped onto the substrate,and a Sigmacotecoated cover-glass was used to seal the chamber with wax sealant. External electrical wires were welded with the gold wires. Electric field was measured by inserting a second pair of parallel gold wires, separated by 5 mm, into the sample to obtain the voltage drop between them. The optical path diagram is shown schematically in Fig.S1 in Supplementary material. The electrophoresis chamber containing the particles was mounted on an inverted microscope(Olympus IX81). To minimize the contribution to the optical trapping effect, the 980-nm tracking beam power was two orders of magnitude lower than that of the 1064-nm laser. Movements of the particle, tracked by the 980 nm laser beam, were detected by a quadrant photodiode (QPD, S7479, Hamamatsu). The voltage reading of the QPD was maintained to be within the linear range of the particle displacement from the trapping center.A lock-in amplifier (SR830, Stanford Research) was used to record the phase of particle motion relative to that of the sinusoidal voltage applied to PZT to form an oscillatory optical tweezers. The wide-field images of the particle were captured by a CCD camera for the purpose of optical alignment.

    3. Experimental results

    3.1. Optical trap and AC electric field

    We apply an AC electric field(20 Hz-100 Hz)to a single spherical particle held by a calibrated optical trap as shown in Fig. 1(a). At such frequencies the period of the applied field was much smaller than the time required to screen the electrodes(>1 s for a separation between the electrodes=1 cm),and much larger than the relaxation time of the electrical double layer[34](<1μs).The relevant equation of motion is given by:

    Fig.1. Electrophoretic measurements using the optical tweezers. (a)Schematic illustration of a charged colloidal particle held by an optical tweezer and driven by an oscillating electric field. The electric field direction is indicated for a particular instant of time. The measured γE/γS(b),magnitude of Qeff(c),andμE(d)plotted as a function of applied electric field frequency for a 0.75-μm radius polystyrene particle,held by a fixed optical trap in DI water(λD=96.1 nm)under two different electric field strengths. It is seen that all measured quantities are relatively independent of the frequency,indicating the quasi-static nature of the experiment.

    Here, m is the mass of the particle, x is the displacement of the particle from the center of the harmonic optical trap, and γ is the drag coefficient. Here the drag coefficient is purposely written without a subscript, since in Eq. (2) there are two forces - optical trap force and electric force - acting on the particle and hence the specification of which drag coefficient applies would depend on the situation to be analyzed below. The Brownian noise term is excluded in the above equation because its contributions will be filtered out by the phase-sensitive, lock-in detection technique employed in our experiments. Here Qeffis the effective charge that provides particle’s electrical coupling to the applied electric field. It should be noted that the deformation of the counterion cloud is very small,so that the resulting electric“retardation”force can be neglected.[35]Our simulations,based on the full numerical solution of the relevant governing equations, have confirmed the negligible effect of retardation.

    By approximating the left-hand side of Eq.(1)to be zero,which is accurate considering the fact that mω2is 5 orders of magnitude lower than ktrap,we have:

    where E∞=E(0)∞exp(-iωt).We write the displacement of the electrophoretic particle as x(t)=D(ω)e-i[ωt-δ(ω)], with the displacement amplitude denoted by D(ω)and the phase δ(ω)defined relative to the applied AC electric field. In our experiments,the phase shift was measured with a lock-in amplifier with an accuracy of a couple of degrees. This is in contrast to the phase shift measurements in a previous paper[15]that has errors in the range of+/-11 degrees. Through the precise measurement of the phase shift, our method enabled the extract of drag coefficient value in an accurate and robust manner.

    3.2. An analysis of the AC experimental approach

    Substituting particle displacement expression for x(t)and the associated displacement velocity ˙x(t) back into Eq. (3),and taking the real part of every term,we get

    Let us define a t0such that ωt0=δ(ω)+π/2, and let ωt =ω(t0+t′)=ωt0+ωt′. Equation (4) can be re-written in an illuminating form as

    At t′=0,the left-hand side vanishes,i.e.,optical trap force is zero,and we only have the right-hand side,from which we obtain

    We label the drag coefficient as γEbecause there is only the electric field force present in Eq.(6). Now let us expand

    in Eq.(5). Then equation(5)can be re-organized in a physically clear manner as

    In Eq.(7),the right-and left-hand sides represent the in-and out-of-phase components, respectively. There is no approximation made in Eq. (7), thus it has to be valid for all values of t′. The only way equation(7)can be true is that both sides must separately be zero, since the time variations on the two sides are orthogonal to each other.

    Physically, equation (7) states that the electrical force comprises two components. One component,E(0)∞Qeffsinδ(ω)cos(ωt′),gives rise to the time-varying electrophoretic velocity ωD(ω)cos(ωt′)with a constant drag coefficient γE. The other component,E(0)∞Qeffcosδ(ω)sin(ωt′),counter-balances the optical trap force ktrapD(ω)sin(ωt′).This optical trap response force is noted to be out of phase with the electrophoretic velocity, in contrast to the external force in Long’s equation(Eq.(1)),which is in-phase with the electrophoretic velocity.

    Equation(7)essentially expresses the fact that 0=0,from which we obtain two independent equations:

    Equation(8b)is identical to Eq.(6). And if we divide Eq.(8b)by Eq. (8a), we obtain the consistency condition imposed by the force balance of the in- and out-of-phase components of the AC experiment:

    Equation(8a)can be re-written for Qeffin terms of the physically measured quantities as

    From Eqs.(9a)and(9b)the mobility is obtained as

    Not limited by the constraint imposed by typical DC measurements that any external non-electrical force would alter the intrinsic electrophoretic motion and yield a Stokes-like drag,the most essential point of the above analysis is that the external force imposed by the optical trap is always 90 degrees out of phase with the electrophoretic velocity. Through this feature of the AC electric field-driven particle in an optical harmonic trap,we can measure electrophoresis drag coefficient γEaccurately without the influences of an in-phase,non-electrical external force.

    3.3. Measured drag coefficient, effective charge and mobility

    Using the measured values D(ω)and δ(ω)in Eqs.9(a)and 9(b), respectively, we have evaluated γEand Qeff,which are plotted in Figs. 1(b) and 1(c) as a function of applied electric field frequency for a 0.75-μm radius polystyrene particle, held by a fixed optical trap in de-ionized (DI) water(λD= 96.1 nm) under two different electric field strengths.In Fig. 1(d) we show the measured values of μEplotted as a function of applied electric field frequency, also for two different electric field strengths. It is seen that the measured γE,Qeff,andμEwere relatively independent of the frequency and the field strength,indicating the quasi-static nature of our measurements, we take as measured values the averages over the measured frequency range.

    It should be noted here that although the electrophoretic mobility determined from Eq. (9c) was determined from the AC measurements, we have also measured it using the traditional approach by using a DC field. Experimentally,this was done by turning off the optical trap and measuring the speed of the same particle in the presence of a DC electric field. The mobility determined by the AC measurements agreed with that by the DC measurements, as expected. It should be noted that the more precise phase shift values in our measurements clearly show the phoretic drag to be non-Stokes, while simultaneously the electrophoretic mobilities obtained by our AC and DC measurements agree within the experimental error. This is in contrast to the earlier experiment by Semenov et al.,[15]in which the phase shift measurement had a much larger error bar. Moreover, by assuming the drag coefficient to be Stokes,Semenov’s study yielded different mobilities between their AC and DC experiments.

    In Fig. 2(a) the measured γE/γSis plotted as a function of κa with open symbols,where κ =1/λD. Non-Stokes drag coefficients values are clearly observed. As κa increases towards 110, however, γEgradually approaches γSwhich is in direct contrast to the 4πηa2/λDvalue of the scaling argument in the thin Debye layer limit.[14,21]

    Fig.2.Measured drag coefficient and effective charge plotted as a function of salt concentration, expressed here as κa, where κ =1/λD. (a)Open symbols denote the values of measured γE/γS plotted as a function of κa. Here γS denotes the calculated Stokes drag coefficient,and γE is obtained by using measured δ in Eq.(9). Error bars on the open circles arise from the errors in amplitude and phase measurements. (b)Measured effective charge magnitude |Qeff| is plotted as a function of κa with open symbols. Filled symbols indicate the surface charge magnitude |QS| obtained from simulations with measured mobility μE as the only input from the experiment. Error bars are from the mobility measurements. All results are for PS spheres with a=0.75μm.

    Values of the measured magnitude of the effective charge|Qeff| are plotted with open symbols as a function of ionic strength in Fig. 2(b). They are in the range of 1700 to 4800 electronic charges, which are one to two orders of magnitude smaller than the generally accepted values deduced from known PS surface charge density[26,36,37]that is in the range of ~-6×104e/μm2. For a 0.75-μm radius PS spheres, one obtains |QS|~1.9×105e, where e denotes the magnitude of the electronic charge. In Fig. 2(b), the magnitude of the surface charge|QS|obtained from our numerical simulations,ob-tained with the experimentally measured mobility values as inputs,are shown by solid symbols,these values are noted to be in the expected range of|QS|for PS surfaces, thus further confirming the huge differences with the magnitude of the effective charge. Also, as QSis always directly related to the mobility,the resulting zeta potential/mobility relation is compared to the literature in Section S5 in Supplementary material.Excellent agreement is obtained.

    4. Simulations

    4.1. Equations and boundary conditions

    Our starting point is the coupled incompressible NS equation and the PNP equations for an electrolyte that is symmetric between the positive and negative ions. The equations are as follows:

    Here ψ stands for electrostatic potential, p(n) stands for the positive (negative) ion density, Jp(n)stands for ionic current density, with the diffusive, electric convective, and flow convective components,η =1×10-3Pa·s is the liquid viscosity,taken to be that of water, ρ?the liquid density of water, dielectric constant of water, ε0εf, is taken to be εf=80, with ε0=8.85×10-12F/m,u is the liquid velocity,and P denotes pressure. The valence of all ions is taken to be one. Here we use the value of Dion=9.32×10-9m2/s.[38]The boundary conditions for Eqs. 10(a)-10(f) are ψ(r)→-E∞·r as r →∞Jn,p|r=a=0, where ?n stands for the unit outward surface normal, and E∞is the applied electric field. Here infinity denotes the simulation domain boundary. In Eqs.(10g)and (10h), E =-?ψ and Γ stands for the particle surface,FE(H)stands for the electric (hydrodynamic viscous) force,and TE(H)denotes the electric(hydrodynamic viscous)stress tensor. Equation (10i) is the Newton’s equation for particle motion.

    Static simulations were performed in the co-moving particle frame with time derivatives in Eqs.(10a)-(10i)set to zero.Constant surface charge boundary condition was used for the boundary condition on the solid particle surface,with its value adjusted according to the criterion that FE+FH=0 on the sphere surface at steady state, subject to the velocity constraint u(r)|r→∞=-μEE∞at the simulation domain boundary, where μEis the measured mobility. The ion densities at simulation domain boundary are given by n(r)|r=∞=n∞,p(r)|r=∞= p∞,with equality between the two as required by overall charge neutrality. The boundary condition at the fluid particle interface is non-slip. The pressure at the simulation boundary is set as constant.For dynamic simulations as shown below, the initial velocity of the particle is set to be zero. A time varying electric field is applied to drive the particle motion. More details on static and dynamic simulations are presented in Section S5 in Supplementary material.

    4.2. Physical picture and discussion

    We would like to describe the salient features of the electrophoretic flow pattern,and to contrast them with those of the Stokes flow field. For a spherical particle acted on by an external force along the direction of unit vectorthe Stokes flow field in the lab frame[39]is given by

    where u0denotes the particle velocity,r denotes the radial coordinate,and θ the polar angle defined relativeIn contrast, the Smoluchowski fluid velocity field for particle electrophoresis in the lab frame can be written as(see Section S2 in Supplementary material):

    Fig. 3. Projected velocity and streamlines in the laboratory frame for same-sized particles under stationery condition (force balance with no acceleration). The white bars in the figures denote the length scale of 750 nm. (a)P2 projection of the electrophoretic velocity component along the electric field direction normalized by electrophoresis velocity v∞=μEE shows very good far-field 1/r3 behavior of the simulated velocity field,indicated by the red dashed line. Yellow dashed line shows the position of particle surface. Blue line shows the P2 projection of the Smoluchowski potential flow velocity field expressed by Eq.(11b). It is clear that the projected 1/r3 behavior stops at a small distance away from the solid surface.The peak position of the project result is delineated by the black dashed line,corresponding to the white dashed curve in panel(b). (b)Streamlines for the negatively charged PS sphere with a=0.75μm,(with λD=96.1 nm,κa=7.8,σ =7000 e/μm2)plotted in the lab frame. Here the colors indicate the magnitude of the velocity, the interface between the inner and outer flow fields (the reference surface) is denoted by white dashed curve. (c)Streamlines for a particle of the same radius as(a),acted on by a constant external body force to moves at same velocity asμEE∞. This represents the Stokes flow field in the lab frame. The contrast with the electrophoretic flow field shown in panel(b)is clear. (d)Streamline for the flow field obtained from a charged particle under the same electric field strength as in panel(b)but fixed by a reverse external(non-electrical)force acting on the particle. This flow field can also be obtained approximately by the superposition of the electrophoretic flow field in panel(b)and the Stokes flow field in panel(c)with a reversed velocity. It should be noted that the external(non-electrical)force required to immobilize the particle is not the same as the electrical force,thereby explaining the difference in the drag coefficients. The difference in the two forces accounts for the non-zero flow field close to the particle.

    In Fig. 3(b) we show the full simulated electrophoretic flow field in the laboratory frame. It is noted that the streamlines of the electrophoretic flow pattern in Fig. 3(b) exhibits a belt of stationary vortices around the equatorial plane of the particle,with the electric field direction defining the north and south poles. The interface between the inner and outer flow fields is delineated here by the white dashed line. In contrast,in Fig.3(c)we show the Stokes flow field with the same-sized particle moving at the same velocity as in panel(b). There are no vortices in the flow field. The difference in the magnitude of the far field velocity is also apparent,as indicated by color,owing to the much slower 1/r velocity decay of the Stokes flow. Since the Reynolds number is negligibly small in electrophoretic flow,here the existence of vortices seen in Fig.3(b)is due to the large net charge and the associated strong local electric field in the Debye layer. It is known mathematically that the conditions for the existence of vortices are(i)the existence of a source term for vorticity (ω=?×v) as derived from the NS equation, and (ii) velocity reversal since the velocity on two sides of a vortex must be opposite. The latter condition is guaranteed in the laboratory frame of the electrophoretic flow close to the solid boundary and in the vicinity of the equatorial plane, because the electroosmotic flow induced by the net (positive) charge within the Debye layer is opposite to that of the solid particle.Hence both conditions are satisfied for the electrophoretic motion in the laboratory frame,but not in the co-moving frame due to the lack of velocity reversal. Hence vortices are not apparent for the electrophoretic motion in the co-moving frame. For the Stokes flow none of the two conditions is satisfied in the laboratory frame. Mathematical details for the non-inertial generation of vorticity,as in the present case,are shown in Section S6 in Supplementary material.

    The physical manifestations of the inner flow field,as afforded by the full numerical solution of the relevant governing equations,represent a realization of the Smoluchowski picture in which the electroosmotic flow near the particle surface is the dominant mechanism of electrophoresis. In the original Smoluchowski solution(see Section S2 in Supplementary material)the inner flow field region was compressed to infinitesimal thickness for the sake of mathematical simplicity. Here we restore the inner flow field region with all its rich electrohydrodynamic flow behavior.

    To further accentuate the difference between the electrophoretic and Stokes flows and their respective drag coefficients, in Fig. 3(d) we show the simulated flow field of a charged particle under the same applied electric field strength as in Fig.3(b)but fixed by an external non-electrical body force acting on the particle. This flow field can also be obtained approximately by the superposition of the electrophoretic flow field in Fig.3(b)and the Stokes flow field in Fig. 3(c) with a reversed velocity. The fact - that the flow field in panel 3(d) is not identically zero - emphasizes our point that the external non-electrical body force required to immobilize the particle is not the same as the electrical force.It also provides another, and perhaps more direct, explanation for why the drag coefficient in Eq.(1)is not the same as the electrophoretic drag coefficient. The difference between the Stokes versus the electrophoretic drag forces accounts for the non-zero flow field close to the particle in Fig.3(d), even though the net force is evaluated to be zero on the solid surface,as required by force balance on an immobile particle.

    5. Closure between theory and experiment

    5.1. Drag coefficient under differentκavalues

    In Fig.4 we summarize the overall results of this work.In Fig.4(a)the ratio between the measured drag coefficients and the Stokes drag coefficient is plotted as open squares onto the solid simulation curves with different κa values, each evaluated on a surface at distance r outside the solid particle’s surface. The vertical dashed line indicates the inner/outer interface position,whose intersections with the different curves yield the drag coefficients at that interface. The measured results cluster around, and are close to, the values at the inner/outer interface. Their values are significantly smaller than the simulated values at the solid surface. The solid squares are calculated from the scaling argument[(4πa2)η/λD]/γS.[14,21]It is not surprising that these values at the solid surface are much larger than γS, since a/λDis in the range of 7-200 in our experimental measurements.

    In Fig. 4(b) we plot the position of the inner/outer flow field interface,in terms of its dimensionless distance κΔ away from the particle surface, as a function of the salt concentration. Solid symbols indicate the inner/outer interface position as determined by the peak values in U/v∞as shown in Fig. 3(a). The open symbols represent the values of Δ determined by using measured drag coefficients as input to find the spherical surfaces on which the calculated drag forces exactly matches the measured values. Both the magnitudes and the trend of the empty symbols track the inner/outer interface position very well. Note that Δ is independent of surface charge density under a given particle radius.

    The fact that the experimentally measured drag coefficient yields a surface position whose salt concentration dependence is identical to that of the inner/outer interface, ties the latter rather uniquely to its role as the reference surface for the electrophoretic drag and effective charge. Below we further reinforce this point by linking this inner/outer interface to the theoretical framework of O’Brien and White,[9]and Long et al.[21]

    5.2. Dynamic simulation and comparison with measured data

    To further substantiate the physical picture that the measured hydrodynamic drag coefficient should be evaluated by using the inner/outer interface (dashed white curve in Fig.3(b))as the reference(slip)surface/plane,we have carried out dynamic simulations by using the moving mesh approach.Details of the simulations are given in Section S6 in Supplementary material. In Fig.5 the open circles are the results of the dynamic simulation plotted as a function of time. The applied electric field is a step function with a very short rising time, shown in the inset. The dashed red curve is the fitting by using Eq.(12)below. In the absence of an optical trap,the dynamic equation of motion can be written as:

    Fig.5. Momentum time derivative evaluated inside the inner flow field(including the solid particle) plotted as a function of time. Black circles are from dynamic simulation and red dashed curve is fitted from Eq.(12). Inset shows the corresponding electric field as a function of time.It shows a linear increase within 5×10-8 s to reach the saturation level.

    For the left-hand side of Eq.(12),we sum up the momenta of the fluid and the solid at each finite element mesh within the inner/outer interface at each time step, and then evaluate its time derivative by using the data from successive time steps.The evaluation of the momenta is necessary since there are relative motions between the fluid and the solid. On the right hand side,the applied electric field is known,and for ˙x(t)we use the center of mass velocity of the evaluated unit. The two parameters γEand Qeffare then varied as the fitting parameters to best reproduce the simulation data. For a=0.75 μm and λD=96.1 nm, the best fit yields values of γE/γS=1.7 and |Qeff|=9.3×10-16C, they are in excellent agreement with the experimentally measured values of 1.69±0.27 and 9.1±1.4×10-16C, respectively. This indicates a consistent picture of drag coefficients evaluated from dynamic and static simulations,as well as consistency with the experimental measurements.

    The fact that the inner flow field is carried along by the charged particle also means that the interface that separates the inner and outer flow fields represents a slip surface/plane.

    5.3. Consistency between the framework of O’Brien and White, and Long et al., with the present theory and experimental results

    We return to Eq. (1) and show that the inner flow field and its well-defined interface with the outer flow field constitute the missing piece in the theoretical framework of O’Brien and White,[22]and Long et al.,[21]that can tie together our theoretical predictions and the experimental results.

    In a seminal paper, O’Brien and White simplified the governing equations for electrophoresis by assuming the coupled hydrodynamic and electric forces on a particle undergoing electrophoresis can be decoupled and solved separately.[22]The problem therefore becomes a superposition of two problems: the pure force problem (particle held fixed in a flow field with no applied field),and the pure electrophoresis problem (particle in an electric field in an electrolyte which is at rest at points far from the particle). Based on this analysis,Long et al.[21]showed that the velocity of the sphere is the sum of the two problems described above, based on global force balance:[14,21]

    Here Felectricstands for the electric force, Fviscdenotes the viscous drag force, and Fextdenotes the non-electric external force. It should be reminded that v stands for the relative velocity between the particle and the far field fluid. Here equation(13a)is noted to be the same as Eq.(1).

    Since we have demonstrated by using dynamic simulation that the inner flow field is carried along by the charged sphere, here we wish to use the inner/outer interface to evaluate all the relevant forces,so as to test Long et al.’s analysis while simultaneously also demonstrate the fact that γE/=γS.The latter can be seen from Eq. (13a) that, if Fext=0, then we must have v=μEE∞. When that happens, the force balance, equation (13b), becomes Feletric+Fvisc=0, i.e., a test of the electrophoretic drag also becomes possible since we can write Feletric=QeffE and Fvisc=γEv,both evaluated at the inner/outer interface.

    All simulations were performed with surface charge density σ =-7×103e/μm2, particle radius a=0.75 μm, and ionic strength 0.01 mM(λD=100 nm).

    To test Eqs. (13a) and (13b) at the interface that divides the inner and outer flow fields,Felectricin Eq.(13b)is the sum of electric forces inside the interface that divides the inner and outer flow fields. This force can be expressed as:

    where ρ(r)=e[p(r)-n(r)]. On the other hand,Fviscstands for the drag force produced by the flow over the same inner/outer interface. It can be written as:

    where σrrand σrθare the normal and tangential elements of the hydrodynamic stress tensor. From Eq. (13b), Fextis the negative of the sum of these two forces.

    Fig.6. The difference between γS and γE,using Eq.(13). (a)External force Fext at interfacial region shows linear dependence with sphere velocity v. The open circles indicate five different values of v. The red dashed line denotes the best fitting with a slope 1.45×10-8 kg·s-1,which is very close to the Stokes drag coefficient 1.41×10-8 kg·s-1(=6πηa, 3% lower than the slope). The external force becomes zero when v is exactly at the electrophoretic velocity. The external electric field strength was maintained at E∞=500 V/m. (b)Viscous force Fvisc at interfacial region with sphere velocity v=μEE∞. The external electric field strength is varied from 100 V/m to 500 V/m in 100 V/m steps.The red dashed line indicates fitting with slope 2.9×10-8 kg·s-1(=2.1×6πηa). All simulations were performed with surface charge density σ =-7×103 e/μm2, radius a=0.75 μm, and ionic strength 0.01 mM(λD=100 nm).

    By maintaining the solid particle stationary and varying the boundary value of the far field fluid velocity v,we evaluate the magnitude of Fextas a function of v,with the externally applied electric field maintained at E∞=500 V/m. We show in Fig.6(a)that the external force Fextdisplays a linear dependence on the relative velocity between the solid particle and the bulk fluid, with a slope given by 1.03γS. That is, equation(13a)is verified:

    It should be noted that even though the net Fextis evaluated at r=a+Δ, the fact that the inner flow field is carried along by the solid particle implies the net force Fextto be the same as that evaluated at the center of mass of the particle,as well as at r=a.

    In Fig.6(b), we focus on the case when Fext=0 so that Feletric=QeffE=Fvisc=γEv.By evaluating the viscous drag at the inner/outer interface and varying the external electric field strength from 100 V/m to 500 V/m in 100 V/m steps,we plot the resulting Fviscas a function of the electrophoretic velocity μEE∞. The slope gives a value that is 2.1 times the Stokes drag coefficient, which is somewhat higher than the experimentally measured value of γE/γS=1.7±0.3. However, both show unmistakable deviation from the Stokes drag coefficient. Therefore, we have explicitly used Long et al.’s analysis, based on the work of O’Brien and White, to show both the correctness of the equation Fext-γsΔv=0, as well as the necessary difference between the Stokes drag coefficient and the electrophoretic drag coefficient, both evaluated at the inner/outer interface.

    6. Conclusions

    We have measured electrophoresis drag coefficient γEby using an AC electric field to drive a charged particle in an optical trap. As far as we are aware, this is the first time such non-Stokes electrophoretic drag coefficients are determined experimentally. Our measurements were possible because our experimental arrangement permitted a time-varying electrophoretic motion of the charged particle with a speed directly proportional to that of the corresponding electric field strength.The non-electric, out-of-phase force from optical trap, fully counter-balanced by the out-of-phase electric force on the particle,does not play a role in hindering the electrophoretic motion. From both the experimental measurements and the direct numerical solutions, we show that the observed non-Stokes electrophoretic drag coefficient can be quantitatively described by the uniquely distinctive inner/outer interface as the reference plane for the evaluation of the hydrodynamic drag. The measured effective charge, Qeff, is consistent with the measured drag coefficient as required by force balance at the inner/outer interface flow fields which is at a distance a few λD’s away from the solid surface. The present experimental-theory study provides new,microscopic insights into the classic problem of electrophoresis. Such understanding can shed light on new applications of electrophoresis to investigate biological nanoparticles. More broadly, the experimental approach taken by this study and the microscopic view of an interface between the inner-outer flows might shed light to the outstanding problems of hydrodynamic forces and drags of other types of phoretic particles or swimming microorganisms.[41]

    Acknowledgment

    Ping Sheng wishes to acknowledge the funding support by Hong Kong RGC grants GRF16303415 and C6004-14G for this work.H.Daniel Ou-Yang wishes to acknowledge the support by US NSF grant 0923299 and the Emulsion Polymers Institute. Ping Sheng,Maijia Liao,and Shixin Xu wish to thank Chun Liu for advice regarding the necessary condition(s) for the appearance of hydrodynamic vortices. H.Daniel Ou-Yang and Ming-Tzo Wei would like to thank Kathryn Reddy,a former undergraduate student of Fordham University who participated in the initial experiments of this study during her REU at Lehigh. H.Daniel Ou-Yang would also like to thank Prof.Joel A.Cohen,Prof.Daan Frenkel,and Dr.Alois W¨urger for their critical comments and helpful suggestions.

    Author contributions

    Ping Sheng and H. Daniel Ou-Yang initiated and supervised the research. Maijia Liao and Ming-Tzo Wei carried out the experiments. Maijia Liao and Ping Sheng contributed to theory and simulations with help from Shixin Xu. Maijia Liao, H. Daniel Ou-Yang, and Ping Sheng analyzed the data with help from Ming-Tzo Wei. Maijia Liao,Ping Sheng,and H. Daniel Ou-Yang wrote the draft manuscript. All participated in revising the manuscript to its final form. Maijia Liao and Ming-Tzo Wei contributed equally to this work.

    精品国产一区二区久久| a级毛片在线看网站| 亚洲国产精品国产精品| 啦啦啦视频在线资源免费观看| 日韩一卡2卡3卡4卡2021年| 国产熟女欧美一区二区| 捣出白浆h1v1| 国产xxxxx性猛交| 秋霞在线观看毛片| 久久影院123| 日韩大片免费观看网站| 国产熟女欧美一区二区| 亚洲第一av免费看| 欧美国产精品va在线观看不卡| av片东京热男人的天堂| 制服诱惑二区| 90打野战视频偷拍视频| 巨乳人妻的诱惑在线观看| 国产日韩欧美亚洲二区| 99久国产av精品国产电影| 欧美少妇被猛烈插入视频| 日韩中文字幕视频在线看片| 一区二区av电影网| 边亲边吃奶的免费视频| 在线天堂最新版资源| 色网站视频免费| 亚洲久久久国产精品| 在线观看免费视频网站a站| 丝袜喷水一区| 日韩制服丝袜自拍偷拍| 成人二区视频| 99热全是精品| 老司机影院成人| 不卡av一区二区三区| 国产精品久久久久久av不卡| 国产成人av激情在线播放| 亚洲av日韩在线播放| 亚洲,欧美,日韩| 成年动漫av网址| 伊人久久大香线蕉亚洲五| 女人被躁到高潮嗷嗷叫费观| 一二三四在线观看免费中文在| 亚洲精品av麻豆狂野| 精品第一国产精品| 一本—道久久a久久精品蜜桃钙片| 久久精品久久久久久噜噜老黄| 欧美日韩一区二区视频在线观看视频在线| 国产在视频线精品| 国产视频首页在线观看| 久久久精品国产亚洲av高清涩受| 大香蕉久久成人网| 高清黄色对白视频在线免费看| 老司机影院毛片| 久久久久久久久久久免费av| 另类亚洲欧美激情| 啦啦啦啦在线视频资源| 亚洲五月色婷婷综合| 在线 av 中文字幕| av电影中文网址| 久久国产亚洲av麻豆专区| 色视频在线一区二区三区| av网站在线播放免费| 久久久久视频综合| 日本vs欧美在线观看视频| 欧美精品一区二区免费开放| 永久免费av网站大全| 国产熟女欧美一区二区| 日韩,欧美,国产一区二区三区| 国产精品久久久av美女十八| 亚洲av在线观看美女高潮| 国产精品蜜桃在线观看| 一区福利在线观看| 免费黄频网站在线观看国产| 老熟女久久久| 在线免费观看不下载黄p国产| 免费观看无遮挡的男女| 青青草视频在线视频观看| 成人影院久久| 日韩制服骚丝袜av| 26uuu在线亚洲综合色| 亚洲av中文av极速乱| 边亲边吃奶的免费视频| 黑丝袜美女国产一区| 亚洲人成网站在线观看播放| 国产亚洲av片在线观看秒播厂| 国产深夜福利视频在线观看| 美女视频免费永久观看网站| 亚洲欧美日韩另类电影网站| 黄片播放在线免费| 亚洲国产精品一区三区| 大码成人一级视频| 精品一区二区三区四区五区乱码 | 少妇熟女欧美另类| 不卡av一区二区三区| av电影中文网址| 一级毛片我不卡| 久久ye,这里只有精品| 亚洲欧洲日产国产| 婷婷色综合大香蕉| 亚洲一区二区三区欧美精品| 男女啪啪激烈高潮av片| 亚洲欧美一区二区三区黑人 | 日本猛色少妇xxxxx猛交久久| 色婷婷av一区二区三区视频| 国产精品一区二区在线观看99| 国产成人精品福利久久| 高清视频免费观看一区二区| 国产亚洲精品第一综合不卡| 国产成人aa在线观看| av天堂久久9| 亚洲欧美成人综合另类久久久| 欧美xxⅹ黑人| 美女脱内裤让男人舔精品视频| 日韩熟女老妇一区二区性免费视频| 亚洲av电影在线进入| 一区二区日韩欧美中文字幕| 日韩欧美一区视频在线观看| videos熟女内射| 十八禁高潮呻吟视频| 精品亚洲乱码少妇综合久久| 男女午夜视频在线观看| 老女人水多毛片| 久久久久久久亚洲中文字幕| 欧美av亚洲av综合av国产av | 亚洲成色77777| 成人免费观看视频高清| 91精品三级在线观看| 女性生殖器流出的白浆| 高清视频免费观看一区二区| 亚洲经典国产精华液单| 久久亚洲国产成人精品v| 午夜老司机福利剧场| 水蜜桃什么品种好| 9热在线视频观看99| 久热久热在线精品观看| 搡女人真爽免费视频火全软件| 最新的欧美精品一区二区| 精品国产露脸久久av麻豆| 伊人亚洲综合成人网| 亚洲欧美精品自产自拍| 亚洲久久久国产精品| 亚洲中文av在线| 午夜福利在线免费观看网站| 欧美xxⅹ黑人| 1024视频免费在线观看| 日本-黄色视频高清免费观看| 高清黄色对白视频在线免费看| 国产免费现黄频在线看| 黄片小视频在线播放| 色94色欧美一区二区| www.av在线官网国产| 女性生殖器流出的白浆| 久久国产精品大桥未久av| 丝袜美腿诱惑在线| 叶爱在线成人免费视频播放| 亚洲欧美一区二区三区黑人 | 日韩中文字幕视频在线看片| 国产又爽黄色视频| 亚洲精品国产av蜜桃| 国产一区有黄有色的免费视频| 久久国产精品男人的天堂亚洲| 国产色婷婷99| 欧美变态另类bdsm刘玥| 精品人妻在线不人妻| kizo精华| 欧美老熟妇乱子伦牲交| 捣出白浆h1v1| 亚洲国产日韩一区二区| 亚洲视频免费观看视频| 欧美日本中文国产一区发布| 国产白丝娇喘喷水9色精品| 丝袜美腿诱惑在线| 久久久亚洲精品成人影院| 欧美精品一区二区大全| 丝袜人妻中文字幕| 伊人亚洲综合成人网| 男女啪啪激烈高潮av片| 国产有黄有色有爽视频| 在线看a的网站| 欧美日韩成人在线一区二区| 一级,二级,三级黄色视频| av在线观看视频网站免费| 国产极品天堂在线| av在线观看视频网站免费| 欧美bdsm另类| 日本vs欧美在线观看视频| www.精华液| 色婷婷久久久亚洲欧美| 91成人精品电影| 国产一区二区三区综合在线观看| 最近最新中文字幕大全免费视频 | 亚洲精品第二区| 欧美av亚洲av综合av国产av | 热99久久久久精品小说推荐| 午夜福利一区二区在线看| 久久久久精品人妻al黑| 九草在线视频观看| 欧美老熟妇乱子伦牲交| 亚洲成色77777| 精品福利永久在线观看| 成人毛片60女人毛片免费| 免费女性裸体啪啪无遮挡网站| 母亲3免费完整高清在线观看 | 黄频高清免费视频| 亚洲一码二码三码区别大吗| 视频在线观看一区二区三区| 看非洲黑人一级黄片| 欧美国产精品一级二级三级| 啦啦啦视频在线资源免费观看| xxxhd国产人妻xxx| 少妇猛男粗大的猛烈进出视频| 我要看黄色一级片免费的| 青春草国产在线视频| 成年av动漫网址| 人妻 亚洲 视频| 九草在线视频观看| 国产成人精品无人区| 日本-黄色视频高清免费观看| 国产成人精品无人区| 免费人妻精品一区二区三区视频| 欧美亚洲 丝袜 人妻 在线| 视频区图区小说| 亚洲精品国产一区二区精华液| 激情五月婷婷亚洲| 欧美中文综合在线视频| av片东京热男人的天堂| 亚洲精品久久久久久婷婷小说| 老司机影院毛片| 熟女电影av网| 中文字幕人妻丝袜制服| 午夜av观看不卡| 亚洲av.av天堂| √禁漫天堂资源中文www| 亚洲欧洲精品一区二区精品久久久 | 看免费成人av毛片| av有码第一页| 国产福利在线免费观看视频| 男女下面插进去视频免费观看| 最新的欧美精品一区二区| 成人国产麻豆网| 日韩av免费高清视频| 色视频在线一区二区三区| 国产综合精华液| 性色avwww在线观看| 建设人人有责人人尽责人人享有的| 青春草国产在线视频| 亚洲国产成人一精品久久久| 国产日韩欧美在线精品| 一本一本久久a久久精品综合妖精 国产伦在线观看视频一区 | 不卡av一区二区三区| 欧美精品亚洲一区二区| 久久久久久久国产电影| 精品视频人人做人人爽| 国产精品二区激情视频| 日本免费在线观看一区| 91国产中文字幕| 18+在线观看网站| 精品人妻偷拍中文字幕| 亚洲综合色惰| 男女边吃奶边做爰视频| 制服人妻中文乱码| 日本爱情动作片www.在线观看| 国产综合精华液| 丝袜脚勾引网站| 丝袜在线中文字幕| 国产不卡av网站在线观看| 欧美变态另类bdsm刘玥| 成年人免费黄色播放视频| 国产女主播在线喷水免费视频网站| 亚洲色图 男人天堂 中文字幕| 高清av免费在线| 日日摸夜夜添夜夜爱| 亚洲四区av| 亚洲五月色婷婷综合| 国产在线一区二区三区精| 亚洲人成网站在线观看播放| 人体艺术视频欧美日本| 免费高清在线观看视频在线观看| 久久久久久伊人网av| 97在线视频观看| av国产久精品久网站免费入址| 国产精品99久久99久久久不卡 | 精品久久蜜臀av无| 最近中文字幕高清免费大全6| 国产黄色视频一区二区在线观看| 精品人妻一区二区三区麻豆| 青春草视频在线免费观看| 五月开心婷婷网| 叶爱在线成人免费视频播放| 91aial.com中文字幕在线观看| 一级片'在线观看视频| 校园人妻丝袜中文字幕| 高清欧美精品videossex| 美女视频免费永久观看网站| 久久久国产欧美日韩av| 成人午夜精彩视频在线观看| 肉色欧美久久久久久久蜜桃| 欧美日韩精品成人综合77777| 99久久中文字幕三级久久日本| 欧美 亚洲 国产 日韩一| 久久这里有精品视频免费| 国产精品无大码| 香蕉国产在线看| 日韩一本色道免费dvd| 伦理电影免费视频| 亚洲精品中文字幕在线视频| 欧美激情极品国产一区二区三区| 黄片无遮挡物在线观看| 亚洲情色 制服丝袜| 中文字幕制服av| 日本免费在线观看一区| 午夜久久久在线观看| 建设人人有责人人尽责人人享有的| 波野结衣二区三区在线| 中文精品一卡2卡3卡4更新| 欧美成人午夜精品| 久久久国产精品麻豆| 国产精品亚洲av一区麻豆 | 在线观看国产h片| 丰满饥渴人妻一区二区三| 成年人免费黄色播放视频| 国产一区二区在线观看av| 欧美国产精品va在线观看不卡| 精品视频人人做人人爽| 视频在线观看一区二区三区| av国产精品久久久久影院| 丝瓜视频免费看黄片| 久久免费观看电影| 成人漫画全彩无遮挡| 国产爽快片一区二区三区| 免费女性裸体啪啪无遮挡网站| 午夜久久久在线观看| 欧美日韩av久久| 亚洲精品一二三| 亚洲色图 男人天堂 中文字幕| 亚洲国产最新在线播放| 久热这里只有精品99| 香蕉国产在线看| 免费在线观看视频国产中文字幕亚洲 | 亚洲精品日本国产第一区| 亚洲人成网站在线观看播放| 国产黄频视频在线观看| 日本av手机在线免费观看| 久久久久久人妻| 大话2 男鬼变身卡| 男女国产视频网站| 久久精品久久久久久久性| 国产成人精品福利久久| 美女国产视频在线观看| av电影中文网址| 亚洲三级黄色毛片| 久久99蜜桃精品久久| 精品少妇内射三级| 欧美精品高潮呻吟av久久| 丝袜人妻中文字幕| 午夜免费鲁丝| 午夜福利视频在线观看免费| 电影成人av| av国产精品久久久久影院| 丝袜在线中文字幕| 亚洲人成77777在线视频| 欧美中文综合在线视频| 在线亚洲精品国产二区图片欧美| 欧美日韩一区二区视频在线观看视频在线| 91久久精品国产一区二区三区| 纯流量卡能插随身wifi吗| 丝袜美腿诱惑在线| 91久久精品国产一区二区三区| 成年av动漫网址| 高清黄色对白视频在线免费看| 久久久久久久亚洲中文字幕| 中文天堂在线官网| 日韩成人av中文字幕在线观看| av在线老鸭窝| 日产精品乱码卡一卡2卡三| 少妇猛男粗大的猛烈进出视频| 在线观看免费视频网站a站| 亚洲精品,欧美精品| 天天躁日日躁夜夜躁夜夜| 国产av一区二区精品久久| 波多野结衣一区麻豆| 一级毛片我不卡| 国产不卡av网站在线观看| 精品国产一区二区久久| 亚洲一区中文字幕在线| 伦精品一区二区三区| 一区二区av电影网| 久久精品夜色国产| 成人亚洲精品一区在线观看| 国产男女内射视频| 欧美成人午夜免费资源| 麻豆av在线久日| 国产免费一区二区三区四区乱码| 男人添女人高潮全过程视频| 大香蕉久久成人网| 成人影院久久| 丁香六月天网| 老鸭窝网址在线观看| 久久国产亚洲av麻豆专区| 亚洲三级黄色毛片| 久久精品久久久久久噜噜老黄| 在线 av 中文字幕| 十八禁网站网址无遮挡| 亚洲成国产人片在线观看| 一二三四在线观看免费中文在| 成年人免费黄色播放视频| 午夜激情久久久久久久| 欧美精品亚洲一区二区| 欧美日韩综合久久久久久| 美女xxoo啪啪120秒动态图| 老女人水多毛片| 国产精品成人在线| 久久97久久精品| 久久狼人影院| 久久精品国产自在天天线| 久久久久精品性色| 国产在线免费精品| 午夜久久久在线观看| 老熟女久久久| 91午夜精品亚洲一区二区三区| 国产成人av激情在线播放| 巨乳人妻的诱惑在线观看| 五月伊人婷婷丁香| 美国免费a级毛片| 亚洲久久久国产精品| 欧美精品国产亚洲| 18禁观看日本| 汤姆久久久久久久影院中文字幕| 一本大道久久a久久精品| 黄色一级大片看看| 免费av中文字幕在线| 999久久久国产精品视频| 这个男人来自地球电影免费观看 | 美女福利国产在线| 亚洲人成网站在线观看播放| 成人国产av品久久久| 极品人妻少妇av视频| 热re99久久精品国产66热6| 交换朋友夫妻互换小说| 亚洲国产欧美网| 亚洲成av片中文字幕在线观看 | 曰老女人黄片| 99精国产麻豆久久婷婷| 日韩成人av中文字幕在线观看| 久久国产精品男人的天堂亚洲| 免费播放大片免费观看视频在线观看| 国产精品99久久99久久久不卡 | 久久精品国产亚洲av高清一级| 一本—道久久a久久精品蜜桃钙片| 妹子高潮喷水视频| 最近最新中文字幕大全免费视频 | 欧美最新免费一区二区三区| 亚洲情色 制服丝袜| 久久国产精品男人的天堂亚洲| 日本色播在线视频| 男人爽女人下面视频在线观看| 黑丝袜美女国产一区| 天美传媒精品一区二区| 国产xxxxx性猛交| 精品国产一区二区久久| 999精品在线视频| 国产成人精品久久久久久| 国产视频首页在线观看| 日本欧美视频一区| 欧美bdsm另类| 亚洲五月色婷婷综合| 黄频高清免费视频| 亚洲,一卡二卡三卡| 免费看不卡的av| 成年女人毛片免费观看观看9 | 水蜜桃什么品种好| 叶爱在线成人免费视频播放| 久久久久人妻精品一区果冻| 欧美在线黄色| 菩萨蛮人人尽说江南好唐韦庄| 香蕉精品网在线| 男女边摸边吃奶| 91精品伊人久久大香线蕉| 80岁老熟妇乱子伦牲交| 亚洲国产欧美在线一区| 十八禁高潮呻吟视频| 成人黄色视频免费在线看| 99久久精品国产国产毛片| 夫妻性生交免费视频一级片| 欧美精品国产亚洲| 国产欧美日韩一区二区三区在线| 2022亚洲国产成人精品| 亚洲精品日本国产第一区| 亚洲精品美女久久av网站| 久久热在线av| 伊人久久国产一区二区| 如日韩欧美国产精品一区二区三区| 亚洲欧美中文字幕日韩二区| 波野结衣二区三区在线| 啦啦啦中文免费视频观看日本| 人人妻人人爽人人添夜夜欢视频| 有码 亚洲区| 国产色婷婷99| 国产欧美日韩综合在线一区二区| 美女视频免费永久观看网站| 狠狠婷婷综合久久久久久88av| 一二三四中文在线观看免费高清| 久久精品国产亚洲av涩爱| 熟妇人妻不卡中文字幕| 午夜日本视频在线| videossex国产| 亚洲三级黄色毛片| 国产精品人妻久久久影院| 欧美 日韩 精品 国产| 亚洲天堂av无毛| 麻豆精品久久久久久蜜桃| 丝袜人妻中文字幕| 国产在视频线精品| 国产探花极品一区二区| 国产综合精华液| 另类亚洲欧美激情| 亚洲成国产人片在线观看| 少妇精品久久久久久久| 亚洲色图 男人天堂 中文字幕| 在线 av 中文字幕| 精品国产一区二区三区久久久樱花| 国产1区2区3区精品| 成人免费观看视频高清| 中国国产av一级| 欧美老熟妇乱子伦牲交| 精品一区在线观看国产| 久久精品aⅴ一区二区三区四区 | 春色校园在线视频观看| 大香蕉久久网| 国产成人a∨麻豆精品| av有码第一页| 卡戴珊不雅视频在线播放| 国产精品二区激情视频| 青春草亚洲视频在线观看| 欧美黄色片欧美黄色片| 久久影院123| 亚洲精华国产精华液的使用体验| 一边亲一边摸免费视频| 中文字幕亚洲精品专区| 久久精品国产自在天天线| 男男h啪啪无遮挡| 性高湖久久久久久久久免费观看| 日韩成人av中文字幕在线观看| 国产又色又爽无遮挡免| 精品一区在线观看国产| 婷婷成人精品国产| 2022亚洲国产成人精品| www日本在线高清视频| 男女免费视频国产| 日本欧美国产在线视频| 午夜福利在线观看免费完整高清在| 黄片播放在线免费| 国产成人午夜福利电影在线观看| 国产精品亚洲av一区麻豆 | 国产亚洲av片在线观看秒播厂| 视频区图区小说| 久久精品国产亚洲av天美| 欧美bdsm另类| 自线自在国产av| 亚洲成国产人片在线观看| 久久狼人影院| 国产亚洲午夜精品一区二区久久| 美女视频免费永久观看网站| 寂寞人妻少妇视频99o| 王馨瑶露胸无遮挡在线观看| 成人亚洲精品一区在线观看| 99久久综合免费| 美女脱内裤让男人舔精品视频| 人妻人人澡人人爽人人| 成年美女黄网站色视频大全免费| 一本大道久久a久久精品| 狠狠婷婷综合久久久久久88av| 精品卡一卡二卡四卡免费| 99国产综合亚洲精品| 超碰成人久久| 黑丝袜美女国产一区| 亚洲国产精品一区三区| 久久鲁丝午夜福利片| 观看美女的网站| 亚洲精品美女久久久久99蜜臀 | 日韩在线高清观看一区二区三区| 久久久久人妻精品一区果冻| 亚洲男人天堂网一区| 免费观看无遮挡的男女| 亚洲成色77777| 欧美日韩国产mv在线观看视频| 大陆偷拍与自拍| 超碰成人久久| 国产精品 国内视频| 一二三四在线观看免费中文在| 亚洲一区二区三区欧美精品| 日韩一区二区三区影片| av福利片在线| 国产在线一区二区三区精| 十分钟在线观看高清视频www| 丝袜喷水一区| 可以免费在线观看a视频的电影网站 | 久久人人爽av亚洲精品天堂| 99re6热这里在线精品视频| 国产成人免费无遮挡视频| 99国产综合亚洲精品| 在现免费观看毛片| 国产精品 国内视频| 超碰97精品在线观看| 一区二区av电影网| 日韩制服丝袜自拍偷拍| 久久久久久伊人网av| 国产又爽黄色视频| 久久久精品94久久精品| 尾随美女入室| 国产免费视频播放在线视频| 精品国产一区二区三区久久久樱花| 女性被躁到高潮视频| 伊人亚洲综合成人网| 国产成人欧美|