首页 | 本学科首页   官方微博 | 高级检索  
相似文献
 共查询到20条相似文献,搜索用时 31 毫秒
1.
Hong Xie  Xuan Zhu 《水文研究》2013,27(25):3685-3693
Evapotranspiration is an important component of the water and energy balance. It is dependent on climate. Precipitation, solar radiation, temperature, humidity, and wind all contribute to the rate of evapotranspiration. In this study, the temporal trends of reference evapotranspiration (ETref) and four main ETref drivers, namely, mean air temperature (Ta), wind speed (u2), net radiation (Rn) and actual vapour pressure (ea) from 1970 to 2009, were calculated based on 75 meteorological stations on the Tibetan Plateau. The results showed that the ETref on the Tibetan Plateau decreased on average by 0.6909 mm a‐1a‐1 from 1970 to 2009. Ta and ea showed an increasing trend, whereas u2 and Rn exhibited a decreasing trend. To explore the underlying causes of the ETref variation, an attribution analysis was performed to quantify the contribution of Ta, u2, Rn and ea, which showed that the changes in u2, Rn and ea produced the negative effect, whereas Ta produced the positive effect on ETref rates. The changes in u2 were found to produce the largest decrease (?0.7 mm) in ETref, followed by ea (?0.4 mm) and Rn (?0.1 mm). Although the significant increase in Ta had a large positive effect (0.51 mm) on ETref rates, changes in the other three variables each reduced ETref rates, resulting in an overall negative trend in ETref. Copyright © 2012 John Wiley & Sons, Ltd.  相似文献   

2.
Assessment of potential climate change impacts on stream water temperature (Ts) across large scales remains challenging for resource managers because energy exchange processes between the atmosphere and the stream environment are complex and uncertain, and few long‐term datasets are available to evaluate changes over time. In this study, we demonstrate how simple monthly linear regression models based on short‐term historical Ts observations and readily available interpolated air temperature (Ta) estimates can be used for rapid assessment of historical and future changes in Ts. Models were developed for 61 sites in the southeastern USA using ≥18 months of observations and were validated at sites with longer periods of record. The Ts models were then used to estimate temporal changes in Ts at each site using both historical estimates and future Ta projections. Results suggested that the linear regression models adequately explained the variability in Ts across sites, and the relationships between Ts and Ta remained consistent over 37 years. We estimated that most sites had increases in historical annual mean Ts between 1961 and 2010 (mean of +0.11 °C decade?1). All 61 sites were projected to experience increases in Ts from 2011 to 2060 under the three climate projections evaluated (mean of +0.41 °C decade?1). Several of the sites with the largest historical and future Ts changes were located in ecoregions home to temperature‐sensitive fish species. This methodology can be used by resource managers for rapid assessment of potential climate change impacts on stream water temperature. Copyright © 2014 John Wiley & Sons, Ltd.  相似文献   

3.
Flow regulation is widely known to modify the thermal regime of rivers. Here, we examine the sensitivity of an empirical approach, the Equilibrium Temperature Concept (ETC), to detect both the effects of hydraulic infrastructures on the annual thermal cycle and the recovery of the thermal equilibrium with the atmosphere. Analysis was undertaken in a Pyrenean river (the Noguera Pallaresa, Ebro basin) affected by a series of reservoirs and hydropower plants. Equilibrium temperature (Te) is defined as the water temperature (Tw) at which the sum of all heat fluxes is zero. Based on the assumption of a linear relationship between Te and Tw, we identified changes in the TeTw regression slope, used as an indicator of a thermal alteration in river flow. We also assessed the magnitude of the alteration by examining the regression slope and its statistical significance. Variations in the regression parameters were used as indicators of the influence of factors other than atmospheric conditions on water temperature. Observed Tw showed a linear relationship with Te at all river stations. However, the slopes of the TeTw relationship appeared to be lower in the reaches downstream from hydraulic infrastructures, particularly below large dams. A seasonal analysis indicated that TeTw relationships had higher slopes and lower p‐values during autumn, while no significant differences were found at other seasons. Although thermal characteristics did not strongly depend on atmospheric conditions downstream of hydraulic infrastructures, the river recovered to pre‐alteration conditions with distance downstream, indicating the natural tendency of water to attain thermal equilibrium with the atmosphere. Accepting associated uncertainties, mostly because of the quality of the data and the lack of consideration of other factors influencing the thermal regime (e.g. discharge), ETC appears to be a simple and effective method to identify thermal alterations in regulated rivers. Copyright © 2015 John Wiley & Sons, Ltd.  相似文献   

