首页 | 本学科首页   官方微博 | 高级检索  
相似文献
 共查询到20条相似文献,搜索用时 26 毫秒
1.
Momentum and turbulent kinetic energy (TKE) budgets across a forest edge have been investigated using large-eddy simulation (LES). Edge effects are observed in the rapid variation of a number of budget terms across this vegetation transition. The enhanced drag force at the forest edge is largely balanced by the pressure gradient force and by streamwise advection of upstream momentum, while vertical turbulent diffusion is relatively insignificant. For variance and TKE budgets, the most important processes at the forest edge are production due to the convergence (or divergence) of the mean flow, streamwise advection, pressure diffusion and enhanced dissipation by canopy drag. Turbulent diffusion, pressure redistribution and vertical shear production, which are characteristic processes in homogeneous canopy flow, are less important at the forest transition. We demonstrate that, in the equilibrated canopy flow, a substantial amount of TKE produced in the streamwise direction by the vertical shear of the mean flow is redistributed in the vertical direction by pressure fluctuations. This redistribution process occurs in the upper canopy layers. Part of the TKE in the vertical velocity component is transferred by turbulent and pressure diffusion to the lower canopy levels, where pressure redistribution takes place again and feeds TKE back to the streamwise direction. In this TKE cycle, the primary source terms are vertical shear production for streamwise velocity variance and pressure redistribution for vertical velocity variance. The evolution of these primary source terms downwind of the forest edge largely controls the adjustment rates of velocity variances.  相似文献   

2.
Mean Flow Near Edges and Within Cavities Situated Inside Dense Canopies   总被引:3,自引:3,他引:0  
A streamfunction-vorticity formulation is used to explore the extent to which turbulent and turbulently inviscid solutions to the mean momentum balance explain the mean flow across forest edges and within cavities situated inside dense forested canopies. The turbulent solution is based on the mean momentum balance where first-order closure principles are used to model turbulent stresses. The turbulently inviscid solution retains all the key terms in the mean momentum balance but for the turbulent stress gradients. Both exit and entry versions of the forest edge problem are explored. The turbulent solution is found to describe sufficiently the bulk spatial patterns of the mean flow near the edge including signatures of different length scales reported in canopy transition studies. Next, the ‘clearing inside canopy’ or the so-called ‘cavity’ problem is solved for the inviscid and turbulent solutions and then compared against flume experiments. The inviscid solution describes the bulk flow dynamics in much of the zones within the cavity. In particular, the solution can capture the correct position of the bulk recirculation zone within the cavity, although with a weaker magnitude. The inviscid solution cannot capture the large vertical heterogeneity in the mean velocity above the canopy, as expected. These features are better captured via the first-order closure representation of the turbulent solution. Given the ability of this vorticity formulation to capture the mean pressure variations and the mean advective acceleration terms, it is ideal for exploring the distributions of scalars and roughness-induced flow adjustments on complex topography.  相似文献   

3.
Most of our knowledge on forest-edge flows comes from numerical and wind-tunnel experiments where canopies are horizontally homogeneous. To investigate the impact of tree-scale heterogeneities (\({>}1\) m) on the edge-flow dynamics, the flow in an inhomogeneous forest edge on Falster island in Denmark is investigated using large-eddy simulation. The three-dimensional forest structure is prescribed in the model using high resolution helicopter-based lidar scans. After evaluating the simulation against wind measurements upwind and downwind of the forest leading edge, the flow dynamics are compared between the scanned forest and an equivalent homogeneous forest. The simulations reveal that forest inhomogeneities facilitate flow penetration into the canopy from the edge, inducing important dispersive fluxes in the edge region as a consequence of the flow spatial variability. Further downstream from the edge, the forest inhomogeneities accentuate the canopy-top turbulence and the skewness of the wind-velocity components while the momentum flux remains unchanged. This leads to a lower efficiency in the turbulent transport of momentum within the canopy. Dispersive fluxes are only significant in the upper canopy. Above the canopy, the mean flow is less affected by the forest inhomogeneities. The inhomogeneities induce an increase in the mean wind speed that was found to be equivalent to a decrease in the aerodynamic height of the canopy. Overall, these results highlight the importance of forest inhomogeneities when looking at canopy–atmosphere exchanges in forest-edge regions.  相似文献   

