首页 | 本学科首页   官方微博 | 高级检索  
相似文献
 共查询到20条相似文献,搜索用时 46 毫秒
1.
The rheological behaviors of a compound droplet in a confined geometry are of importance in many industrial and natural processes. However, a detailed numerical simulation of the finite deformation and its transition to the breakup of the multi-core compound droplet in an axisymmetric T-junction channel is still lacking. The present study is to fill this gap through the numerical simulations of a two-core compound droplet that deforms and breaks up in this channel configuration. The numerical results are obtained by the axisymmetric front-tracking method. Our new finding is that the compound droplet in the channel can experience the finite deformation or the breakup depending on the flow condition or the configuration of the channel. In the finite deformation mode (i.e. non-breakup mode), the droplet rapidly reaches the maximum deformation before approaching the perpendicular rigid wall. The most deformation occurs with the outer droplet, and the inner droplet closer to the wall is less deformed than the other inner core droplet. In the breakup mode, three breakup patterns are recognized: (i) breakup type I occurring when breaking up only the outer droplet; (ii) breakup type 2 occurring when breaking up only the inner droplets; (iii) breakup type 3 occurring when breaking up both inner and outer droplets. The transition from the non-breakup mode to the breakup mode is available when increasing the Reynolds number Re (from 0.16 to 40.0), the capillary number Ca (from 0.04 to 4.0), the size Ro of the outer droplet and the middle-to-outer fluid viscosity ratio μ21, or decreasing the size Ri of the inner droplets, the radial size C2 of the channel (normalized by the channel axial size C1) and the interfacial tension ratio of the inner to the outer droplets. The transition diagrams based on some of these parameters are also proposed to provide a more complete picture of the two-core compound droplet behaving in the axisymmetric T-junction channel..  相似文献   

2.
The determination of the critical Weg number separating the different breakup regimes has been extensively studied in several experimental and numerical works, while empirical and semi-analytical approaches have been proposed to relate the critical Weg number with the Ohl number. Nevertheless, under certain conditions, the Reg number and the density ratio ε may become important. The present work provides a simple but reliable enough methodology to determine the critical Weg number as a function of the aforementioned parameters in an effort to fill this gap in knowledge. It considers the main forces acting on the droplet (aerodynamic, surface tension and viscous) and provides a general criterion for breakup to occur but also for the transition among the different breakup regimes. In this light, the present work proposes the introduction of a new set of parameters named as Weg,eff and Cal monitored in a new breakup plane. This plane provides a direct relation between gas inertia and liquid viscosity forces, while the secondary effects of Reg number and density ratio have been embedded inside the effective Weg number (Weg,eff)  相似文献   

3.
The breakup mechanism and instability of a power law liquid jet are investigated in this study. The power law model is used to account for the non-Newtonian behavior of the liquid fluid. A new theoretical model is established to explain the breakup of a power law liquid jet with axisymmetric and asymmetric disturbances, which moves in a swirling gas. The corresponding dispersion relation is derived by a linear approximation, and it is applicable for both shear-thinning and shear-thickening liquid jets. Analysis results are calculated based on the temporal mode. The analysis includes the effects of the generalized Reynolds number, the Weber number, the power law exponent, and the air swirl strength on the breakup of the jet. Results show that the shear-thickening liquid jet is more unstable than its Newtonian and shear-thinning counterparts when the effect of the air swirl is taken into account. The axisymmetric mode can be the dominant mode on the power law jet breakup when the air swirl strength is strong enough, while the non-axisymmetric mode is the domination on the instability of the power liquid jet with a high We and a low Re n . It is also found that the air swirl is a stabilizing factor on the breakup of the power law liquid jet. Furthermore, the instability characteristics are different for different power law exponents. The amplitude of the power law liquid jet surface on the temporal mode is also discussed under different air swirl strengths.  相似文献   