4.
Hydrologic models that rely on site specific linear and non‐linear regression water temperature (Tw) subroutines forced solely with observed air temperature (Ta) may not accurately estimate Tw in mixed‐use urbanizing watersheds where hydrogeological and land use complexity may confound common Tw regime assumptions. A nested‐scale experimental watershed study design was used to test Tw model predictions in a representative mixed‐use urbanizing watershed of the central USA. The linear regression Tw model used in the Soil and Water Assessment Tool (SWAT), a non‐linear regression Tw model, and a process‐based Tw model that accounts for watershed hydrology were evaluated. The non‐linear regression Tw model tested at a daily time step performed significantly (P < 0.01) better than the linear Tw model currently used in SWAT. Both regression Tw models overestimated Tw in lower temperature ranges (Tw < 10.0 °C) with percent bias (PBIAS) values ranging from ?28.2% (non‐linear Tw model) to ?66.1% (linear regression Tw model) and underestimated Tw in the higher temperature range (Tw > 25.0 °C) by 3.2%, and 7.2%, respectively. Conversely, the process‐based Tw model closely estimated Tw in lower temperature ranges (PBIAS = 4.5%) and only slightly underestimated Tw in the higher temperature range (PBIAS = 1.7%). Findings illustrate the benefit of integrating process‐based Tw models with hydrologic models to improve model transferability and Tw predictive confidence in urban mixed‐land use watersheds. The findings in this work are distinct geographically and in terms of mixed‐land use complexity and are therefore of immediate value to land‐use managers in similarly urbanizing watersheds globally. Copyright © 2015 John Wiley & Sons, Ltd.  相似文献   

5.
Sustainable water management in semi-arid agriculture practices requires quantitative knowledge of water fluxes within the soil-vegetation-atmosphere system. Therefore, we used stable-isotope approaches to evaluate evaporation (Ea), transpiration (Ta), and groundwater recharge (R) at sites in Senegal's Groundnut basin and Ferlo Valley pasture region during the pre-monsoon, monsoon, and post-monsoon seasons of 2021. The approaches were based upon (i) the isothermal evaporation model (for quantifying Ea); (ii) water and isotope mass balances (to partition Ea and Ta for groundnut and pasture); and (iii) the piston displacement method (for estimating R). Ea losses derived from the isothermal evaporation model corresponded primarily to Stage II evaporation, and ranged from 0.02 to 0.09 mm d−1 in the Groundnut basin, versus 0.02–0.11 mm d−1 in Ferlo. At the groundnut site, Ea rates ranged from 0.01 to 0.69 mm d−1; Ta was in the range 0.55–2.29 mm d−1; and the Ta/ETa ratio was 74%–90%. At the pasture site, the ranges were 0.02–0.39 mm d−1 for Ea; 0.9–1.69 mm d−1 for Ta; and 62–90% for Ta/ETa. The ETa value derived for the groundnut site via the isotope approach was similar to those from eddy covariance measurements, and also to the results from the previous validated HYDRUS-1D model. However, the HYDRUS-1D model gave a lower Ta/ETa ratio (23.2%). The computed groundwater recharge for the groundnut site amounted to less than 2% of the local annual precipitation. Recommendations are made regarding protocols for preventing changes to isotopic compositions of water in samples that are collected in remote arid regions, but must be analysed days later. The article ends with suggestions for studies to follow up on evidence that local aquifers are being recharged via preferential pathways.  相似文献   