4.
The structure of turbulent flows along a transition between tall-forested canopies and forest clearings continues to be an active research topic in canopy turbulence. The difficulties in describing the turbulent flow along these transitions stem from the fact that the vertical structure of the canopy and its leaf area distribution cannot be ignored or represented by an effective roughness length. Large-eddy simulation (LES) runs were performed to explore the effect of a homogeneous variation in the forest leaf area index (LAI) on the turbulent flow across forest edges. A nested grid numerical method was used to ensure the development of a deep boundary layer above the forest while maintaining a sufficiently high resolution in the region close to the ground. It was demonstrated that the LES here predicted first-order and second-order mean velocity statistics within the canopy that agree with reported Reynolds-Averaged Navier–Stokes (RANS) model results, field and laboratory experiments. In the simulations reported here, the LAI was varied between 2 and 8 spanning a broad range of observed LAI in terrestrial ecosystems. By increasing the forest LAI, the mean flow properties both within the forest and in the clearing near the forest edge were altered in two fundamental ways: near the forest edge and into the clearing, the flow statistical properties resembled the so-called back-facing step (BFS) flow with a mean recirculation zone near the edge. Another recirculation zone sets up downstream of the clearing as the flow enters the tall forest canopy. The genesis of this within-forest recirculation zone can be primarily described using the interplay between the mean pressure gradients (forcing the flow) and the drag force (opposing the flow). Using the LES results, a simplified analytical model was also proposed to explain the location of the recirculation zone inside the canopy and its dependence on the forest LAI. Furthermore, a simplified scaling argument that decomposes the mean velocity at the outflow edge into a superposition of ‘exit flow’ and BFS-like flow with their relative importance determined by LAI was explored.  相似文献   

5.
The turbulent flow in and above plant canopies is of fundamental importance to the understanding oftransport processes of momentum,heat and mass between plant canopies and atmosphere,and to microme-teorology.The Reynolds stress equation model(RSM)has been applied to calculate the turbulence in cano-pies in this paper.The calculated mean wind velocity profiles,Reynolds stress,turbulent kinetic energy andviscous dissipation rate in a corn canopy and a spruce forest are compared with field observed data and withWilson's and Shaw's model.The velocity profiles and Rynolds stress calculated by both models are in goodagreement,and the length scale of turbulence appears to be similar.  相似文献   

6.
植被内部及其上方湍流场的数值模拟   总被引:3,自引:0,他引:3       下载免费PDF全文
尹协远  J.D.Lin 《气象学报》1988,46(2):194-201
植被内部及其上方的湍流流场对于了解植被与大气之间的动量、热量和质量交换过程极其重要。本文把高阶湍流封闭模型的Reynolds应力方程模型(RSM)应用于植被湍流的计算,得到了风速、湍流动能、Reynolds应力及能量耗散率的垂直分布,与现场观测数据比较,甚为满意。  相似文献   

7.
We use a nonhydrostatic shelterbelt boundary-layer turbulence model with Mellor–Yamada second-order closure to evaluate quantitatively the dynamic processes of surface boundary-layer flow perturbed by shelterbelts of different densities and to understand the shelter mechanism. We first analyze the drag exerted on air by shelterbelts of different densities, a root cause of any shelter function, and the resulting wind reduction. The results show that the effectiveness of a shelter is determined not only by its total drag but also by the distribution of the drag-generated momentum deficit in the sheltered area, and that medium-dense shelterbelts have the maximum shelter effect. We also analyze the horizontal momentum budget and find that the shelter mechanism is the product of several processes. The results reveal that strong vertical mean transport and the pressure gradient also play important roles in shelter efficiency. The pressure perturbation caused by the shelter extends far downstream of the shelter, and combines with advective transport to provide the larger shelter efficiency of medium-dense shelterbelts. We finally analyze the changes of perturbed pressure, turbulence, and vertical velocity with shelterbelt density to further clarify the shelter mechanism.  相似文献   