4.
Compound droplet collision has been found in various industrial and academic applications. Such colliding phenomena of two droplets in simple shear flow are numerically resolved by a front-tracking technique. Initially, each compound droplet, assumed cylindrical, contains one concentric inner droplet. They are separated at the lateral and vertical intervals denoted by Δx0 and Δy0. Because of the shear flow, the compound droplets interact with each other and exhibit three collision modes: passing-over, merging (i.e. coalescence) and reversing. These modes and their transition are affected by many parameters including the Reynolds Re and Capillary Ca numbers (based on the properties of the outer fluid), the viscosity ratios μ13 and μ23, the interfacial tension ratio σ12 of the inner to outer interfaces, the ratio of the radii of the inner to outer droplets R12 and the initial distance between them (in terms of Δx0 and Δy0). It is found that from a merging mode, decreasing Re from 2.51 to a value less than or equal to 1.0 induces a transition to a reversing mode, whereas, increasing Ca from 0.005 to a value greater than or equal to 0.04 leads to a transition to a passing-over mode. A transition from the merging to passing-over modes also appears when varying μ23 in the range of 0.1–10.0 or varying R12 in the range of 0.2–0.8. A transition from a passing-over mode to a reversing one is available when increasing Δx0 or decreasing Δy0. Three modes of collision all occur when μ13 is varied in the range of 0.1–10.0. However, the variation of σ12 does not induce any transition between different modes. Several phase diagrams in terms of Re versus Ca, or Δx0 versus Δy0 are also proposed to show the transitions between these modes.  相似文献   

5.
The dynamics of compound droplets is more and more attractive because of their applications in a wide range of industrial and natural processes. This study aims to improve the understanding of dynamical rheological behaviors of a compound droplet moving in a nozzle with a conical shape in the downstream region via front-tracking-based simulations. The numerical results show that the compound droplet experiences three stages of deformation: the entrance stage (in front of the conical region), the transit stage (within the conical region), and the exit stage (in the exit of the nozzle). The droplet receives the maximum deformation in the axial direction during the transit stage, and the radially maximum deformation occurs during the exit stage. Because of the acceleration induced by the conical region, the inner droplet of the compound droplet can break up into smaller droplets during the exit stage. To reveal the transition between the finite deformation and the breakup, many parameters including the Capillary number Ca (varied in the range of 0.0125–1.6), the droplet size relative to the nozzle size R1/R0 (varied in the range of 0.2–0.9), the droplet radius ratio R21 (varied in the range of 0.3–0.8), the viscosity ratios μ21 and μ31 (varied in the range of 0.05–3.2), the interfacial tension ratio σ21 (varied in the range of 0.125–8.0), the conical angle α (varied in the range of 4°–34°) and the initial location of the inner droplet (i.e. the droplet eccentricity) are considered. From the finite deformation mode, the transition to the breakup mode of the inner droplet occurs when increasing any of Ca, R1/R0, R21 and α, or decreasing any of μ21 and σ21. The breakup mode is also enhanced when the inner droplet is initially located closer to the leading side of the outer droplet. However, varying μ31 induces no transition between these modes. The regime diagrams of these modes, based on these parameters, are also proposed.  相似文献   

6.
Dynamical behaviors of compound droplets including those interacting with rigid walls in an axisymmetric channel appear in various industrial and natural processes. However, so far, no detailed investigation has been carried out for such interactions of compound droplets. Motivating from this missing gap, we here numerically study the finite deformation and breakup of an initially concentric compound droplet when it moves toward a rigid wall at the bottom of an axisymmetric vertical channel. The method used is a finite difference-based front-tracking method. The numerical results reveal that when the compound droplet is delivered toward the wall, it is deformed and can break up into smaller droplets. For the cases of finite deformation (i.e. non-breakup), while the outer droplet is radially stretched, the inner droplet first moves downward in the direction of the outer flow but then gets back. Thereby, a thin film is created between the outer and inner interfaces at the droplet top and thus prevents the outer droplet further deforming and breaking up. In contrast, if breakup happens, the outer droplet is further stretched, and most of the middle fluid moves outward toward the outer droplet edge to form a blob. Breakup can be available in one of three patterns: off-axis breakup, on-axis breakup, and inner breakup. The off-axis breakup mode only happens with the outer droplet while the inner breakup mode is only for the inner droplet. Various parameters are investigated to show the transition between a non-breakup mode to a mode of breakup. Such parameters contributing the transition include the Capillary number Ca (varied in the range of 0.01–2.5), the channel aspect ratio (varied in the range of 0.4–2.0), the ratio of the inner to outer droplet radii (varied in the range of 0.3–0.8), the droplet size relative to the channel size (varied in the range of 0.2–0.9), the interfacial tension ratio of the inner to outer interfaces (varied in the range of 0.1–4.0), and the viscosity ratio of the middle to outer fluids (varied in the range of 0.16–6.3). In contrast, some others, e.g. the Reynolds number, the viscosity ratio of the inner to the outer, do not induce any transition. From the numerical results, regime diagrams of breakup and non-breakup based on these parameters are proposed.  相似文献   