6.
We have determined the elastic thicknessTe of the oceanic lithosphere along two volcanic chains of the South Central Pacific: Cook-Austral and Society islands. We used a three-dimensional spatial method to model the lithospheric flexure assuming a continuous elastic plate. The model was constrained by geoid height data from the SEASAT satellite.Along the Cook-Austral chain the elastic thickness increases westward, from 2–4 km at McDonald hot spot to 14 km at Rarotonga. At McDonald seamount, however, the data are better explained by a local compensation model. The observed trend shows an increase ofTe with age of plate at loading time. However, the elastic layer under the Cook-Austral appears systematically thinner by several kilometers than expected for “normal” seafloor, suggesting that substantial thermal thinning has taken place in this region. Considering the apparent thermal age of the plate instead of crustal age improves noticeably the results. Along the Society chainTe varies from 20 km under Tahiti to 13 km under Maupiti which is located 500 km westward. When plotting together the Society and Cook-AustralTe results versus age of load, we notice that within the first five million years after loading,Te decreases significantly while tending rapidly to an equilibrium value. This may be interpreted as the effect of initial stress relaxation which occurs just after loading inside the lower lithosphere and suggests that the presently measured elastic thickness under the very young Tahiti load ( 0.8 Ma) is not yet the equilibrium thickness.  相似文献   

7.
Changes in rainfall erosivity are an expected consequence of climate change. Long‐term series of the single storm erosion index, EI, may be analysed to detect trends in rainfall erosivity. An indirect approach has to be applied for estimating EI, given that long series of rainfall intensities are seldom available. In this paper, a method for estimating EI from the corresponding rainfall amount, he, was developed for Sicily. This method was then applied at 17 Sicilian locations, representative of different climatic zones of the region, to generate a long series (i.e. from 1916 to 1999 in most cases) of EI values. Linear and step (step located at 1970) trends in annual and seasonal erosivity were detected by both classical approaches (Mann–Kendall test, Wilcoxon‐Mann‐Whitney rank‐sum test) and a new empirical approach (quantile approach, QA), based on the determination of the erosivity values corresponding to selected probability levels. A power relationship between EI and he with a space‐ and time‐variable scale factor and a time‐variable process parameter yielded the most accurate predictions of EI. However, a simpler model, using a time‐variable scale factor and a constant process parameter, yielded reasonably accurate EI estimates. Annual erosivity did not increase in Sicily during the twentieth century. At the most, it decreased at a few locations (three of the 17 considered locations). Significant trends were observed more frequently for winter erosivity (six locations) than for summer erosivity (two locations), suggesting that the erosive storms of winter determined the occasional occurrence of a negative trend in annual erosivity. In general, the QA compared reasonably well with more classical approaches. The QA appears promising since step trends for different return periods may be detected but efforts are needed to statistically formalize the proposed approach. Copyright © 2007 John Wiley & Sons, Ltd.  相似文献   

8.
ABSTRACT

From a series of 83 daily heat balances for the lake la Godivelle (France), daily estimates of the rate of heating S (or variation in the heat content) of the lake are analysed by the methods of Edinger et al. (1968), Keijman (1974) for equilibrium temperature. This enables the calculation procedure and the formulation of the exchange coefficient K, and the equilibrium temperature Te, to be distinguished which give the best estimate of the rate of heating. In addition the data are used to describe the daily and seasonal behaviour of K and Te for the summer and autumn stratification phases of the lake.  相似文献   

9.
In spite of important differences in structural response to near‐fault and far‐fault ground motions, this paper aims at extending well‐known concepts and results, based on elastic and inelastic response spectra for far‐fault motions, to near‐fault motions. Compared are certain aspects of the response of elastic and inelastic SDF systems to the two types of motions in the context of the acceleration‐, velocity‐, and displacement‐sensitive regions of the response spectrum, leading to the following conclusions. (1) The velocity‐sensitive region for near‐fault motions is much narrower, and the acceleration‐sensitive and displacement‐sensitive regions are much wider, compared to far‐fault motions; the narrower velocity‐sensitive region is shifted to longer periods. (2) Although, for the same ductility factor, near‐fault ground motions impose a larger strength demand than far‐fault motions—both demands expressed as a fraction of their respective elastic demands—the strength reduction factors Ry for the two types of motions are similar over corresponding spectral regions. (3) Similarly, the ratio um/u0 of deformations of inelastic and elastic systems are similar for the two types of motions over corresponding spectral regions. (4) Design equations for Ry (and for um/u0) should explicitly recognize spectral regions so that the same equations apply to various classes of ground motions as long as the appropriate values of Ta, Tb and Tc are used. (5) The Veletsos–Newmark design equations with Ta=0.04 s, Tb=0.35 s, and Tc=0.79 s are equally valid for the fault‐normal component of near‐fault ground motions. Copyright © 2001 John Wiley & Sons, Ltd.  相似文献   