8.
In this paper, we lay the foundations of a systematic mathematical formulation for the governing equations for flow through an urban canopy (e.g., coarse-scaled building array) where the effects of the unresolved obstacles on the flow are represented through a distributed mean-momentum sink. This, in turn, implies additional corresponding terms in the transport equations for the turbulence quantities. More specifically, a modified k-- model is derived for the simulation of the mean wind speed and turbulence for a neutrally stratified flow through and over a building array, where groups of buildings in the array are aggregated and treated as a porous medium. This model is based on time averaging the spatially averaged Navier--Stokes equations, in which the effects of the obstacle--atmosphere interaction are included through the introduction of a volumetric momentum sink (representing drag on the unresolved buildings in the array).The k-- turbulence closure model requires two additional prognostic equations, namely one for the time-averaged resolved-scale kinetic energy of turbulence,, and another for the dissipation rate, , of . The transport equation for is derived directly from the transport equation for the spatially averaged velocity, and explicitly includes additional sources and sinks that arise from time averaging the product of the spatially averaged velocity fluctuations and the distributed drag force fluctuations. We show how these additional source/sink terms in the transport equation for can be obtained in a self-consistent manner from a parameterization of the sink term in the spatially averaged momentum equation. Towards this objective, the time-averaged product of the spatially averaged velocity fluctuations and the distributed drag force fluctuations can be approximated systematically using a Taylor series expansion. A high-order approximation is derived to represent this source/sink term in the transport equation for . The dissipation rate () equation is simply obtained as a dimensionally consistent analogue of the equation. The relationship between the proposed mathematical formulation of the equations for turbulent flow within an urban canopy (where the latter is treated as a porous medium) and an earlier heuristic two-band spectral decomposition for parameterizing turbulence in a plant canopy is explored in detail.  相似文献   

9.
An analyticalone-dimensional model of momentum transferby vegetation with variable foliage distribution,sheltering and drag coefficientis developed. The model relies on a simpleparameterization of the ratio of theabove-canopy friction velocity, u*, to thewind speed at the top of the canopy,u(h), to predict vegetation roughness length(z0) and displacement height(d) as functions of canopy height (h) and dragarea index. Model predictionsof d/h and z/h compare very favorably withobserved values.A model sensitivity analysis suggests that shelteringeffects for momentum transfertend to make canopies with non-uniform foliagedistribution resemble canopies withmore uniform foliage distribution and that anyinfluence wind speed has on d/hand z0/h is more likely to be related to theinfluence that wind speed may haveon u*/u(h) rather than the influence windspeed may have on the foliage dragcoefficient. Model results indicate that z0/hand d/h are sensitive to uncertaintiesin the numerical values of the model parameters,foliage density and distribution,sheltering effects and variations in drag coefficientwithin the canopy. In additionz0/h is also shown to be sensitive to thepresence or absence of the roughnesssublayer. Given the simplicity of the model it issuggested that it may be of usefor land surface parameterizations in large scalemodels.  相似文献   

10.
Edge Flow and Canopy Structure: A Large-Eddy Simulation Study   总被引:4,自引:4,他引:0  
Sharp heterogeneities in forest structure, such as edges, are often responsible for wind damage. In order to better understand the behaviour of turbulent flow through canopy edges, large-eddy simulations (LES) have been performed at very fine scale (2 m) within and above heterogeneous vegetation canopies. A modified version of the Advanced Regional Prediction System (ARPS), previously validated in homogeneous conditions against field and wind-tunnel measurements, has been used for this purpose. Here it is validated in a simple forest-clearing-forest configuration. The model is shown to be able to reproduce accurately the main features observed in turbulent edge flow, especially the “enhanced gust zone” (EGZ) present around the canopy top at a few canopy heights downwind from the edge, and the turbulent region that develops further downstream. The EGZ is characterized by a peak in streamwise velocity skewness, which reflects the presence of intense intermittent wind gusts. A sensitivity study of the edge flow to the forest morphology shows that with increasing canopy density the flow adjusts faster and turbulent features such as the EGZ become more marked. When the canopy is characterized by a sparse trunk space the length of the adjustment region increases significantly due to the formation of a sub-canopy wind jet from the leading edge. It is shown that the position and magnitude of the EGZ are related to the mean upward motion formed around canopy top behind the leading edge, caused by the deceleration in the sub-canopy. Indeed, this mean upward motion advects low turbulence levels from the bottom of the canopy; this emphasises the passage of sudden strong wind gusts from the clearing, thereby increasing the skewness in streamwise velocity as compared with locations further downstream where ambient turbulence is stronger.  相似文献   