7.
In this paper we investigate experimentally the injection of a negatively buoyant jet into a homogenous immiscible ambient fluid. Experiments are carried out by injecting a jet of dyed fresh water through a nozzle in the base of a cylindrical tank containing rapeseed oil. The fountain inlet flow rate and nozzle diameter were varied to cover a wide range of Richardson Ri (8 × 10−4 < Ri < 1.98), Reynolds Re (467 < Re < 5,928) and Weber We (2.40 < We < 308.56) numbers. Based on the Re, Ri and We values for the experiments, we have determined a regime map to define how these values may control the occurrence of the observed flow types. Whereas Ri plays a stronger role when determining the maximum penetration height, the effect of the Reynolds number is stronger predicting the flow behaviour for a specific nozzle diameter and injection velocity.  相似文献   

8.
The paper is mainly focused to the vast number of researchers who work within direct injection (DI) engine fuel spray simulations. The most common simulation framework today within the community is the Reynolds Averaged Navier Stokes (RANS) approach together with the Lagrangian Particle Tracking (LPT) method. In fact, this study is one of the first studies where high resolution LES/LPT diesel spray modeling is considered. The potential of LES to deepen the present day multidimensional LPT fuel spray simulations is discussed. Spray evolution is studied far from an injector by modeling a spray as a particle laden jet (PLJ). The effect of d on mixing in non-atomizing and atomizing sprays is thoroughly investigated at jet inlet Reynolds number Re?=?104 and Mach number Ma?=?0.3. Based on and justified by rather recent and also quite old ideas, novel and compact views on droplet breakup in turbulent flows are pointed out from the literature. We use LES/LPT to illustrate that even in a low Weber number flow (We?<?13) the droplet breakup modeling may need considerable attention in contrast to what is typically assumed in the present-day breakup models. LES and LPT techniques are first applied to essentially confirm certain expected droplet size effects on spray shape in non-atomizing monodisperse sprays. In the simulations LES e.g. produces an expected turbulent dispersion pattern that depends on droplet diameter (d) without a droplet dispersion model in contrast to RANS. A new compact droplet breakup model is formulated and tested for droplets that break with a natural resonance time rate according to the Poisson process. As a result of the study: 1) the analysis gives a rigorous and enriching proof of currently existing views on droplet size effects on mixing, and 2) the presented a priori analysis points out the importance of modeling the resonance breakup even at a low We.  相似文献   

9.
The change in flow characteristics downstream of a circular cylinder (inner cylinder) surrounded by an outer permeable cylinder was investigated in shallow water using particle image velocimetry technique. The diameter of the inner cylinder and the water height were kept constant during the experiments as d?=?50?mm and h w ?=?25?mm, respectively. The depth-averaged free-stream velocity was also kept constant as U?=?170?mm/s which corresponded to a Reynolds number of Red?=?8,500 based on the inner cylinder diameter. In order to examine the effect of diameter and porosity of the outer cylinder on flow characteristics of the inner cylinder, five different outer cylinder diameters (D?=?60, 70, 80, 90 and 100?mm) and four different porosities (???=?0.4, 0.5, 0.6 and 0.7) were used. It was shown that both porosity and outer cylinder diameter had a substantial effect on the flow characteristics downstream of the circular cylinder. Turbulent statistics clearly demonstrated that in comparison with the bare cylinder (natural case), turbulent kinetic energy and Reynolds stresses decreased remarkably when an outer cylinder was placed around the inner cylinder. Thereby, the interaction of shear layers of the inner cylinder has been successfully prevented by the presence of outer cylinder. It was suggested by referring to the results that the outer cylinder having 1.6????D/d????2.0 and 0.4????D/d????0.6 should be preferred to have a better flow control in the near wake since the peak magnitude of turbulent kinetic energy was considerably low in comparison with the natural case and it was nearly constant for these mentioned porosities ??, and outer cylinder to inner cylinder diameter ratios D/d.  相似文献   