10.
This paper focuses on examining the effects of frequency content of the ground motion on the inelastic demands imposed on both single degree of freedom (SDF) and multi degree of freedom (MDF) steel‐framed systems. A detailed literature review is conducted to identify the indicator that best represents the frequency content of ground motion. The mean period (Tm) of ground motion is selected owing to its ability to distinguish between various spectral shapes of ground motion, and its relationship with magnitude, distance and site characteristics. Inelastic displacement demands on SDF systems for target ductility levels are first studied in the light of Tm, using a suite of 128 ground motion records. The study is then extended to MDF systems with the help of incremental dynamic analysis by employing the same ground motion ensemble to assess the influence of Tm on various engineering demand parameters. The results obtained indicate that, for SDF systems, the amplification of displacements occurs when the period ratio between elastic period (Te) and Tm is lower than unity. For MDF systems, the results demonstrate that the influence of higher modes on the base shear and maximum storey drift profile becomes more pronounced, as Tm approaches the higher mode periods of the structure. These observations, for both SDF and MDF systems, tend to be more evident for higher levels of inelasticity. The significance of the results, with particular reference to European seismic design procedures, is highlighted. Copyright © 2010 John Wiley & Sons, Ltd.  相似文献   

11.
A procedure based on rigorous non‐linear analysis is presented that estimates the peak deformation among all isolators in an asymmetric building due to strong ground motion. The governing equations are reduced to a form such that the median normalized deformation due to an ensemble of ground motions with given corner period Td depends primarily on four global parameters of the isolation system: the isolation period Tb, the normalized strength η, the torsional‐to‐lateral frequency ratio Ωθ, and the normalized stiffness eccentricity eb/r. The median ratio of the deformations of the asymmetric and corresponding symmetric systems is shown to depend only weakly on Tb, η, and Ωθ, but increases with eb/r. The equation developed to estimate the largest ratio among all isolators depends only on the stiffness eccentricity and the distance from the center of mass to the outlying isolator. This equation, multiplied by an earlier equation for the deformation of the corresponding symmetric system, provides a design equation to estimate the deformations of asymmetric systems. This design equation conservatively estimates the peak deformation among all isolators, but is generally within 10% of the ‘exact’ value. Relative to the non‐linear procedure presented, the peak isolator deformation is shown to be significantly underestimated by the U.S. building code procedures. Copyright © 2003 John Wiley & Sons, Ltd.  相似文献   

12.
Three practical schemes for computing the snow surface temperature Ts, i.e. the force–restore method (FRM), the surface conductance method (SCM), and the Kondo and Yamazaki method (KYM), were assessed with respect to Ts retrieved from cloud‐free, NOAA‐AVHRR satellite data for three land‐cover types of the Paddle River basin of central Alberta. In terms of R2, the mean Ts, the t‐test and F‐test, the FRM generally simulated more accurate Ts than the SCM and KYM. The bias in simulated Ts is usually within several degrees Celsius of the NOAA‐AVHRR Ts for both the calibration and validation periods, but larger errors are encountered occasionally, especially when Ts is substantially above 0 °C. Results show that the simulated Ts of the FRM is more consistent than that of the SCM, which in turn was more consistent than that of the KYM. This is partly because the FRM considers two aspects of heat conduction into snow, a stationary‐mean diurnal (sinusoidal) temperature variation at the surface coupled to a near steady‐state ground heat flux, whereas the SCM assumes a near steady‐state, simple heat conduction, and other simplifying assumptions, and the KYM does not balance the snowpack heat fluxes by assuming the snowpack having a vertical temperature profile that is linear. Copyright © 2005 John Wiley & Sons, Ltd.  相似文献   