11.
A sequence of numerical calculations has been made for the equilibrium balances of eddies and mean currents in open and partially blocked, periodic channels. The physical model employed is a two-layer, quasigeostrophic, wind-driven one, with important bottom friction and weak lateral friction. The resolved eddies provide the interior fluxes of momentum and potential vorticity which allow the mean state to be a balanced one. The set of calculations does not provide a parameter study as such, but does provide examples of the influences of alternative physical processes and geometrical constraints. These alternatives include the presence or absence of a partial barrier across the channel, the length of the channel, the addition of a transient component to the wind-driving, and the addition of a topographic sill across the channel gap. Particular attention is focused upon the steadily driven general circulation of a β-plane channel, because of the structural simplicity of the solution. The results may be broadly summarized as follows. The eddies are generated by a baroclinic instability of the mean flow. They act to intensity the upper layer mean jet and mean cross-jet potential vorticity gradient (through eddy horizontal Reynolds stress and relative vorticity flux divergence, respectively) and to transfer downwards mean zonal momentum, energy, and potential vorticity gradient (through eddy interfacial pressure drag, vertical pressure work, and vortex stretching flux divergence, respectively). In the case of a zonally uniform channel, the meridional heat flux is found not to conform closely to previously proposed parameterizations. The presence of a partial meridional barrier and a topographic obstacle are found to strongly influence the equilibrium solution, while neither a change in the basin length nor the presence of a transient wind component appear to importantly alter the solution.  相似文献   

12.
The Regional Atmospheric Modelling System (RAMS)-based Forest Large-Eddy Simulation (RAFLES) model is used to simulate the effects of large rectangular prism-shaped semi-porous barriers of varying densities under neutrally buoyant conditions. RAFLES model resolves flows inside and above forested canopies and other semi-porous barriers, and it accounts for barrier-induced drag on the flow and surface flux exchange between the barrier and the air. Unlike most other models, RAFLES model also accounts for the barrier-induced volume and aperture restriction via a modified version of the cut-cell coordinate system. We explicitly tested the effects of the numerical representation of volume restriction, independent of the effects of the drag, by comparing drag-only simulations (where we prescribed neither volume nor aperture restrictions to the flow), restriction-only simulations (where we prescribed no drag), and control simulations where both drag and volume plus aperture restrictions were included. Previous modelling and empirical work have revealed the development of important areas of increased uplift upwind of forward-facing steps, and recirculation zones downwind of backward-facing steps. Our simulations show that representation of the effects of the volume and aperture restriction due to the presence of semi-porous barriers leads to differences in the strengths and locations of increased-updraft and recirculation zones, and the length and strength of impact and adjustment zones when compared to simulation solutions with a drag-only representation. These are mostly driven by differences to the momentum budget of the streamwise wind velocity by resolved turbulence and pressure gradient fields around the front and back edges of the barrier. We propose that volume plus aperture restriction is an important component of the flow system in semi-porous environments such as forests and cities and should be considered by large-eddy simulation (LES).  相似文献   

13.
Dispersive Stresses at the Canopy Upstream Edge   总被引:1,自引:0,他引:1  
The derivation of flow and mass transfer models in canopy and porous media environments involves the spatial-averaging of the flow properties and their subscale equations. The averaging of the momentum equation generates the dispersive stress terms that represent the subscale spatial variations of the unresolved velocity field. While previous studies ignored the dispersive stresses in their flow models, recent evidence indicates that the dispersive stresses may be important. Here we focus our attention on the magnitude of the normal dispersive stresses in the entry region of a ‘forest patch’, where the in-canopy velocities are large and the longitudinal derivatives do not cancel out. Highly detailed particle image velocimetry measurements, at a temporal and spatial resolution of 5 Hz and 1.4 mm, are obtained inside and around a 1-m long model canopy which consists of transparent vertical cylinders 6 mm in diameter and 74.3 mm high (h). The cylinders are randomly distributed to form a relatively sparse forest patch with a leaf area density of 7.56 m−1 and a fluid volume fraction (porosity) of 0.965. We present results of the double averaged flow properties at three different regions of the forest patch; the upstream edge (x ≈ 0), the fully-developed interior region (x ≈ 10h) and the downstream edge (x ≈ 13h). We find that the normal dispersive stresses around the entry region of the forest patch are significantly larger than the normal Reynolds stresses. An order of magnitude analysis of the relevant terms in the momentum equation indicates that the longitudinal derivatives of the dispersive stresses are of the same order of magnitude as that of the drag force and similar to that of the horizontal convection term. The longitudinal derivatives of the Reynolds stresses are smaller, though cannot be ignored. Comparing these results with the characteristic profiles measured in the fully-developed region indicates that the dispersive stresses, which are generated at the forest patch entrance, decrease along an adjustment region while maintaining their profile shape. We find that the dispersive stresses influence the rate at which momentum penetrates into the canopy. These observations suggest that under certain flow conditions, dispersive stresses may dominate the momentum balance and therefore must be considered in future canopy and porous media flow models.  相似文献   