10.
Fully-developed turbulent flow in a concentric annulus, r1/r2 = 0.5, Reh = 12,500, with the outer wall rotating at a range of rotation rates N = Uθ,wall/Ub from 0.5 up to 4 is studied by large-eddy simulations. The focus is on the effects of moderate to very high rotation rates on the mean flow, turbulence statistics and eddy structure. For N up to ∼2, an increase in the rotation rate dampens progressively the turbulence near the rotating outer wall, while affecting only mildly the inner-wall region. At higher rotation rates this trend is reversed: for N = 2.8 close to the inner wall turbulence is dramatically reduced while the outer wall region remains turbulent with discernible helical vortices as the dominant turbulent structure. The turbulence parameters and eddy structures differ significantly for N = 2 and 2.8. This switch is attributed to the centrifuged turbulence (generated near the inner wall) prevailing over the axial inertial force as well as over the counteracting laminarizing effects of the rotating outer wall. At still higher rotation, N = 4, the flow gets laminarized but with distinct spiralling vortices akin to the Taylor–Couette rolls found between the two counter-rotating cylinders without axial flow, which is the limiting case when N approaches to infinity. The ratio of the centrifugal to axial inertial forces, Ta/Re2  N2 (where Ta is the Taylor number) is considered as a possible criterion for defining the conditions for the above regime change.  相似文献   

11.
We have recently discovered a new type of self-excited flapping jets due to a flexible film whose leading edge is fixed at the nozzle exit [Exp Ther Fluid Sci, 106, 226-233]. This paper is to report the experimental investigation on mixing characteristics of the jet induced by a rectangular FEP film. Hot wire anemometry and flow visualization are used to examine the flapping jet flow versus the non-flapping counterpart. Experiments are conducted under the following conditions: i.e., L/D = 1.0 (fixed), W/D = 0.03 ~ 1.0 (varying) and Re = 10000 ~ 45000 (varying); where W and L are the film's width and length, D is the nozzle-exit diameter, and Re is the Reynolds number defined by Re UoD/ν with Uo and ν being the jet-exit velocity and fluid viscosity.It is found that the jet-flapping frequency fF varies with W in a complex fashion while it grows roughly linearly with increasing Uo for W/D ≥ 0.5. The flapping Strouhal number StF fFD/Uo ranges in 0.13 ≤ StF ≤ 0.23 for Re = 15,000 ~ 45,000. These Strouhal numbers are substantially lower than that (≈ 0.45 ~ 0.7) for the primary vortex generation in the free jet, but one to two orders of magnitude higher than those from the conventional self-exciting fluidic devices. In general, the flapping jet decays and spreads more rapidly than does the free jet. As W increases, the decaying and spreading rates both grow. Of significance, the centerline evolutions of Taylor and Kolmogorov scales versus the integral scale are examined to characterize the small scales of turbulence against the large-scale motion.  相似文献   

12.
When a spherical liquid drop is subjected to a step change in relative gas velocity, aerodynamic forces lead to drop deformation and possible breakup into a number of secondary fragments. To investigate this flow, a digital in-line holography (DIH) diagnostic is proposed which enables rapid quantification of spatial statistics with limited experimental repetition. To overcome the high uncertainty in the depth direction experienced in previous applications of DIH, a crossed-beam, two-view configuration is introduced. With appropriate calibration, this diagnostic is shown to provide accurate quantification of fragment sizes, three-dimensional positions and three-component velocities in a large measurement volume. These capabilities are applied to investigate the aerodynamic breakup of drops at two non-dimensional Weber numbers, We, corresponding to the bag (We  = 14) and sheet-thinning (We = 55) regimes. Ensemble average results show the evolution of fragment size and velocity statistics during the course of breakup. Results indicate that mean fragment sizes increase throughout the course of breakup. For the bag breakup case, the evolution of a multi-mode fragment size probability density is observed. This is attributed to separate fragmentation mechanisms for the bag and rim structures. In contrast, for the sheet-thinning case, the fragment size probability density shows only one distinct peak indicating a single fragmentation mechanism. Compared to previous related investigations of this flow, many orders of magnitude more fragments are measured per condition, resulting in a significant improvement in data fidelity. For this reason, this experimental dataset is likely to provide new opportunities for detailed validation of analytic and computational models of this flow.  相似文献   

