999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

How Are El Ni?o and La Ni?a Events Improved in an Eddy-Resolving Ocean General Circulation Model?

2015-11-24 09:49:50HUALiJuanandYUYongQiang
關鍵詞:血漿

HUA Li-Juanand YU Yong-Qiang

1State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics, Institute of Atmospheric Physics, Chinese Academy of Sciences, Beijing 100029, China

2College of Earth Science, University of Chinese Academy of Sciences, Beijing 100049, China

How Are El Ni?o and La Ni?a Events Improved in an Eddy-Resolving Ocean General Circulation Model?

HUA Li-Juan1,2and YU Yong-Qiang1*

1State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics, Institute of Atmospheric Physics, Chinese Academy of Sciences, Beijing 100029, China

2College of Earth Science, University of Chinese Academy of Sciences, Beijing 100049, China

The present study compares the performance of two versions of the LASG/IAP (State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics/Institute of Atmospheric Physics) Climate System Ocean Model (LICOM) in reproducing the interannual variability associated with El Ni?o and La Ni?a events in the tropical Pacific. Both versions are forced with the identical boundary conditions from observed or reanalysis data, in which one version has a finer spatial resolution of (1/10)° in the horizontal domain and 55 vertical layers, and the other version has a coarse resolution of 1° in the horizontal domain and 30 vertical layers. ENSO simulations form the two versions are compared with observations and, in particular, the improvements with regard to ENSO by the finer resolution ocean model are emphasized. As a result of the finer spatial resolution, both the vertical temperature gradient and vertical velocity are better represented in the equatorial Pacific than they are by the coarse resolution model; and thus, the corresponding vertical advections of temperature are more reasonable. Besides the mean climatology, simulated ENSO events and relevant feedbacks are much improved in the finer resolution model. A heat budget analysis suggests that both thermocline feedback and Ekman feedback are mainly responsible for the rapid increase in temperature anomalies during the developing and mature phases of ENSO events.

El Ni?o, eddy resolving, OGCM, ENSO feedback

1 Introduction

An OGCM with a horizontal grid distance of less than (1/10)° is called an eddy-resolving OGCM, because it directly represents mesoscale eddies in most regions of the global ocean. An eddy-resolving OGCM can better capture the complex topography of the sea floor and the land-sea distribution. In addition, it can better describe the western boundary currents and other narrow currents. Consequently, it is important to develop eddy-resolving OGCMs, to better understand oceanic dynamics, as welltheir effects on climate change.

Eddy-resolving OGCMs have been widely used to investigate mesoscale eddies and the corresponding effects on climate in the North Atlantic basin (Smith et al., 2000; Oschlies, 2002). To date, several studies have diagnosed the ability of global eddy-resolving models in reproducing mesoscale eddies (Shriver et al., 2007; Thoppil et al., 2011). However, simulations of equatorial Pacific temperature and currents associated with ENSO in eddy-resolving OGCMs have rarely been analyzed in detail. In particular, we do not know how increased model resolution improves ENSO simulation. Meanwhile, despite being forced with observed wind stress and heat flux, many coarse resolution OGCMs exhibit very good performance in reproducing some basic ENSO characteristics, but few studies have explored the dynamic feedbacks in these models.

In the ocean-atmosphere system of the Pacific, many studies have used the Bjerknes index to evaluate the simulation ability of ENSO (Kim and Jin, 2011a, b). In particular, there are five linear feedback processes associated with the Bjerknes index: zonal advection feedback, thermocline feedback, Ekman feedback, mean advection feedback, and thermodynamic feedback. Thermocline feedback is the advection of anomalous temperature by mean vertical current, while Ekman feedback is the advection of mean temperature by anomalous vertical current. In the present study, we use the observed wind stress to force two versions of an ocean model. Thus, the ocean dynamic feedbacks are mainly discussed in this paper. The zonal advection feedback, thermocline feedback, and also the Ekman feedback play a positive role. However, the mean advection feedback and thermodynamic feedback have a damping effect on the Bjerknes index. Each feedback depends on the climatological mean state as well as the responses of the ocean (atmosphere) to the atmosphere (ocean). The present study analyzes the above important linear feedback processes related to ENSO.

The present paper attempts to answer two important questions: Does the simulation of equatorial Pacific temperature and currents improve as the model resolution is increased? And if so, what internal mechanism is for the impr responsible ovement?

2 Model and data

