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

    Shape and diffusion instabilities of two non-spherical gas bubbles under ultrasonic conditions

    2024-01-25 07:30:14WurihanBao包烏日汗andDeXinWang王德鑫
    Chinese Physics B 2024年1期
    關(guān)鍵詞:烏日

    Wurihan Bao(包烏日汗) and De-Xin Wang(王德鑫)

    College of Physics and Electronics,Inner Mongolia Minzu University,Tongliao 028043,China

    Keywords: non-spherical bubble,shape instability,diffusive instability

    1.Introduction

    When sound waves propagate through a liquid, they induce a sequence of cavitation phenomena within minute bubbles immersed in the liquid due to the influence of a driving acoustic field.[1]As the driving acoustic pressure attains the inertial cavitation threshold, the bubbles collapse, accompanied by concurrent light radiation under specific conditions,a phenomenon termed sonoluminescence.[2,3]In recent years,owing to the rapid progress in science and technology,the extensive application of ultrasonic cavitation technology has become prevalent across various domains, encompassing ultrasonic cleaning, pulverization, sterilization, ultrasonic extraction,and ultrasonic therapy.

    The presence of multiple bubbles is primarily responsible for the acoustic cavitation effect in liquid.However,the intricate interplay between these bubbles and environmental conditions poses challenges for accurately simulating diverse experimental scenarios using numerical modeling.In contrast,simulating the cavitation of two bubbles in liquid is considerably more straightforward compared to the complexities associated with multiple-bubble acoustic cavitation.When subjected to a driving acoustic field, the oscillation of bubbles give rise to mutual attractive or repulsive forces, known as the secondary Bjerknes forces.Under the influence of the secondary Bjerknes forces, distinct acoustic phenomena occur on the bubbles immersed in the liquid,deviating from those observed in isolated bubble systems.[4–6]

    Acoustic cavitation bubbles, renowned for their intricate dynamics, have instigated a comprehensive exploration spanning both theoretical and experimental domains.In this captivating realm, the attention has particularly focused on the dynamics of two gas bubbles due to their versatile applications.Noteworthy is the work of Luet al.,[7]who extended the theoretical framework by employing the velocity potential superposition theory to elucidate the complexities inherent in the ultrasonic dynamics of two gas bubbles with distinct frequencies.Their research not only sheds light on the intricate interactions intrinsic to bubble behavior but also establishes the groundwork for a deeper comprehension of their dynamics.Empirical insights further enrich this landscape.Maet al.[8–10]delved into the dynamics of bubble volume under secondary Bjerknes forces,utilizing the Lagrange equation to reveal the interplay of parameters influencing oscillation amplitude and initial phase.This empirical depth not only complements theoretical insights but also provides tangible validation through experimental observations.Extending the inquiry to non-spherical bubbles,Wanget al.[11]harnessed the modified Keller–Miksis equation to quantify fluctuations in the radius of two gas bubbles,secondary Bjerknes forces,and temperature during acoustic cavitation.These experimental findings seamlessly align with theoretical constructs,fortifying the foundational principles governing bubble dynamics.Moreover,the pioneering study conducted by Zilonovaet al.[12]has delved into the intriguing dynamics of bubble-bubble interactions within viscous media, revealing their pronounced influence on the behavior of individual bubbles.That study notably highlights a distance-dependent relationship, wherein mutual influence diminishes as the inter-bubble distance increases,becoming negligible at distances of a few millimeters.This insight assumes even greater significance when extended to the realm of viscous media,warranting further exploration.In the broader context of the cavitation of two bubbles, it is noteworthy to acknowledge the significant contributions of other researchers.Studies conducted in Refs.[13–19] have illuminated various facets of cavitation involving two bubbles,collectively enhancing our understanding of this intricate phenomenon.

    In all of these investigations, Prosperettiet al.have predominantly focused on the stability characteristics of spherical bubbles, primarily considering the influence of fluid viscosity.However, the stability characteristics of non-spherical bubbles have not been addressed in their research.[20]Studies on bubble dynamics often rely on the assumption of ideal spherical bubbles.However, real-world experiments involve the presence of numerous non-spherical bubbles,rendering the spherical model inapplicable.Wanget al.[21]have proposed a non-spherical model for sonoluminescence in a non-spherical acoustic field and have determined a stable range for sonoluminescence.They have discovered that bubbles exhibit nonspherical oscillations within specific parameter intervals.Wuet al.have investigated the translation and non-spherical oscillation of a single bubble, concluding the work that, under a constant initial bubble radius and driving acoustic pressure amplitude, the non-spherical oscillation becomes more pronounced with an increasing initial speed of the bubble center.Moreover, the corresponding instability region in theR0–Paphase diagram gradually expands.[22]

    In the present study,based upon literature,[16]we investigate the impact of various initial radius, distance, and perturbation parameters on the shape instability and diffusive equilibrium curve properties of two non-spherical gas bubbles.To accomplish this,we employ the modified Keller–Miksis equation.Our findings contribute theoretical evidence for the examination of instability characteristics in non-spherical multibubbles and bubble clouds.

    2.Theoretical model for two non-spherical gas bubbles

    To facilitate the analysis, we initially developed a simplified theoretical model of two non-spherical gas bubbles,as illustrated in Fig.1.Although this is a schematic diagram of a three-dimensional model,the bubble dynamics model under the three-dimensional model is very complex.In this study,only two non-spherical gas bubbles under two-dimensional conditions are considered, and their long radii of ellipses are represented asR1andR2respectively.It can be observed from Fig.1 thatr2=r1?d,whered=(0,0,d)represents the coordinate vector of the displacementdbetween bubbles.Assuming thatR1+R2?d, we consider a random pointPlocated outside the two bubbles,r1andr2represent the displacements of the centers of bubbles 1 and 2 from the pointP.Under the assumptions of an incompressible and irrotational ideal fluid,the particle velocity potential function at any point within the liquid satisfies the Laplace equation

    If we take the central point of bubble 1 as the origin coordinate, the velocity potential of the liquid due solely to radial pulsations of bubble 1 is then represented as[16]

    Similarly,potential generated by radial pulsations of bubble 2 would be

    whereΦ1(r1)is the velocity potentials of bubble 1,Φ2(r2)is the velocity potentials of bubble 2, respectively.To carry out the calculation further we need to be able to transform coordinates from(r1,θ1)to(r2,θ2),and vice versa(see Fig.1).Herepnis thenthLegendre polynomial.We can expandΦ2in the local coordinate of bubble 1 as

    Similarly,we haveΦ1in the local coordinate of bubble 2,

    Fig.1.A three-dimensional schematic diagram of two non-spherical gas bubbles in an acoustic field.

    To describe the oscillation of two non-spherical gas bubbles in the acoustic field, we introduce a perturbation in the driving acoustic pressure.The distance between the center point of any non-spherical bubble and any pointPon the surface of that bubble can be expressed as follows:[17]

    whereRi(t) is the initial radius of the bubble with no perturbation.Ymn(θ,Φ) represents the spherical harmonics ofnorder, wheren ≥2, andai(t) is the amplitude of the surface distortion.During the oscillation process, the bubbles are rotationally symmetric.Then,we only need to consider the bubbles’radial oscillation,and their velocity potentials can be expressed as

    We further view the part containingθiin Eq.(7)as a perturbation item and neglect its high order small quantities,then we obtain the velocity potentials inside and outside the nonspherical bubble interfaces as follows:

    Then substituting Eq.(8) into the Bernoulli equation yields

    whereP(ri)andP(∞)denote the pressures at the positionsPand infinity in the liquid,εis a small perturbation parameter,

    whereμis the viscosity coefficient of the liquid, andσis the surface tension coefficient.P0is the ambient pressure,Pd=?Pasin(ωt)is the driving acoustic pressure on the bubble surface,ω=2π fwithfbeing the frequency of acoustic field,andγdenotes the polytropic exponent.

    WhenPdis not too large, because of its symmetry, nonspherical bubbles from persisting for an extended period under uniform acoustic driving.Therefore, it becomes necessary to introduce a small non-spherical symmetrical acoustic pressure to correct the driving pressure:[23]

    whereδis the boundary layer thickness,(i=1,2), onlyn=2 is considered in this study, quadruple distortion.By substituting the two-non-spherical-gas-bubble potentialΦ(ri,θi) into Eq.(9), whenr=R1, we can obtain the shape instability equation of bubble 1 as follows:

    Similarly,whenr=R2,we can obtain the shape instability equation of bubble 2 as follows:

    The right sides of Eqs.(14) and (16) describe the spherical asymmetric perturbation of sound waves,[23]whereδ pnrepresents the deviation of the driving sound pressure from the spherical symmetry.

    Considering the compressibility of the liquid in a highly intense acoustic field,the Rayleigh–Plesset model cannot satisfy the requirement that the radial oscillation velocity of bubbles in the liquid is lower than the speed of sound.Consequently, in this study, we adopt the assumption of incompressible liquid.To meet this condition, we used the Keller–Miksis equation to describe the variation in the radius of two non-spherical gas bubbles, instead of the Rayleigh–Plesset equations (13) and (15).The Keller–Miksis equation is expressed as

    3.Factors influencing the instability characteristics of two non-spherical gas bubbles

    Compared to a single bubble, calculating the shape instabilityR0–Paphase diagram and diffusive equilibrium curve of two non-spherical gas bubbles requires more extensive effort and entails increased complexity.Therefore, our primary focus in this paper is to examine the influence of bubble radius, distance and perturbation parameter on the instability characteristics of cavitation bubbles.Based on previous research,[23,24]the instability of non-spherical bubbles can be classified into three main types: shape instability,diffusion instability, and chemical instability.The shape instability type can be further classified into three categories: (a) Rayleigh–Taylor(RT)instability,(b)rebound instability,and(c)parameter instability.Studies have confirmed that RT and parameter instabilities play a dominant role in the stability characteristics of bubbles, with their respectiveR0–Paphase diagrams exhibiting comparable features.Therefore, our investigation primarily focuses on the RT instability, for which the corresponding instability criterion can be expressed as[25]

    To provide a more intuitive description of the instability characteristics of non-spherical bubbles,we conducted numerical simulations to observe the instability states over five acoustic cycles.We incorporated the parameter conditions for the three instabilities of a single bubble as provided in the literature[24]into our modified models for two non-spherical gas bubbles.Subsequently, we investigated the instability characteristics of the bubbles by considering the parameters of perturbation amplitude,which changes over time,and the bubble radii.These parameters are as follows: (a)R0=2.5 μm,Pa=1.5 atm;(b)R0=4.0μm,Pa=1.3 atm;(c)R0=5.2μm,Pa=1.0 atm.We performed calculations and generated diagrams illustrating the normalized perturbation amplitude as it varies with time and the bubbles’radii,presented in Figs.2(a)–2(c).Using Eqs.(13)–(16) and employing the Runge–Kutta method,we conducted numerical simulations for the three instabilities of non-spherical bubble.In our proposed model,we assigned the following physical parameters:f=25 kHz,ρ=1000 kg/m3,c=1485 m/s,γ=1.4,P0=1.013×105Pa,P0=P(∞),σ=0.0725 N/m,andμ=0.001 kg/m·s.

    Figure 2(solid lines)illustrates the behavior of cavitation bubbles over five complete acoustic cycles.Within each cycle,the bubbles initially undergo slow expansion.During the first 0.5 cycles of acoustic pressure, their radii reach a peak and then sharply decrease.Subsequently, within 0.6 cycles of acoustic pressure, the radii reach their minima.Following this, the bubbles oscillate in a stable manner.The corresponding maximal expansion ratios of the cavitation bubbles for the three instabilities are 1.7, 8.0, and 20, respectively.It can be observed that the expansion ratios of the bubbles differ depending on their initial radius.Under the same amplitude of the driving pressure,it is observed that cavitation bubbles with larger initial radius tend to exhibit relatively higher expansion ratios.This phenomenon could potentially be attributed to the proximity of larger initial bubbles to the resonant radius within the examined parameters.Notably,previous studies[26–28]have demonstrated that the maximum expansion ratio is achieved when bubbles are driven under resonant conditions.The oscillation of the bubble radii can be roughly divided into three periods: expansion,collapse,and rebound.

    In order to explore the impact of the initial radius, distance,and perturbation parameters of two bubbles on their instability characteristics, we conducted numerical simulations to obtain the shape instabilityR0–Paphase diagrams and diffusive equilibrium curves.The diffusive equilibrium condition can be mathematically expressed as[25]

    Here,C∞is the gas concentration far away from the bubble,C0is the saturation(mass)concentration of the gas in liquid,and〈〉denotes the time average.

    Fig.2.The dashed lines correspond to the normalized perturbation amplitude a2/R0 varying with acoustic cycles; the solid lines correspond to the radius R/R0 varying with acoustic cycles: (a) Rayleigh–Taylor instability(R0 =2.5μm, Pa =1.5 atm), (b)rebound instability(R0 =4.0 μm, Pa =1.3 atm), (c) parameter instability (R0 =5.2 μm,Pa=1.0 atm).

    The shape instability caused by the interaction between two bubbles can be observed through the secondary Bjerknes force.[29]The secondary Bjerknes force between two gas bubbles can be expressed as

    whereV1andV2denote the volumes of bubbles 1 and 2, respectively,V1=4/3πR31andV2=4/3πR32;eris the radial unit vector.In order to present the results,the variation of the secondary Bjerknes force is often portrayed using the secondary Bjerknes force coefficient,denoted asfB.It is defined as

    Therefore,by observing the sign offB,we can discern the interaction between the two bubbles: whenfB>0, the bubbles attract each other;whenfB<0,the bubbles repel each other.

    3.1.Influence of the radius between bubbles

    Figures 3(b)–3(c)illustrate the influence of different initial radius of bubble 2 on the shape instability and the diffusion equilibrium curve of bubble 1 inR0–Paphase diagrams for the two-non-spherical-gas-bubble model presented in this paper.Figure 3(a) shows the results of a single bubble under the same conditions, in comparison with the two-bubble model.Bubbles are driven by sound waves with a frequency of 25 kHz in water at 293 K and 1 atm ambient pressure,without considering perturbation.In the figures,black areas represent regions of shape instability,whereas white areas represent regions of shape stability.The stability threshold corresponds to the boundary between the black and white areas of the figures.The curves represent the diffusive equilibrium of the bubble under different relative concentrations ofC∞/C0,with a positive slope indicating a stable diffusive equilibrium and a negative slope indicating an unstable diffusive equilibrium.

    Fig.3.Shape instability and diffusive equilibrium curves in R0–Pa phase diagrams, under a constant distance d =20 mm without perturbation for various relative concentrations C∞/C0: (a)single bubble,[(b),(c)]the initial radius of bubble 2 with R20=1μm and R20=5μm,respectively.

    It is known that in the case of single-bubble cavitation,after undergoing several oscillation cycles,due to rectification effects,the interior of the bubble eventually contains only argon gas.[30]The results presented in Figs.3(b) and 3(c) reveal important findings regarding the stability of two bubbles.When the distance between bubbles is kept constant, it becomes evident that the shape instability region of bubble 1 in theR0–Paphase diagram expands gradually as the initial radius of bubble 2 increases.Compared to the case of a single bubble,the unstable region has significantly increased,but its range is different.This could be due to the interaction forces between the two bubbles, specifically the secondary Bjerknes force, coming into play.In air and air-saturated water,the relative concentration of argon gas,C∞/C0=0.01.From Figs.3(b) and 3(c), it can be observed that with an increase in the initial radius of bubble 2, the curve ofC∞/C0=0.01 with a positive slope covers a decreasingParegion.This suggests that in terms of diffusive stability,due to the interaction forces between the two bubbles, the two-bubble model is not as stable as the single-bubble model under the same cavitation conditions.

    Additionally, from Figs.3(b) and 3(c), it is evident that some continuous stable regions emerge within the instability region.The stability theory of Mathieu(or Hill)equation can provide a good understanding of the general characteristics of instability.[31]Bubble oscillations are highly nonlinear,the Mathieu (or Hill) equation replacesan=bn/R3/2and retains only the linear terms in the spherical asymmetric perturbation.Calculations using the Mathieu(or Hill)equation can explain the emergence of continuous stable regions within the instability region of theR0–Paphase diagram.[32]

    3.2.Influence of the distance between bubbles

    The distance between bubbles plays a significant role in the secondary Bjerknes forces.Unlike the behavior of a single bubble, where stability and dynamics are relatively wellunderstood,the presence of two bubbles introduces a new level of complexity.Due to the influence of secondary Bjerknes forces and variations in the inherent resonance frequencies of the two bubbles,varying initial distances between bubbles lead to distinct interactions.To investigate the influence of the distance between bubbles on the shape instability of bubble 1 in the presence of two non-spherical gas bubbles,we conducted simulations and generated the shape instabilityR0–Paphase diagrams and diffusive equilibrium curves without perturbation.When the initial radius of bubble 2 is 5 μm, the initial distances between bubbles shown in Figs.4(a) and 4(b) are 2 mm and 0.2 mm,respectively.

    Fig.4.Shape instability and the diffusive equilibrium curves of bubble 1 in R0–Pa phase diagrams for various relative concentrations C∞/C0.The distance between bubbles is (a) d =2 mm and (b) d =0.2 mm,respectively.

    Comparing Fig.3(c) with Fig.4, one can conclude that the unstable region gradually increases with the reducing distance between bubbles, which may due to the frequent fluctuations in secondary Bjerknes forces.Therefore, we calculated the variation of the secondary Bjerknes force coefficientfBin theR10–R20plane with different distance between bubbles,as shown in Fig.5.According to the region boundary in Fig.4,we selected the driving sound pressure as 1.15 atm.The white regions represent repulsive force,(i.e,fB<0),while the gray scale regions represent attractive forces(i.e,fB>0).The darker the color, the higher the absolute value offB.For better visualization, the complete symmetric data are shown.In Fig.5,as the distance between the bubbles reduces(i.e.,they approach each other),the repulsion region gradually decreases and the boundaries of these regions shift towards larger bubble radius,while the changes in the attraction region are minimal.This phenomenon may be attributed to the natural resonance of bubbles increasing as they approach each other.[29]Consequently, when bubbles approach each other, there is a possibility that the sign offBfor the cases near the boundaries changes, potentially leading to an inversion of the secondary Bjerknes force in some bubble pairs.A reduction in the distance between bubbles results in the widening of pressure ranges associated with the positive slope regions in the diffusion curves, particularly when considering a relative concentration ofC∞/C0at 0.01.This observation aligns harmoniously with the outcomes documented in earlier research studies.[25]

    Fig.5.The variations of the secondary Bjerknes force coefficient fB in the R10–R20 plane for driving pressure Pa =1.15 atm under different distance between bubbles.The planes are given for distances: (a)d=20 mm,(b)d=2 mm,and(c)d=0.2 mm.

    3.3.Influence of the spherical asymmetric perturbation parameter

    The spherical asymmetric perturbation parameterδ p2in Eq.(12) is challenging to evaluate directly.In this paper,we have chosen to simplify by settingδ p2= 1×10?5Pa,1×10?7Pa,values that were previously used for a single bubble under similar conditions[23]and were demonstrated with the best fit to the experimental data in Ref.[33],and the parameters used in the instability analysis of other bubble models.[25]

    In Figs.6(a)–6(c), we can see that, in the bubble-1 shape instabilityR0–Paphase diagram and diffusive equilibrium curve under a varying initial radius with driving pressure when the distance between bubbles is 20 mm,the initial radius of bubble 2 is 5μm and the spherical asymmetric perturbation parameterδ p2is 0Pa,1×10?5Pa,1×10?7Pa,respectively.

    From Figs.6(a)and 6(b),we can observe that,under the same perturbation parameter,with the initial radius of bubble 2 in the two-bubble model being constant, the instability region of bubble 1 is larger compared to that of a single bubble.This illustrates that perturbation parameter has a substantial influence on the shape instability of bubbles.

    Comparing Fig.3(c) with Fig.6, we can observe that in the presence of perturbation parameter in the two-bubble model, larger perturbation parameters lead to an increase in the instability region of bubble 1.

    Fig.6.Shape instability and the diffusive equilibrium curves of bubble 1 in R0–Pa phase diagrams,for d=20 mm,and the spherical asymmetric perturbation parameter: (a) the single bubble model, δ p2 =1×10?5Pa, (b) the two-nonspherical-bubble model, R20 =5μm,δ p2=1×10?5Pa,(c)the two-nonspherical-bubble model,R20=5μm,δ p2=1×10?7Pa.

    Conversely,smaller perturbation parameter tends to facilitate the formation of stable oscillations, resulting in a reduction of the instability region.However,it is important to note that the relationship between perturbation parameter and shape instability is likely nonlinear,and the specific mechanisms involved require further investigation.

    From Figs.3(c), 6(b), and 6(c), we selected the region boundary as 1.3 atm and calculated the secondary Bjerknes force coefficientfBin theR10–R20plane corresponding to the three figures, as shown in Fig.7.The perturbation parameter has minimal impact on the overall trend of the secondary Bjerknes force.It only induces slight variations at the boundaries between the repulsion and attraction regions in the twononspherical-bubble model.

    Fig.7.The secondary Bjerknes force coefficient fB in the R10–R20 plane for driving presure Pa=1.3 atm with different spherical asymmetric perturbation parameter.The planes are given for perturbation: (a)δ p2=0,(b)δ p2=1×10?5Pa,and(c)δ p2=1×10?7Pa.

    4.Conclusion

    We have investigated the instability characteristics of two-nonspherical-gas-bubble cavitation and explored the influence of various factors.Using numerical simulations and the Keller–Miksis equation, we establish a theoretical model considering initial bubble radius, distance, and perturbation parameter.With this model,we find the unique stability characteristics arising from interactions between two bubbles as opposed to a single bubble under similar cavitation conditions.The pivotal role played by the secondary Bjerknes force and resonance frequencies is highlighted.Furthermore,alterations in the initial distance of two bubbles and perturbation parameters exert substantial influence on the extent of the instability region.Notably,the nonlinear nature of the relationship between perturbation parameters and shape stability calls for further in-depth investigation to elucidate its intricacies.These findings contribute to understanding two-nonsphericalgas-bubble cavitation and provide a foundation for achieving stable bubble sonoluminescence.

    Acknowledgment

    Project supported by the Scientific Research Project of Higher Education in the Inner Mongolia Autonomous Region(Grant No.NJZY23100).

    猜你喜歡
    烏日
    一葉輕舟承載祖孫情
    馬頭琴奏出世界風:難忘父母甩賣了半個牛群
    蘇日娜、周念祺、于婧、烏日罕作品
    烏日根 對表演的熱忱從未改變
    時尚北京(2020年11期)2020-11-16 02:08:18
    烏日嘎和他的烏潤合爾
    滿族文學(2020年4期)2020-09-03 04:29:23
    烏日嘎和他的烏潤合爾
    烏日嘎的蒙古馬
    安徽文學(2020年1期)2020-01-15 04:27:42
    烏日更達賴:播綠還“心債”
    額爾敦-烏日勒對AS家兔肝低密度脂蛋白受體(RDLR)表達的影響
    蒙藥薩烏日勒為主方治療薩病臨床療效分析
    精品久久国产蜜桃| 人体艺术视频欧美日本| 美女高潮的动态| 偷拍熟女少妇极品色| 爱豆传媒免费全集在线观看| 久久久亚洲精品成人影院| 亚洲精品成人久久久久久| 又粗又硬又长又爽又黄的视频| 亚洲最大成人av| 一二三四中文在线观看免费高清| 精品久久久久久久久av| 日韩强制内射视频| 卡戴珊不雅视频在线播放| 免费搜索国产男女视频| 亚洲五月天丁香| 久久精品夜色国产| 国产精品一二三区在线看| 亚洲av不卡在线观看| 久久久久久久久久久丰满| 欧美性猛交黑人性爽| 男女国产视频网站| 国产亚洲最大av| 一本—道久久a久久精品蜜桃钙片 精品乱码久久久久久99久播 | 国产精品一区二区在线观看99 | 黄片无遮挡物在线观看| 丰满人妻一区二区三区视频av| 久久久久精品久久久久真实原创| 国产成人免费观看mmmm| 激情 狠狠 欧美| 亚洲国产精品成人综合色| 午夜激情福利司机影院| 亚洲av电影不卡..在线观看| 欧美一级a爱片免费观看看| 中文字幕久久专区| 成人美女网站在线观看视频| 美女被艹到高潮喷水动态| 尤物成人国产欧美一区二区三区| 成人三级黄色视频| 国产三级中文精品| 不卡视频在线观看欧美| 欧美成人a在线观看| 久久国内精品自在自线图片| 99热6这里只有精品| 久久久久九九精品影院| eeuss影院久久| АⅤ资源中文在线天堂| 一二三四中文在线观看免费高清| 久久国内精品自在自线图片| 国产精品久久视频播放| 人人妻人人澡人人爽人人夜夜 | 97在线视频观看| 七月丁香在线播放| 欧美一级a爱片免费观看看| 免费黄网站久久成人精品| av线在线观看网站| 午夜福利成人在线免费观看| 在线播放国产精品三级| 最近2019中文字幕mv第一页| 97超视频在线观看视频| 99在线人妻在线中文字幕| 久久6这里有精品| 国产精品一区www在线观看| 一级av片app| 成人综合一区亚洲| 特级一级黄色大片| 国产精品久久久久久久电影| 美女大奶头视频| 色综合亚洲欧美另类图片| 国产一区二区亚洲精品在线观看| 国产成人免费观看mmmm| 狂野欧美激情性xxxx在线观看| 国产精华一区二区三区| 国产精品蜜桃在线观看| 夜夜爽夜夜爽视频| 国产av码专区亚洲av| 91av网一区二区| 国产精品一区二区三区四区久久| 18禁在线无遮挡免费观看视频| 精品无人区乱码1区二区| 久久精品国产鲁丝片午夜精品| 日本一本二区三区精品| 亚洲伊人久久精品综合 | 国产视频首页在线观看| 最近视频中文字幕2019在线8| 狠狠狠狠99中文字幕| 狠狠狠狠99中文字幕| 在现免费观看毛片| 中文字幕免费在线视频6| 天堂av国产一区二区熟女人妻| 波多野结衣巨乳人妻| 夜夜看夜夜爽夜夜摸| 亚洲精品乱码久久久久久按摩| 亚洲精品自拍成人| 国产女主播在线喷水免费视频网站 | 熟妇人妻久久中文字幕3abv| 国语自产精品视频在线第100页| 免费无遮挡裸体视频| 午夜福利网站1000一区二区三区| 一区二区三区高清视频在线| 中文资源天堂在线| 欧美+日韩+精品| 插逼视频在线观看| 蜜臀久久99精品久久宅男| 99久久精品热视频| 高清av免费在线| 狠狠狠狠99中文字幕| 欧美xxxx性猛交bbbb| 日本三级黄在线观看| 免费人成在线观看视频色| 国产精品99久久久久久久久| 亚洲中文字幕日韩| 免费无遮挡裸体视频| 国产精品蜜桃在线观看| 国产免费又黄又爽又色| 精品酒店卫生间| 久久精品人妻少妇| 我的女老师完整版在线观看| 搡老妇女老女人老熟妇| 最近的中文字幕免费完整| 国产高清三级在线| 中文字幕av在线有码专区| 黄色配什么色好看| 日本免费在线观看一区| 男的添女的下面高潮视频| 日韩中字成人| 中文字幕熟女人妻在线| 亚洲国产精品成人综合色| 免费黄色在线免费观看| 内地一区二区视频在线| 日韩在线高清观看一区二区三区| 国产精品.久久久| 91久久精品电影网| 国产成人精品婷婷| 免费搜索国产男女视频| 国产单亲对白刺激| 看片在线看免费视频| 欧美激情在线99| 欧美成人a在线观看| 久久99热6这里只有精品| 91av网一区二区| 欧美日本亚洲视频在线播放| 2022亚洲国产成人精品| 午夜视频国产福利| 亚洲电影在线观看av| 搡女人真爽免费视频火全软件| 一级毛片电影观看 | 国产日韩欧美在线精品| 两个人的视频大全免费| 亚洲自偷自拍三级| 校园人妻丝袜中文字幕| 在线免费观看不下载黄p国产| 精品久久久久久久久久久久久| av在线亚洲专区| 亚洲国产欧美在线一区| 丰满乱子伦码专区| 女人十人毛片免费观看3o分钟| 一边摸一边抽搐一进一小说| 22中文网久久字幕| 五月伊人婷婷丁香| 最近视频中文字幕2019在线8| 日本欧美国产在线视频| 亚洲国产色片| 色综合亚洲欧美另类图片| 一级毛片电影观看 | 成人鲁丝片一二三区免费| 黄色欧美视频在线观看| 日日干狠狠操夜夜爽| 黄片无遮挡物在线观看| 人妻系列 视频| 久久久国产成人精品二区| 观看免费一级毛片| 伊人久久精品亚洲午夜| 国产三级中文精品| 欧美一区二区国产精品久久精品| 亚洲国产精品成人综合色| 99热6这里只有精品| 又黄又爽又刺激的免费视频.| 亚洲在久久综合| 亚洲欧美精品专区久久| 久久精品夜夜夜夜夜久久蜜豆| 久久99热这里只频精品6学生 | 久久久久久国产a免费观看| 亚洲高清免费不卡视频| 一卡2卡三卡四卡精品乱码亚洲| 91狼人影院| 你懂的网址亚洲精品在线观看 | 少妇猛男粗大的猛烈进出视频 | 最近中文字幕2019免费版| 大话2 男鬼变身卡| 高清午夜精品一区二区三区| 亚洲一级一片aⅴ在线观看| 亚洲av电影在线观看一区二区三区 | 亚洲av成人精品一二三区| 高清午夜精品一区二区三区| 午夜日本视频在线| 色吧在线观看| 一二三四中文在线观看免费高清| 免费在线观看成人毛片| 秋霞伦理黄片| 夫妻性生交免费视频一级片| 汤姆久久久久久久影院中文字幕 | 国产黄片美女视频| 人人妻人人澡欧美一区二区| 黄片无遮挡物在线观看| 亚洲熟妇中文字幕五十中出| 午夜亚洲福利在线播放| 好男人在线观看高清免费视频| 午夜福利在线观看吧| 中文亚洲av片在线观看爽| 亚洲欧美精品自产自拍| 亚洲精品色激情综合| 亚洲综合精品二区| 精品久久久久久电影网 | 成人无遮挡网站| 免费看av在线观看网站| 欧美变态另类bdsm刘玥| 在线观看66精品国产| 亚洲精品aⅴ在线观看| 床上黄色一级片| 国产中年淑女户外野战色| 免费一级毛片在线播放高清视频| 熟女电影av网| 亚洲欧洲国产日韩| 少妇猛男粗大的猛烈进出视频 | 天天躁夜夜躁狠狠久久av| 久久久久久久久久黄片| 插阴视频在线观看视频| 久久草成人影院| 中文亚洲av片在线观看爽| 哪个播放器可以免费观看大片| 国产成人精品婷婷| 一区二区三区乱码不卡18| 日本三级黄在线观看| 大香蕉久久网| av国产免费在线观看| 99九九线精品视频在线观看视频| 久久久国产成人精品二区| 久久久欧美国产精品| 亚洲精品自拍成人| 日本猛色少妇xxxxx猛交久久| 建设人人有责人人尽责人人享有的 | 亚州av有码| 最近手机中文字幕大全| 直男gayav资源| 少妇熟女aⅴ在线视频| 亚洲精品乱久久久久久| 国产精品久久久久久精品电影小说 | 秋霞在线观看毛片| 欧美另类亚洲清纯唯美| 深爱激情五月婷婷| 国产69精品久久久久777片| 亚洲,欧美,日韩| 婷婷色综合大香蕉| 日本免费a在线| 欧美高清性xxxxhd video| 身体一侧抽搐| 一边亲一边摸免费视频| 亚洲欧美一区二区三区国产| a级一级毛片免费在线观看| 欧美一区二区精品小视频在线| 日韩三级伦理在线观看| 女人十人毛片免费观看3o分钟| 婷婷六月久久综合丁香| 久久精品久久精品一区二区三区| 18禁动态无遮挡网站| 网址你懂的国产日韩在线| www.色视频.com| 午夜免费激情av| 18禁动态无遮挡网站| 女人久久www免费人成看片 | 性色avwww在线观看| 在线播放国产精品三级| 麻豆久久精品国产亚洲av| 日日撸夜夜添| av视频在线观看入口| 国产一级毛片在线| 成人一区二区视频在线观看| 夜夜看夜夜爽夜夜摸| 亚洲怡红院男人天堂| 国产一区二区亚洲精品在线观看| 有码 亚洲区| 久久久久久大精品| 三级国产精品片| 国产精品国产高清国产av| 色综合亚洲欧美另类图片| 午夜激情欧美在线| 国产精品一区二区在线观看99 | 2022亚洲国产成人精品| 天天躁日日操中文字幕| 99热精品在线国产| 日韩欧美精品v在线| 丰满少妇做爰视频| 美女被艹到高潮喷水动态| 男女下面进入的视频免费午夜| 在线免费十八禁| 国语自产精品视频在线第100页| 日韩欧美三级三区| 男人和女人高潮做爰伦理| 麻豆成人av视频| 亚洲怡红院男人天堂| 国产淫片久久久久久久久| 色5月婷婷丁香| 日韩精品有码人妻一区| 毛片一级片免费看久久久久| 精华霜和精华液先用哪个| 中文资源天堂在线| 日韩人妻高清精品专区| 免费av毛片视频| 久久久久网色| 色5月婷婷丁香| 亚洲激情五月婷婷啪啪| 日韩强制内射视频| 69av精品久久久久久| 亚洲欧美精品自产自拍| 亚洲欧美日韩东京热| av在线老鸭窝| 久久久久久久亚洲中文字幕| 国产亚洲精品av在线| 高清午夜精品一区二区三区| av免费在线看不卡| 一个人观看的视频www高清免费观看| 久久精品国产99精品国产亚洲性色| 国产伦精品一区二区三区四那| 精品人妻偷拍中文字幕| 99九九线精品视频在线观看视频| 欧美成人精品欧美一级黄| 国产探花极品一区二区| 哪个播放器可以免费观看大片| 淫秽高清视频在线观看| 麻豆乱淫一区二区| 色5月婷婷丁香| 亚洲欧美日韩无卡精品| 亚洲欧美日韩卡通动漫| 国产在线一区二区三区精 | 日韩av在线免费看完整版不卡| 又粗又硬又长又爽又黄的视频| 狂野欧美白嫩少妇大欣赏| av女优亚洲男人天堂| 看免费成人av毛片| 级片在线观看| 国产一级毛片在线| 真实男女啪啪啪动态图| 欧美潮喷喷水| 欧美日韩一区二区视频在线观看视频在线 | 国产三级中文精品| 日韩大片免费观看网站 | 1024手机看黄色片| 99久久人妻综合| 国产精品人妻久久久久久| 精品熟女少妇av免费看| 少妇丰满av| 欧美成人a在线观看| 国产爱豆传媒在线观看| 国产成人精品婷婷| 日日干狠狠操夜夜爽| 国产精品久久久久久精品电影小说 | 久久精品夜夜夜夜夜久久蜜豆| 床上黄色一级片| 最近手机中文字幕大全| 免费黄网站久久成人精品| 国产免费男女视频| 成人毛片60女人毛片免费| 精品久久久久久电影网 | 麻豆一二三区av精品| 国产免费福利视频在线观看| 久久99蜜桃精品久久| 亚洲国产欧美在线一区| 一级毛片我不卡| 欧美极品一区二区三区四区| 伦精品一区二区三区| 中文欧美无线码| 色综合站精品国产| 国产淫片久久久久久久久| 一级爰片在线观看| 国产精品一区二区在线观看99 | 久久久久久久久中文| 亚洲第一区二区三区不卡| 精品一区二区免费观看| 七月丁香在线播放| 舔av片在线| 两个人视频免费观看高清| 禁无遮挡网站| 国产亚洲av片在线观看秒播厂 | 99热这里只有是精品在线观看| 免费人成在线观看视频色| 国产精品一区二区在线观看99 | 久久精品国产亚洲av涩爱| 天堂中文最新版在线下载 | 我要看日韩黄色一级片| 99久国产av精品国产电影| av线在线观看网站| 少妇熟女aⅴ在线视频| 欧美又色又爽又黄视频| 精品久久久久久久久av| 大香蕉97超碰在线| videos熟女内射| 中文亚洲av片在线观看爽| 亚洲怡红院男人天堂| 卡戴珊不雅视频在线播放| 欧美成人精品欧美一级黄| 国产精品乱码一区二三区的特点| 蜜臀久久99精品久久宅男| 九九热线精品视视频播放| 亚洲欧美中文字幕日韩二区| 国产探花在线观看一区二区| 日本免费在线观看一区| 国产精品一及| 亚洲va在线va天堂va国产| 久久久久久久久久久丰满| h日本视频在线播放| 在线观看一区二区三区| 99国产精品一区二区蜜桃av| 亚洲欧美精品综合久久99| 国产伦一二天堂av在线观看| 2021天堂中文幕一二区在线观| 狂野欧美激情性xxxx在线观看| 国产亚洲91精品色在线| 国产白丝娇喘喷水9色精品| 国产一区二区在线观看日韩| 一区二区三区四区激情视频| 蜜桃亚洲精品一区二区三区| 天堂中文最新版在线下载 | 亚洲国产欧美人成| 熟女人妻精品中文字幕| 免费人成在线观看视频色| 神马国产精品三级电影在线观看| 最近视频中文字幕2019在线8| 国产美女午夜福利| 联通29元200g的流量卡| 国产高潮美女av| 久久精品综合一区二区三区| 日韩一区二区三区影片| av在线亚洲专区| 免费看a级黄色片| 亚洲av一区综合| 女人被狂操c到高潮| 精品酒店卫生间| 美女国产视频在线观看| 亚洲精品国产av成人精品| 国产色婷婷99| 国内少妇人妻偷人精品xxx网站| 综合色av麻豆| 啦啦啦啦在线视频资源| 国产毛片a区久久久久| 99久久精品热视频| 春色校园在线视频观看| 国产成人精品一,二区| 纵有疾风起免费观看全集完整版 | 精品国产露脸久久av麻豆 | 亚洲色图av天堂| 日韩欧美三级三区| 国产人妻一区二区三区在| 久久这里有精品视频免费| 欧美日本亚洲视频在线播放| 色播亚洲综合网| 久久99热6这里只有精品| 精品国内亚洲2022精品成人| 永久网站在线| 少妇被粗大猛烈的视频| 国产色婷婷99| 欧美又色又爽又黄视频| 男人的好看免费观看在线视频| 最近中文字幕2019免费版| 欧美一区二区亚洲| 夜夜爽夜夜爽视频| 在线观看一区二区三区| 99热这里只有是精品在线观看| 精品人妻偷拍中文字幕| 亚洲综合精品二区| 国产伦一二天堂av在线观看| 亚洲自偷自拍三级| 99久久人妻综合| 日韩制服骚丝袜av| 国产单亲对白刺激| 精品酒店卫生间| 中国国产av一级| 成人综合一区亚洲| 亚洲av免费在线观看| 国产精品野战在线观看| 国产av一区在线观看免费| 在线观看av片永久免费下载| 成人漫画全彩无遮挡| 99久久人妻综合| 国产乱人视频| 亚洲婷婷狠狠爱综合网| 日韩视频在线欧美| 草草在线视频免费看| 亚洲欧美日韩高清专用| 日韩国内少妇激情av| 高清日韩中文字幕在线| 成人鲁丝片一二三区免费| 黄色日韩在线| 高清av免费在线| 久久久久久久国产电影| av.在线天堂| 69av精品久久久久久| 又黄又爽又刺激的免费视频.| 日韩欧美精品v在线| 免费搜索国产男女视频| 亚洲在线观看片| 精品人妻视频免费看| 国语对白做爰xxxⅹ性视频网站| 国产精品乱码一区二三区的特点| 我的女老师完整版在线观看| 一级av片app| 老司机影院毛片| 欧美性猛交╳xxx乱大交人| 高清在线视频一区二区三区 | 少妇丰满av| 亚洲av不卡在线观看| 精华霜和精华液先用哪个| 在线免费观看不下载黄p国产| 插逼视频在线观看| 最近视频中文字幕2019在线8| 十八禁国产超污无遮挡网站| 91午夜精品亚洲一区二区三区| 日本一二三区视频观看| 中文字幕av在线有码专区| 美女黄网站色视频| 干丝袜人妻中文字幕| 国产免费福利视频在线观看| 99久久无色码亚洲精品果冻| 欧美日本视频| 舔av片在线| 国产成人一区二区在线| 日本爱情动作片www.在线观看| 天天躁日日操中文字幕| 精品久久久噜噜| 免费不卡的大黄色大毛片视频在线观看 | 一本—道久久a久久精品蜜桃钙片 精品乱码久久久久久99久播 | 国产精品久久久久久精品电影小说 | 99热精品在线国产| 欧美精品国产亚洲| 又黄又爽又刺激的免费视频.| 少妇的逼好多水| 成人午夜精彩视频在线观看| 国产视频首页在线观看| 尾随美女入室| 免费观看人在逋| h日本视频在线播放| 亚洲精品国产av成人精品| 一区二区三区高清视频在线| 18禁在线无遮挡免费观看视频| 亚洲av二区三区四区| 一个人观看的视频www高清免费观看| 国产激情偷乱视频一区二区| 国产伦在线观看视频一区| 久久精品综合一区二区三区| 欧美变态另类bdsm刘玥| 色综合亚洲欧美另类图片| 女的被弄到高潮叫床怎么办| 亚洲精品日韩av片在线观看| 大香蕉97超碰在线| 你懂的网址亚洲精品在线观看 | 久久久久久久国产电影| 亚洲精品日韩在线中文字幕| 午夜精品在线福利| 高清午夜精品一区二区三区| 日本免费a在线| 免费av观看视频| 日韩三级伦理在线观看| 亚洲天堂国产精品一区在线| 免费看a级黄色片| 99热这里只有是精品在线观看| 精品一区二区三区人妻视频| 禁无遮挡网站| 亚洲欧美日韩东京热| 亚洲av电影在线观看一区二区三区 | a级一级毛片免费在线观看| 美女黄网站色视频| 成年av动漫网址| 国产午夜精品久久久久久一区二区三区| 国产精品久久久久久精品电影小说 | 三级毛片av免费| 美女脱内裤让男人舔精品视频| 国产男人的电影天堂91| 91精品国产九色| 国产爱豆传媒在线观看| 中文乱码字字幕精品一区二区三区 | 长腿黑丝高跟| 国产大屁股一区二区在线视频| 人体艺术视频欧美日本| 91在线精品国自产拍蜜月| 国产成人午夜福利电影在线观看| 日本猛色少妇xxxxx猛交久久| 中文在线观看免费www的网站| 精品久久久久久电影网 | 国产三级中文精品| 国内精品宾馆在线| 一级黄色大片毛片| av国产久精品久网站免费入址| 国内精品一区二区在线观看| 一区二区三区四区激情视频| 午夜免费男女啪啪视频观看| 精品久久国产蜜桃| 久久久久久久亚洲中文字幕| 国产成年人精品一区二区| 蜜桃久久精品国产亚洲av| 久久久久精品久久久久真实原创| 麻豆精品久久久久久蜜桃| 国产av在哪里看| 最近的中文字幕免费完整| 亚洲精品乱码久久久久久按摩| 国产精品熟女久久久久浪| 1000部很黄的大片| 久久久久九九精品影院| 大话2 男鬼变身卡| 99久久精品国产国产毛片| 大话2 男鬼变身卡| 2021少妇久久久久久久久久久| 寂寞人妻少妇视频99o| 亚洲av免费高清在线观看| 亚洲中文字幕一区二区三区有码在线看| 国产精品久久久久久av不卡| 久久久久久久久久久免费av|