14.
Turbulence above and within canopies has characteristics distinct from that over rough surfaces. The vertical transport of momentum and scalars is dominated by coherent structures whose origin is now thought to be the result of the unstable inflexion in the profile of the mean wind speed established by the application of canopy drag. This distinctive property leads to the failure of the standard Monin–Obukhov flux–profile relationships over homogeneous canopies, relationships that are assumed in many surface exchange schemes within numerical weather prediction and general circulation models. A modification of the flux–profile relationships is presented that incorporates the effects of the canopy turbulence. The subsequent impacts on the evolution of the surface energy balance and boundary-layer state are investigated within a simple numerical model for the evolution of the boundary layer and canopy state. By comparing cases with and without the modification it is shown that canopy-generated turbulence can lead, not only to the alteration of the flux–profile relationships above the canopy, but also to a different evolution of the surface energy balance and differences in near-surface conditions that would be significant in numerical weather prediction. More fundamentally, the modifications to the flux–profile relationships imply that parameters such as the roughness length and displacement height for canopies should not be considered as invariant properties, but rather as properties that depend on the flow and hence vary systematically with the diabatic stability of the boundary layer.  相似文献   

15.
A Large Eddy Simulation (LES) model representing the air flow within and above a plant canopy layer has been completed. Using this model, the organized structures of turbulent flow in the early developmental stages of a crop are simulated and discussed in detail.The effect of the drag due to vegetation is expressed by a term added to the three-dimensional Navier-Stokes equation averaged over the grid scale. For the formulation of sub-grid turbulence processes, the equations for the time-dependent SGS (Sub-Grid-Scale) turbulence energy equation is used, which includes the effects of dissipation (both by viscosity and leaf drag), shear production and diffusion.The organized structure of turbulent flow at the air-plant interface, obtained numerically by the model, yields its contribution to momentum transfer. The three-dimensional large eddy structures, which are composed of spanwise vortices (rolls) and streamwise vortices (ribs), are simulated near the air-plant interface. They are induced by the shear instability at inflection points of the velocity profile. The structure clearly has a life cycle. The instantaneous image of the structure is similar to those observed in the field observations of Gaoet al. (1989) and in the laboratory flume experiments of Ikeda and Ota (1992). These organized structures also account for the well known fact that the sweep motion of turbulence dominates momentum transport within and just above a plant canopy, and the motion of ejection prevails in the higher regions.  相似文献   

16.
Forested landscapes often exhibit large spatial variability in vertical and horizontal foliage distributions. This variability may affect canopy-atmosphere exchanges through its action on the development of turbulent structures. Here we investigate in neutral stratification the turbulent structures encountered in a maritime pine forest characterized by a high, dense foliated layer associated with a deep and sparse trunk space. Both stand and edge regions are considered. In situ measurements and the results of large-eddy simulations are used and analyzed together. In stand conditions, far from the edge, canopy-top structures appear strongly damped by the dense crown layer. Turbulent wind fluctuations within the trunk space, where the momentum flux vanishes, are closely related to these canopy-top structures through pressure diffusion. Consequently, autocorrelation and spectral analyses are not quite appropriate to characterize the vertical scale of coherent structures in this type of canopy, as pressure diffusion enhances the actual scale of structures. At frequencies higher than those associated with canopy-top structures, wind fluctuations related to wake structures developing behind tree stems are observed within the trunk space. They manifest themselves in wind velocity spectra as secondary peaks in the inertial subrange region, confirming the hypothesis of spectral short-cuts in vegetation canopies. In the edge region specific turbulent structures develop just below the crown layer, in addition to canopy-top structures. They are generated by the wind shear induced by the sub-canopy wind jet that forms at the edge. These structures provide a momentum exchange mechanism similar to that observed at the canopy top but in the opposite direction and with a lower magnitude. They may develop as in plane mixing-layer flows, with some perturbations induced by canopy-top structures. Wake structures are also observed within the trunk space in the edge region.  相似文献   