The low resolution model (hereafter referred to asLICOM_L) is the LASG/IAP (State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics/Institute of Atmospheric Physics) Climate System Ocean Model, version 2 (LICOM2.0), which has been widely applied in many studies (e.g., Yu et al., 2011; Liu et al., 2012). The model is built using an Arakawa B grid. In particular, its dynamical framework is based on a latitude-longitude grid system with a 1° × 1° horizontal resolution. However, the meridional resolution between 10°S and 10°N increases from 1° to 0.5°, to better resolve the equatorial waveguide with an acceptable computational cost. The grid distance between 10° and 20° varies gradually from 0.5° to 1°. In addition, there are 30 layers with 15 equal-depth levels in the upper 150 m. The model is configured with realistic topography, except that the North Pole is set up as an isolated island. More details can be found in Liu et al. (2012).

The high resolution model (hereafter referred to as LICOM_H) is a quasi-global eddy-resolving OGCM based on LICOM2.0 but with many updates. For instance, the horizontal resolution is higher, at (1/10)° compared to 1° in LICOM2.0; there are more vertical layers (55, compared to 30), the thickness of the first layer is 5 m, there are 36 uneven layers in the upper 300 m, and thus the mean thickness is less than 10 m; and the model domain of LIOCM_H covers 66°N-79°S, meaning the Arctic Ocean is excluded. More details can be found in Yu et al. (2012). In addition, LICOM_H can simulate the observed meridional overturning circulation and meridional heat transport well (Mo and Yu, 2012). The performance regarding the representation of the Indonesian Throughflow, especially the mean vertical structures of the along-strait velocities, is improved in LICOM_H (Feng et al., 2013). Moreover, the outputs from LICOM_H have been successfully used to analyze the deep circulations of the South China Sea (Xie et al., 2013), as well as the eddy energy sources and sinks (Yang et al., 2013). Based on the control experiments in Yu et al. (2012), H. L. Liu (2012, personal communication) recently carried out longer term simulations with the two versions of LICOM (LICOM_H and LICOM_L) and using the same datasets from the Coordinated Ocean-Ice Reference Experiments (Large and Yeager, 2004) for the period 1948-2007. The model data used in the present study are for the period January 1958 to December 2001.

The ocean reanalysis data from the Simple Ocean Data Assimilation (SODA) (Carton and Giese, 2008) are compared with the model simulations. The ocean model of SODA has a horizontal resolution of 0.5° × 0.5° and 40 vertical levels. The monthly dataset used in this study is that of temperature from SODA2.0.2. The climatology fields are defined as the average over the period 1958-2001.

3 Results

Climatological mean temperature in the equatorial upper ocean averaged over the whole 44 years from SODA, LICOM_H, and LICOM_L are shown in Fig. 1. These plots demonstrate that both LICOM_H and LICOM_L can capture the basic structure of the observed temperature, including the depth and zonal tilt of the equatorial thermocline, which is deeper in the west and shallower in the east. However, some discrepancies nevertheless exist, as indicated by the contours, which represent the temperature differences between the model simulations and reanalysis data. It is important to note that the simulated temperature in LICOM_H coincides well with the reanalysis in the eastern equatorial Pacific, especially the region between 100°W and 80°W, but there is an excessive warm bias in this region in LICOM_L. This illustrates that LICOM_H has a better simulation ability than LICOM_L in the eastern basin. In the western equatorial Pacific, both models demonstrate simulation bias. Moreover, it is important to note that the temperature bias in the off-equator region is more obvious than that at the equator (not shown). Thus, the temperature bias in the western equatorial Pacific could be mostly a result of the simulation in the off-equator region—an idea that will be analyzed in a future study. The present study pays attention to the simulation improvement in the eastern equatorial Pacific in LICOM_H.

Figure 1 Time-mean temperature (color scale) in the equatorial upper ocean (5°N-5°S) in (a) Simple Ocean Data Assimilation (SODA), (b) High Resolution Model, and Low Resolution Model, and the corresponding temperature differences (contours) between (b) High Resolution Model and SODA and (c) Low Resolution Model and SODA. Units: °C.

Figure 2 The temperature standard deviation (color scale) in the equatorial upper ocean (5°N-5°S) in (a) SODA, (b) High Resolution Model, and (c) Low Resolution Model (units: °C). The coefficients (contours) denote the response of the temperature anomaly (in 5°N-5°S) to the zonal wind stress anomaly (in (5°N-5°S, 120°E-90°W)) (units: °C (N m-2)-1).