13.
In this paper, downscaling models are developed using various linear regression approaches, namely direct, forward, backward and stepwise regression, for obtaining projections of mean monthly maximum and minimum temperatures (Tmax and Tmin) to lake‐basin scale in an arid region in India. The effectiveness of these regression approaches is evaluated through application to downscale the predictands for the Pichola lake region in the state of Rajasthan in India, which is considered to be a climatically sensitive region. The predictor variables are extracted from (i) the National Centers for Environmental Prediction (NCEP) reanalysis dataset for the period 1948–2000 and (ii) the simulations from the third‐generation Canadian Coupled Global Climate Model (CGCM3) for emission scenarios A1B, A2, B1 and COMMIT for the period 2001–2100. The selection of important predictor variables becomes a crucial issue for developing downscaling models as reanalysis data are based on a wide range of meteorological measurements and observations. A simple multiplicative shift was used for correcting predictand values. Direct regression was found to yield better performance among all other regression techniques for the training data set, while the forward regression technique performed better in the validation data set, explored in the present study. For trend analysis, the Mann–Kendall non‐parametric test was performed. The results of downscaling models show that an increasing trend is observed for Tmax and Tmin for A1B, A2 and B1 scenarios, whereas no trend is discerned with the COMMIT scenario by using predictors. Copyright © 2010 John Wiley & Sons, Ltd.  相似文献   

14.
Gerard Govers  Jan Diels 《水文研究》2013,27(25):3777-3790
Experimental work has clearly shown that the effective hydraulic conductivity (Ke) or effective infiltration rate (fe) on the local scale of a plot cannot be considered as constant but are dependent on water depth and rainfall intensity because non‐random microtopography‐related variations in hydraulic conductivity occur. Rainfall–runoff models generally do not account for this: models assume that excess water is uniformly spread over the soil surface and within‐plot variations are neglected. In the present study, we propose a model that is based on the concepts of microtopography‐related water depth‐dependent infiltration and partial contributing area. Expressions for the plot scale Ke and fe were developed that depend on rainfall intensity and runon from upslope (and thus on water depth). To calibrate and validate the model, steady state infiltration experiments were conducted on maize fields on silt loam soils in Belgium, with different stages and combinations of rainfall intensity and inflow, simulating rainfall and runon. Water depth–discharge and depth–inundation relationships were established and used to estimate the effect of inundation on Ke. Although inflow‐only experiments were found to be unsuitable for calibration, the model was successfully calibrated and validated with the rainfall simulation data and combined rainfall–runon data (R²: 0.43–0.91). Calibrated and validated with steady state infiltration experiments, the model was combined with the Green–Ampt infiltration equation and can be applied within a two‐dimensional distributed rainfall–runoff model. The effect of water depth–dependency and rainfall intensity on infiltration was illustrated for a hillslope. Copyright © 2012 John Wiley & Sons, Ltd.  相似文献   

15.
From stable carbon isotope analysis of tree-rings of Chinese pine (Pinus tabulaeformis) from Mt. Helan, China, we found that high-δ13C values were related to high mean temperatures from June to August (T 68), and Iow-δ13C values corresponded to low T 68. From these data, a transfer function has been used to reconstruct summer temperatures (T 68) for the Mt. Helan region. The explained variance of reconstruction is 34.9% (F=15.01, p<0.001). The time period containing the highest summer temperatures in northern China (late 1920-1930s) was confirmed by our reconstruction. The data indicate that there is a tele-connection between summer temperatures in Mt. Helan area and sea-surface-temperatures in the tropical Pacific. The extreme low temperature periods around the years of 1920 and 1947 for Mt. Helan region correspond well to the cold climate in the tropical Pacific. Along with other analyses, this suggests that climate variations in the Mt. Helan region are driven not only by local events, but also by the global climate. Significant periodicities appearing in the reconstruction are 2.56 and 2.63 years.  相似文献   

16.
A database of the electron temperature (Te) comprising of most of the available LEO satellite measurements is used for studying the solar activity variations of Te. The Te data are grouped for two levels of solar activity (low LSA and high HSA), five altitude ranges between 350 and 2000 km, and day and night. By fitting a theoretical expression to the Te values we obtain variation of Te along magnetic field lines and heat flux for LSA and HSA. We have found that Te increases with increase in solar activity at low and mid-latitudes during nighttime at all altitudes studied. During daytime the Te response to solar activity depends on latitude, altitude, and season. This analysis shows existence of anti-correlation between Te and solar activity at mid-latitudes below 700 km during the equinox and winter day hours. Heat fluxes show small latitudinal dependence for daytime but substantial for nighttime.  相似文献   

17.
Abstract