13.
Nonlinear instability and breakup of an annular liquid sheet has been modeled in this paper. The liquid sheet is considered to move axially and is exposed to co-flowing inner and outer gas streams. Also, the effect of outer gas swirl on sheet breakup has been studied. In the developed model a perturbation expansion method has been used with the initial magnitude of the disturbance as the perturbation parameter. This is a comprehensive model in that other geometries of planar sheet and a coaxial jet can be obtained as limiting cases of very large inner radius and inner radius equal to zero, respectively. In this temporal analysis, the effect of liquid Weber number, initial disturbance amplitude, inner gas-to-liquid velocity ratio, outer gas-to-liquid velocity ratio and outer gas swirl strength on the breakup time is investigated. The model is validated by comparison with earlier analytical studies for the limiting case of a planar sheet as well as with experimental data of sheet breakup length available in literature. It is shown that the linear theory cannot predict breakup of an annular sheet and the developed nonlinear model is necessary to accurately determine the breakup length. In the limiting case of a coaxial jet, results show that gas swirl destabilizes the jet, makes helical modes dominant compared to the axisymmetric mode and decreases jet breakup length. These results contradict earlier linear analyses and agree with experimental observations. For an annular sheet, it is found that gas flow hastens the sheet breakup process and shorter breakup lengths are obtained by increasing the inner and the outer gas velocity. Axially moving inner gas stream is more effective in disintegrating the annular sheet compared to axially moving outer gas stream. When both gas streams are moving axially, the liquid sheet breakup is quicker compared to that with any one gas stream. In the absence of outer gas swirl, the axisymmetric mode is the dominant instability mode. However, when outer gas flow has a swirl component higher helical modes become dominant. With increasing outer gas swirl strength, the maximum disturbance growth rate increases and the most unstable circumferential wave number increases resulting in a highly asymmetric sheet breakup with shorter breakup lengths and thinner ligaments.  相似文献   

14.
The vaporization characteristics of a liquid hexanes jet in a lab-scale test section with a plain orifice-type injector were experimentally investigated. The experimental measurements were carried out on the basis of the infrared laser extinction method using two He–Ne lasers (one at 632.8 nm and the other at 3.39 μm). The momentum flux ratio (qF/A) was varied from 20 to 60 over 20 steps, and the supplying air temperature (TA) was changed from 20 to 260 °C over 120 steps. The objectives of the current study were to assess the vaporization characteristics of a liquid hexanes jet and to derive a correlation between flow conditions and hexanes vapor concentration in a jet-in-crossflow configuration. From the results of the experimental measurement, it was concluded that hexanes vapor concentration increased with the increase of the momentum flux ratio and the supplying air temperature. An experimental correlation between flow conditions and hexanes vapor concentration (ZF) was proposed as a function of the normalized horizontal distance (x/do), the supplying air temperature (TA), the momentum flux ratio (qF/A), the fuel jet Reynolds number (ReF), and the fuel jet Weber number (WeF).  相似文献   

15.
为研究射流在水动力作用下的碎化特性,采用有限体积法对轴对称坐标下Navier-Stokes方程进行了求解,考虑重力和表面张力的影响,并通过Volume-of-Fluid法与Level-Set法成功捕捉到界面的不稳定发展、变形及射流碎化过程,分析了流场内部速度场和压力场分布,结果表明,射流碎化长度随Re/We“5数呈指数型增加,最后探讨了射流速度、直径及周围流体密度、粘性等参量对射流的碎化过程的影响规律.  相似文献   