So, why are the temperature biases much improved in the eastern Pacific in LICOM_H? Because of the identical wind stress forcing used in the two models, the response of temperature to winds must be due to dynamic processes associated with the different spatial resolutions. Figure 2 presents the spatial pattern of temperature standard deviation in the equatorial upper ocean (color scale). In the eastern Pacific, the maximum temperature standard deviation is 2.4°C in the reanalysis and 2.4°C in LICOM_H, but 1.6°C in LICOM_L. This illustrates that the standard deviation in LICOM_H is more reasonable, which confirms again the significant improvement in the eastern ocean in LICOM_H. However, some biases also exist. The value of temperature standard deviation in the central Pacific is 2.8°C, which is a bit larger than that in the reanalysis (2.4°C). For the western Pacific, the standard deviation in LICOM_H around 160°E is larger than the reanalysis too. The overestimated standard deviation in LIOCM_H in the western Pacific may be attributable to the same cause of the bias in the mean state, mentioned above. Meanwhile, the contours in Fig. 2 present the regression coefficients of the temperature anomaly to the zonal wind stress anomaly in SODA, LICOM_H, and LICOM_L. The largest positive regression coefficient is 200°C (N m-2)-1located at (60 m, 110-100°W) in SODA, 180°C (N m-2)-1located at (60 m, 120-100°W) in LICOM_H, and 120°C (N m-2)-1located at (60 m, 100-80°W) in LICOM_L. Both the strength and location of the maximum regression in LICOM_H are very similar to those in the reanalysis. This implies that the response of temperature to wind in the eastern equatorial Pacific is weaker in LICOM_L. It is worth noting that the largest negative coefficient in LICOM_H, located at (150 m, 160°E), is larger than in the reanalysis. This is due to the shallower mean thermocline in the western equatorial Pacific in LICOM_H, which would cause a more sensitive response of temperature to the wind stress anomaly. Although a certain degree of bias exists in the western basin, the simulated regression of temperature on wind in the eastern equatorial Pacific in LICOM_H is mostly close that in the reanalysis. This finding is manifested in the simulation improvement in the east, which is in line with the results shown in Fig. 1.

To better understand why the subsurface temperature anomalies are much improved in LICOM_H, we conduct a heat budget analysis using monthly mean temperature and currents in the region (0-100 m, (150-90°W, 5°N-5°S)), where the maximum variability of temperature is observed. We base the analysis on the following temperature equation, as suggested by An and Jin (2004), for describing the heat budget of the ocean subsurface layer, in which the heat flux and subgrid-scale contributions (e.g., small-scale oceanic diffusion, heat flux due to tropical instability waves) are attributed to the residual termR':

where the overbar denotes long-term mean temperatureand velocity, and the prime denotes the anomaly relative to the mean state. The first set of brackets includes the linear advection terms, while the second set presents the non-linear terms. The present study focuses mainly on linear advection. Figure 3 shows the heat budget results, in which the selected El Ni?o years (1965, 1972, 1976, 1982, 1986, 1991, and 1997) are the seven warmest years in the Ni?o3 region averaged in December to February (DJF) and, likewise, the chosen La Ni?a years (1970, 1973, 1975, 1984, 1988, 1998, and 1999) are the seven coldest years. In this way, we have the same number of events for both types, meaning the analysis is more reasonable and reliable. In addition, the data in the ENSO developing phase (from June to December) are used for the analysis. The plots demonstrate that the most important terms areand. The termreflects the effect of zonal advection feedback (Jin and An, 1999). Furthermore, the contribution of the vertical advection of temperature,(representative of the thermocline and Ekman feedback, respectively), is in line with many previous studies (Jin, 1997; Jin and An, 1999; Huang et al., 2012). It is clear that the terms' andpossess the most significant differences between LICOM_H and LICOM_L, and thus both play an important role in inducing the simulation improvement in the eastern equatorial Pacific. The physical mechanism involved is described below.

As discussed above, the vertical advection terms of' andare the key factors in the improved simulation of ENSO. Furthermore, in the upper ocean,?zT' in LICOM_H is closer to the observation than that in LICOM_L, due to the increased vertical resolution—clearly seen in the results presented in Figs. 1 and 2. Thus, the differences of w' andbetween LICOM_H and LICOM_L need to be detected in detail. Given the fact that the currents dataset in SODA suffers from some severe biases, we use an empirical formula to diagnose thezonal and meridional velocity, and then calculate w′ and. The anomalous Ekman currents are calculated based on Chang and Philander (1994), using