17.
Averaging procedures for flow within vegetation canopies   总被引:13,自引:5,他引:13  
Most one-dimensional models of flow within vegetation canopies are based on horizontally averaged flow variables. This paper formalizes the horizontal averaging operation. Two averaging schemes are considered: pure horizontal averaging at a single instant, and time averaging followed by horizontal averaging. These schemes produce different forms for the mean and turbulent kinetic energy balances, and especially for the wake production term describing the transfer of energy from large-scale motion to wake turbulence by form drag. The differences are primarily due to the appearance, in the covariances produced by the second scheme, of dispersive components arising from the spatial correlation of time-averaged flow variables. The two schemes are shown to coincide if these dispersive fluxes vanish.  相似文献   

18.
A model of \(\hbox {CO}_{2}\) atmospheric transport in vegetated canopies is tested against measurements of the flow, as well as \(\hbox {CO}_{2}\) concentrations at the Norunda research station located inside a mixed pine–spruce forest. We present the results of simulations of wind-speed profiles and \(\hbox {CO}_{2}\) concentrations inside and above the forest canopy with a one-dimensional model of profiles of the turbulent diffusion coefficient above the canopy accounting for the influence of the roughness sub-layer on turbulent mixing according to Harman and Finnigan (Boundary-Layer Meteorol 129:323–351, 2008; hereafter HF08). Different modelling approaches are used to define the turbulent exchange coefficients for momentum and concentration inside the canopy: (1) the modified HF08 theory—numerical solution of the momentum and concentration equations with a non-constant distribution of leaf area per unit volume; (2) empirical parametrization of the turbulent diffusion coefficient using empirical data concerning the vertical profiles of the Lagrangian time scale and root-mean-square deviation of the vertical velocity component. For neutral, daytime conditions, the second-order turbulence model is also used. The flexibility of the empirical model enables the best fit of the simulated \(\hbox {CO}_{2}\) concentrations inside the canopy to the observations, with the results of simulations for daytime conditions inside the canopy layer only successful provided the respiration fluxes are properly considered. The application of the developed model for radiocarbon atmospheric transport released in the form of \(^{14}\hbox {CO}_{2}\) is presented and discussed.  相似文献   

19.
Using analyses of data from extant direct numerical simulations and large-eddy simulations of boundary-layer and channel flows over and within urban-type canopies, sectional drag forces, Reynolds and dispersive shear stresses are examined for a range of roughness densities. Using the spatially-averaged mean velocity profiles these quantities allow deduction of the canopy mixing length and sectional drag coefficient. It is shown that the common assumptions about the behaviour of these quantities, needed to produce an analytical model for the canopy velocity profile, are usually invalid, in contrast to what is found in typical vegetative (e.g. forest) canopies. The consequence is that an exponential shape of the spatially-averaged mean velocity profile within the canopy cannot normally be expected, as indeed the data demonstrate. Nonetheless, recent canopy models that allow prediction of the roughness length appropriate for the inertial layer’s logarithmic profile above the canopy do not seem to depend crucially on their (invalid) assumption of an exponential profile within the canopy.  相似文献   

20.
The barotropic processes associated with the development of a precipitation system are investigated through analysis of cloud-resolving model simulations of Mei-yu torrential rainfall events over eastern China in mid-June 2011. During the model integration period, there were three major heavy rainfall events: 9–12, 13–16 and 16–20 June. The kinetic energy is converted from perturbation to mean circulations in the first and second period, whereas it is converted from mean to perturbation circulations in the third period. Further analysis shows that kinetic energy conversion is determined by vertical transport of zonal momentum. Thus, the prognostic equation of vertical transport of zonal momentum is derived, in which its tendency is associated with dynamic, pressure gradient and buoyancy processes. The kinetic energy conversion from perturbation to mean circulations in the first period is mainly associated with the dynamic processes. The kinetic energy conversion from mean to perturbation circulations in the third period is generally related to the pressure gradient processes.  相似文献   

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

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