16.
17.
The paper presents average flow visualizations and measurements, obtained with the Particle Image Velocimetry (PIV) technique, of a submerged rectangular free jet of air in the range of Reynolds numbers from Re = 35,300 to Re = 2200, where the Reynolds number is defined according to the hydraulic diameter of a rectangular slot of height H. According to the literature, just after the exit of the jet there is a zone of flow, called zone of flow establishment, containing the region of mixing fluid, at the border with the stagnant fluid, and the potential core, where velocity on the centerline maintains a value almost equal to the exit one. After this zone is present the zone of established flow or fully developed region. The goal of the paper is to show, with average PIV visualizations and measurements, that, before the zone of flow establishment is present a region of flow, never mentioned by the literature and called undisturbed region of flow, with a length, LU, which decreases with the increase of the Reynolds number. The main characteristics of the undisturbed region is the fact that the velocity profile maintains almost equal to the exit one, and can also be identified by a constant height of the average PIV visualizations, with length, LCH, or by a constant turbulence on the centerline, with length LCT. The average PIV velocity and turbulence measurements are compared to those performed with the Hot Film Anemometry (HFA) technique. The average PIV visualizations show that the region of constant height has a length LCH which increases from LCH = H at Re = 35,300 to LCH = 45H at Re = 2200. The PIV measurements on the centerline of the jet show that turbulence remains constant at the level of the exit for a length, LCT, which increases from LCT = H at Re = 35,300 to LCT = 45H at Re = 2200. The PIV measurements show that velocity remains constant at the exit level for a length, LU, which increases from LU = H at Re = 35,300 to LU = 6H at Re = 2200 and is called undisturbed region of flow. In turbulent flow the length LU is almost equal to the lengths of the regions of constant height, LCH, and constant turbulence, LCT. In laminar flow, Re = 2200, the length of the undisturbed region of flow, LU, is greater than the lengths of the regions of constant height and turbulence, LCT = LCH = 45H. The average PIV and HFA velocity measurements confirm that the length of potential core, LP, increases from LP = 45H at Re = 35,300 to LP = 78H at Re = 2200, and are compared to the previous experimental and theoretical results of the literature in the zone of mixing fluid and in the fully developed region with a good agreement.  相似文献   

18.
The interaction of a particle-laden jet (ReUe=6600)(ReUe=6600) with a single synthetic jet or a continuous control jet located upstream of the main jet exit (i.e., within the main jet nozzle) was examined experimentally using PIV and PTV. A reduction technique was used to calculate 3-D flow fields from multiple 2-D measurement planes to study the complex 3-D interactions. The synthetic jet was shown to influence the particles both directly and indirectly through the manipulation of the carrier fluid’s drag force on the particles. The synthetic jet impulse directly vectors the particles away from the synthetic jet, while the formation of large vortical structures indirectly affects the particles, spreading throughout the measurement domain. By comparison, a continuous control jet only vectors the particles away from itself. The lowest Stokes number particles respond similarly to the carrier fluid, while higher Stokes number particles are less responsive to the control and only follow the strong vortical structures (i.e., higher circulation), which suggests that the preferential concentration concept depends on both the Stokes number as well as the strength of the coherent structures.  相似文献   

19.
The qualitative and quantitative behaviour of double row vortical structures in the near field region of a plane wall jet are studied experimentally by flow visualization and hot-wire measurements. Ensemble averaging is employed to investigate the interaction of vortices with the wall. In the flow visualization study, a double row vortical structure, which includes a primary vortex formed in the outer layer region and a secondary vortex induced in the inner layer region, and the vortex lift-off phenomenon are clearly observed during the development of the wall jet. The phase averaged results of the velocity measurements show that the instability leading to induction of the secondary vortex is stimulated by the primary vortex. In the early stage of wall jet transition, the inflection point of the inner layer velocity profile moves transversely from the wall surface to the inner layer region due to passage of the well-organized primary vortex in the outer layer region. The inner layer instability is thus induced and the instability wave rolls up to form the secondary vortex. Furthermore, the secondary vortex will convect downstream faster than the primary vortex, and this difference in convective speed will lead to the subsequent phenomenon of vortex lift-off from the wall surface.List of symbols A1, A2, . . . primary vortex - B1,B2, . . . secondary vortex - fe forcing frequency - f fundamental frequency - H nozzle exit height - Re Reynolds number,U j H/ - T period of the referred signal (=13.5 ms) - t, t time scale - U streamwise mean velocity - U c convection speed - U j jet exit velocity - U m local maximum velocity - ut' streamwise turbulence intensity - uv turbulent shear stress - V transverse mean velocity - v transverse turbulence intensity - X streamwise coordinate - Y transverse coordinate - X Ai streamwise location of vortexAi - X Bi streamwise location of vortexBi - X ave averaged streamwise location of the vortex - Y m wall jet inner layer width, the distance from wall to whereU=U m - Y 1/2 wall jet half-width, the distance from wall to whereU=1/2U m in outer layer region - t time interval (=0.267 s) - phase averaged value  相似文献   

20.
设为首页 | 免责声明 | 关于勤云 | 加入收藏

Copyright©北京勤云科技发展有限公司  京ICP备09084417号