Figure 3 Heat budget analysisaveraged over the Ni?o3 region (5°N-5°S, 150-90°W) in the upper ocean (0-100 m) for the warm phase (red) and cold phase (blue) in (a) High Resolution Model and (b) Low Resolution Model. Units: °C/month.

where β is the planetary vorticity gradient, ρ is the density of seawater; H1is the mean mixed-layer depth (here, we use 50 m); τxand τyare the zonal and meridional wind stress anomalies, respectively; and rs= 0.5 d-1is the dissipation rate (Su et al., 2010). Figure 4 exhibits the anomalous Ekman currents, together with the w′ averaged over the Ni?o3 region simulated in LICOM_H and LICOM_L. The results show that the temporal evolution and magnitude of w′ in LICOM_H are more similar to the anomalous Ekman currents than they are in LICOM_L. This means that LICOM_H simulates a more reasonable w′ than LICOM_L, which could contribute to the Ekman feedback.

Keepinxgthe oytherparameters as theyare in Fig.4, we replace τ and τ with monthly-meanwind stress in the above equation to obtain the monthly-mean vertical velocity. Figure 5 plots thesimulated in LICOM_H and LICOM_L in comparison with the climatological mean Ekman current. A double-peak pattern exists in LICOM_H, which is identical to the mean Ekman current. Furthermore, the magnitude ofin LICOM_H is more reasonable according to the three time series. Therefore, the improved simulation ofin LICOM_H could also result in the improved thermocline feedback. In addition, the strengthened climatological mean upwelling is associated with the improvement in the climatological mean thermocline, and thus the more reasonable w also contributes to the Ekman feedback. Besides their different spatial resolutions, it is necessary to emphasize two other main differences between LICOM_H and LICOM_L: the horizontal viscosity coefficient is different in the two versions of LICOM, and the parameterization of mesoscale eddies (Gent and McWilliams, 1990) is turned off in the equations of tracers in LICOM_H. Considering the above findings, the weakened horizontal viscosity coefficient is the key reason, as it results in stronger horizontal and vertical currents.

Figure 4 The anomalies of Ekman vertical currents (black line), the anomalies of vertical currents in High Resolution Model (red line), and the anomalies of vertical currents in Low-Resolution Model (green line), averaged over the Ni?o3 region (5°N-5°S, 150-90°W) in the upper ocean (0-50 m). Units: m s-1.

4 Summary

Because the simulation of the equatorial Pacific temperature and currents associated with ENSO in LICOM_H has rarely been explored in previous studies, we conduct two numerical experiments, one using LICOM_H and one using LICOM_L, to help understand how the finer resolution of LICOM_H helps to improve the simulation of ENSO. Our main findings can be summarized as follows:

(1) The finer resolution improves the simulation of the mean state of, for example, the vertical temperature gradient in the upper ocean. Also, the vertical velocity nearthe surface is comparable with the Ekman current estimated from the observed wind stress.

Figure 5 The monthly-mean vertical currents (black line), the monthly-mean vertical currents in High Resolution Model (red line), and the monthly-mean vertical currents in Low-Resolution Model (green line), averaged over the Ni?o3 region (5°N-5°S, 150-90°W) in the upper ocean (0-50 m) (units: m s-1).

(2) The vertical advections of temperature in the temperature equation are improved due to the contributions of both temperature gradient and velocity.

(3) Through heat budget analysis, the differences in vertical advections are much greater than those of zonal and meridional advection between LICOM_H and LICOM_L. This is manifested in the key roles played by the thermocline feedback and Ekman feedback in causing the amplitude of ENSO.

(4) Although the external wind stress forcing is exactly the same in the LICOM_H and LICOM_L experiments, the regression of temperature on wind stress is better reproduced in LICOM_H compared to that in LICOM_L, due to the above-mentioned more reasonable ENSO feedbacks.

It is important to acknowledge that there are some other details related to the high resolution of the model that need be further investigated in future studies. For instance, the role of nonlinear temperature advections, vertical mixing, and the effect of tropical instability waves are not considered in the present study. In addition, an attempt should be made to better understand the temperature bias in the western equatorial Pacific in LICOM_H.