It is shown that magnetic fields generated by flows v r,(r,t)er+vT where vT is an arbitrary toroidal component (er˙vT≡V≡vT≡0), cannot be maintained indefinitely against ohmic dissipation. The poloidal field variable max |r 2 B r| is shown to decay strictly monotonically with an undetermined decay rate. A bound on the growth of the toroidal field norm ∥T1 is established solely dependent on the rate of conversion of poloidal to toroidal field, so that when the poloidal field is negligible then ∥T1 decays strictly monotonically. The main application of these results is to models of stellar evolution based on axisymmetric differential rotation and spherically symmetric contraction. This symmetric velocity theorem overlaps with two already known theorems, namely the toroidal velocity theorem where v r≡0 and the radial velocity theorem where vT≡0. The new theorem does not entirely include the already established ones, principal differences being in the rates of decay and the field variables for which the decay is proven (see Table 1).  相似文献   

18.
Examples of data from DE-2 satellite instruments are presented. These illustrate the behaviour of plasma parameters in the F-region and adjacent topside ionosphere during rapid sub-auroral ion drift (SAID) events. In particular, a variety of behaviours of the electron temperature (Te) is demonstrated, both within and equatorward of the SAID region. The Sheffield University plasmasphere-ionosphere model (SUPIM) is used to perform calculations in which a model SAID is applied to a plasma flux tube. The model results indicate that strongly elevated ion temperature (a recognised signature of SAID events) is on occasion sufficient to raise Te to observed values by ion-electron heat transfer. On other occasions, an additional heat source is required. It is suggested that such a source for the electron gas may be due to interaction between the ring current and the plasmasphere at high altitudes. The magnitude of the downward heat flux is consistent with that necessary to produce sub-auroral red arcs. The resulting strongly heated electron gas causes vibrational excitation of molecular nitrogen in the thermosphere.  相似文献   

19.
This paper deals with the estimation of peak inelastic displacements of SDOF systems, representative of typical steel structures, under constant relative strength scenarios. Mean inelastic deformation demands on bilinear systems (simulating moment resisting frames) are considered as the basis for comparative purposes. Additional SDOF models representing partially‐restrained and concentrically‐braced (CB) frames are introduced and employed to assess the influence of different force‐displacement relationships on peak inelastic displacement ratios. The studies presented in this paper illustrate that the ratio between the overall yield strength and the strength during pinching intervals is the main factor governing the inelastic deformations of partially‐restrained models and leading to significant differences when compared with predictions based on bilinear structures, especially in the short‐period range. It is also shown that the response of CB systems can differ significantly from other pinching models when subjected to low or moderate levels of seismic demand, highlighting the necessity of employing dedicated models for studying the response of CB structures. Particular attention is also given to the influence of a number of scalar parameters that characterise the frequency content of the ground motion on the estimated peak displacement ratios. The relative merits of using the average spectral period Taver, mean period Tm, predominant period Tg, characteristic period Tc and smoothed spectral predominant period To of the earthquake ground motion, are assessed. This paper demonstrates that the predominant period, defined as the period at which the input energy is maximum throughout the period range, is the most suitable frequency content scalar parameter for reducing the variability in displacement estimations. Finally, noniterative equivalent linearisation expressions based on the secant period and equivalent damping ratios are presented and verified for the prediction of peak deformation demands in steel structures. Copyright © 2011 John Wiley & Sons, Ltd.  相似文献   

20.
Adsorption of reactive black 5 (RB5) from aqueous solution onto chitosan was investigated in a batch system. The effects of solution pH, initial dye concentration, and temperature were studied. Adsorption data obtained from different batch experiments were modeled using both pseudo first‐ and second‐order kinetic equations. The equilibrium adsorption data were fitted to the Freundlich, Tempkin, and Langmuir isotherms over a dye concentration range of 45–100 µmol/L. The best results were achieved with the pseudo second‐order kinetic and Langmuir isotherm equilibrium models, respectively. The equilibrium adsorption capacity (qe) was increased with increasing the initial dye concentration and solution temperature, and decreasing solution pH. The chitosan flakes for the adsorption of the dye was regenerated efficiently through the alkaline solution and was then reused for dye removal. The activation energy (Ea) of sorption kinetics was estimated to be 13.88 kJ/mol. Thermodynamic parameters such as changes in free energy (ΔG), enthalpy (ΔH), and entropy (ΔS) were evaluated by applying the van't Hoff equation. The thermodynamics of reactive dye adsorption by chitosan indicates its spontaneous and endothermic nature.  相似文献   

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

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