Acknowledgements. This study was jointly supported by the National Basic Research Program of China (Grant No. 2013CB956204), the “Strategic Priority Research Program—Climate Change: Carbon Budget and Relevant Issues” of the Chinese Academy of Sciences (Grant No. XDA05110302), and the Jiangsu Collaborative Innovation Center for Climate Change.

References

An, S. I., and F. F. Jin, 2004: Nonlinearity and asymmetry of ENSO, J. Climate, 17, 2399-2412.

Carton, J. A., and B. S. Giese, 2008: A reanalysis of ocean climate using Simple Ocean Data Assimilation (SODA), Mon. Wea. Rev., 136, 2999-3017.

Chang, P., and S. G. Philander, 1994: A coupled ocean-atmosphere instability of relevance to the seasonal cycle, J. Atmos. Sci., 51, 3627-3648.

Feng, X., H. L. Liu, F. Wang, et al., 2013: Indonesian throughflow in an eddy-resolving ocean model, Chin. Sci. Bull., 58, 4504-4514.

Gent, P. R., and J. C. McWilliams, 1990: Isopycnal mixing in ocean circulation models, J. Phys. Oceanogr., 20, 150-155.

Huang, B. Y., Y. Xue, H. Wang, et al., 2012: Mixed layer heat budget of the El Ni?o in NCEP climate forecast system, Climate Dyn., 39, 365-381.

Hcy屬于蛋氨酸循環(huán)的中間產(chǎn)物,且臨床上已有不少研究報道證實,血漿Hcy水平的升高可促進動脈粥樣硬化,同時促進動脈以及靜脈血栓形成,從而導致心腦血管疾病的發(fā)生[5-8]。

Jin, F. F., 1997: An equatorial ocean recharge paradigm for ENSO. Part I: Conceptual model, J. Atmos. Sci., 54, 811-829.

Jin, F. F., and S. I. An, 1999: Thermocline and zonal advection feedbacks within the equatorial ocean recharge oscillator model for ENSO, Geophys. Res. Lett., 26, 2989-2992.

Kim, S. T., and F.-F. Jin, 2011a: An ENSO stability analysis. Part I: Results from a hybrid coupled model, Climate Dyn., 36, 1593-1607.

Kim, S. T., and F.-F. Jin, 2011b: An ENSO stability analysis. Part II: Results from the twentieth and twenty-first century simulations of the CMIP3 models, Climate Dyn., 36, 1609-1627.

Large, W. G., and S. Yeager, 2004: Diurnal to Decadal Global Forcing for Ocean and Sea-Ice Models: The Data Sets and Flux Climatologies, Technical Report TN-460+STR, NCAR, 105pp. Liu, H. L., P. F. Lin, Y. Q. Yu, et al., 2012: The baseline evaluation of LASG/IAP Climate system Ocean Model (LICOM) version 2, Acta Meteor. Sinica, 26, 318-329.

Mo, H. E., and Y. Q. Yu, 2012: Simulation of volume and heat transport along 26.5°N in the Atlantic, Atmos. Oceanic Sci. Lett., 5, 373-378.

Oschlies, A., 2002: Improved representation of upper-ocean dynamics and mixed layer depths in a model of the North Atlantic on switching from eddy-permitting to eddy-resolving grid resolution, J. Phys. Oceanogr., 32, 2277-2298.

Shriver, J. F., H. E. Hurlburt, O. M. Smedstad, et al., 2007: 1/32° real-time global ocean prediction and value-added over 1/16° resolution, J. Mar. Syst., 65, 3-26.

Su, J. Z., R. H. Zhang, T. Li, et al., 2010: Causes of the El Ni?o and La Ni?a amplitude asymmetry in the eastern equatorial Pacific, J. Climate, 23, 605-617.

Thoppil, P. G., J. G. Richman, P. J. Hogan, 2011: Energetics of a global ocean circulation model compared to observations, Geophys. Res. Lett., 38, L15607, doi:10.1029/2011GL048347.

Xie, Q., J. G. Xiao, D. X. Wang, et al., 2013: Analysis of deep and bottom circulation in the South China Sea based on eight quasiglobal ocean model outputs, Chin. Sci. Bull., 58, 4000-4011.

Yang, H. Y., L. X. Wu, H. L. Liu, et al., 2013: Eddy energy sources and sinks in the South China Sea, J. Geophys. Res., 118, 4716-4726.

Yu, Y. Q., H. Liu, and P. Lin, 2012: A quasi-global (1/10)° eddy-resolving ocean general circulation model and its preliminary results, Chin. Sci. Bull., 57, 3908-3916.

Yu, Y. Q., W. P. Zheng, B. Wang, et al., 2011: Versions g1.0 and g1.1 of the LASG/IAP Flexible Global Ocean-Atmosphere-Land System Model, Adv. Atmos. Sci., 28, 99-117.

Hua, L.-J., and Y.-Q. Yu, 2015: How are El Ni?o and La Ni?a events improved in an eddy-resolving ocean general circulation model? Atmos. Oceanic Sci. Lett., 8, 245-249,

10.3878/AOSL20150015.

8 January 2015; revised 9 March 2015; accepted 3 April 2015; published 16 September 2015

YU Yong-Qiang, yyq@lasg.iap.ac.cn

猜你喜歡
血漿
糖尿病早期認知功能障礙與血漿P-tau217相關性研究進展
血漿置換加雙重血漿分子吸附對自身免疫性肝炎合并肝衰竭的細胞因子的影響
血漿corin、NEP、BNP與心功能衰竭及左室收縮功能的相關性
血漿B型利鈉肽在慢性心衰診斷中的應用
miRNA-145和miRNA-143在川崎病患兒血漿中的表達及意義
CHF患者血漿NT-proBNP、UA和hs-CRP的變化及其臨床意義
自體富血小板血漿在周圍神經(jīng)損傷修復中的潛在價值
胰腺癌患者血漿中microRNA-100水平的測定及臨床意義
腦卒中后中樞性疼痛相關血漿氨基酸篩選
系統(tǒng)性硬化病患者血漿D-dimer的臨床意義探討
主站蜘蛛池模板: 国产福利一区在线| 手机永久AV在线播放| 凹凸国产分类在线观看| 久久毛片基地| 欧美成a人片在线观看| 中文字幕亚洲无线码一区女同| 国产女同自拍视频| 亚洲一区无码在线| 国产成人超碰无码| 国模沟沟一区二区三区| 91福利在线看| 国产h视频免费观看| 亚洲区视频在线观看| 美女被操黄色视频网站| 97国产精品视频自在拍| 99热国产这里只有精品9九| 97se亚洲综合在线天天| 久久亚洲日本不卡一区二区| 免费a在线观看播放| 国产综合精品日本亚洲777| 欧美成人亚洲综合精品欧美激情| 色婷婷久久| 国产波多野结衣中文在线播放| 国产亚洲精品自在线| 亚洲婷婷丁香| 又猛又黄又爽无遮挡的视频网站| 亚洲国产日韩视频观看| 国产小视频在线高清播放| 99久久国产综合精品2023| 国产区网址| 亚洲第一极品精品无码| 免费看黄片一区二区三区| 99这里精品| 韩国v欧美v亚洲v日本v| 国产精品网拍在线| 欧美精品二区| 成人国内精品久久久久影院| 欧美成人午夜影院| 国产精品久久久久久搜索| 美女无遮挡免费视频网站| 99性视频| 国产va免费精品| aaa国产一级毛片| 扒开粉嫩的小缝隙喷白浆视频| 91综合色区亚洲熟妇p| 国产成人精品视频一区二区电影| 免费av一区二区三区在线| 波多野结衣国产精品| 91福利免费| 亚洲男女天堂| 色婷婷在线影院| 尤物成AV人片在线观看| 无码免费的亚洲视频| 国产精品视频第一专区| 国产永久免费视频m3u8| 婷婷色在线视频| 久久无码av三级| 亚洲色图狠狠干| 久热这里只有精品6| 欧美不卡视频在线| 福利姬国产精品一区在线| 18禁色诱爆乳网站| 欧美三级日韩三级| 欧美一区精品| 午夜老司机永久免费看片| 国产精品蜜芽在线观看| 国产精品第三页在线看| 高清久久精品亚洲日韩Av| 浮力影院国产第一页| 乱码国产乱码精品精在线播放| 熟女日韩精品2区| 国产综合无码一区二区色蜜蜜| 亚洲另类第一页| 国产一线在线| 亚洲成a人在线播放www| 综合人妻久久一区二区精品 | 欧美精品亚洲二区| av午夜福利一片免费看| 99色亚洲国产精品11p| 九九九九热精品视频| 99久久亚洲综合精品TS| 乱系列中文字幕在线视频|