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

    Nonlinear Responses of Oceanic Temperature to Wind Stress Anomalies in Tropical Pacific and Indian Oceans:A Study Based on Numerical Experiments with an OGCM

    2015-01-05 02:01:57HUALijuan華莉娟andYUYongqiang俞永強
    Journal of Meteorological Research 2015年4期

    HUA Lijuan(華莉娟)and YU Yongqiang(俞永強)

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

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

    Nonlinear Responses of Oceanic Temperature to Wind Stress Anomalies in Tropical Pacific and Indian Oceans:A Study Based on Numerical Experiments with an OGCM

    HUA Lijuan1,2(華莉娟)and YU Yongqiang2?(俞永強)

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

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

    As a highly nonlinear dynamic system,oceanic general circulation models(OGCMs)usually exhibit nonlinear responses to prescribed wind stress forcing.To explore mechanisms for these nonlinear responses, we designed and conducted three idealized numerical experiments with an OGCM with modified wind stress forcing.In the experiments,the climatological mean wind stress was identical,and the only differences in external forcing were wind stress anomalies.The wind anomalies were set to zero in a control run,and the observed wind stress anomalies with and without reversed signs were superimposed on the mean climatology in two sensitivity experiments.Forced by the prescribed wind stress anomalies in sensitivity runs,the OGCM well reproduced the El Ni?no-Southern Oscillation(ENSO)and the Pacific and Indian Ocean Dipole(IOD) in the Indian Ocean,as wellas the asymmetry between positive and negative phases of these modes.Relative to the control run,the two sensitivity runs exhibited almost identical changes in the mean climate state, although the wind stress anomalies were reversed in these two experiments.Thus,it was concluded that the asymmetry of wind stress anomalies contributes only slightly to the mean state changes and ocean internal dynamics was the main contributor.Further heat budget analysis suggested that nonlinear temperature advection terms,including both mean advection and perturbed advection,favor the ENSO/IOD rectified effect on the mean state.

    ENSO,IOD,numerical simulation,asymmetry,rectification,nonlinearity

    1.Introduction

    It is well known that the probability density function of the observed tropical Pacific sea temperature does not always follow a standard normal distribution(Burgers and Stephenson,1999;An and Jin, 2004;Zhang et al.,2009).Instead,it shows a significantly skewed distribution in most regions.For example,in the equatorial central and eastern Pacific,both the strength and duration of warm temperature anomalies in El Ni?no events are markedly different from those in their counterparts;this difference is commonly referred to as El Ni?no-Southern Oscillation(ENSO)asymmetry(Zebiak and Cane,1987; Okumura and Deser,2010).Further studies indicated that the asymmetry of the temperature anomalies is closely related to nonlinear dynamical heating(NDH) (An and Jin,2004;An et al.,2005;Su et al.,2010). Furthermore,ENSO asymmetry can be caused by high-frequency winds through a nonlinear rectification process,which increases the skewness of sea surface temperature(SST)anomalies in the eastern equato-rial Pacific(Rong et al.,2011).Other researchers have reported that the atmospheric responses to SST anomalies of the same magnitudes but opposite signs are asymmetric,which might,in turn,reinforce the original SST asymmetry(Ohba and Ueda,2009). The interaction between ENSO and the climatological mean state is another important issue,but the specifics of this interaction remain unknown.Through numerical experiments with simplified or fully coupled models,many studies have emphasized the impact of the mean state on ENSO events,especially in modulations of ENSO period and amplitude(An and Jin,2001;Fedorov and Philander,2001).The climatological basic state of the tropical atmosphere can influence ENSO events by modulating westerly winds over the tropical Pacific(Yan and Zhang,2002).However,recent studies have begun to consider the time-mean effect of ENSO events(Sun et al.,2012;Sun et al.,2014;Hua et al.,2015).Sun(2003)found that stronger ENSO variations attempt to cool the western Pacific and warm the eastern Pacific.Sun and Zhang(2006)further noted that the time-mean effect of ENSO would cool the warm pool while warming the thermocline water and broad region in the surface eastern Pacific.Some of these studies have been empirically based(Rodger et al.,2004;Choi et al.,2012).These studies have suggested that decadal variability in the tropical Pacific SST probably results from the residual effect of ENSO.

    The asymmetry of temperature anomalies is not only observed in the tropical Pacific,but also in the Indian Ocean.The Indian Ocean Dipole(IOD)is a zonal mode of SST interannual variation(Saji et al.,1999).Observational(Saji and Yamagata,2003), modeled(Cai et al.,2005),and theoretical(Zhong et al.,2005)studies have indicated that the IOD is a seasonally dependent mode that develops rapidly in summer and reaches a mature phase in autumn. IOD events are significantly asymmetric,which is similar to ENSO events,and the strength of a positive IOD is greater than that of a negative IOD. Meanwhile,strong negative IOD temperature skewness anomalies exist in the eastern IOD in autumn (Hong et al.,2008).Hong et al.(2008)considered that the asymmetry in wind-advection-SST feedback, SST-cloud-radiation feedback,and wind-evaporation-SST feedback all contribute to this negative skewness. In particular,wind-advection-SST feedback is related to NDH.It shows that the NDH favors cold anomalies, rather than restraining warm anomalies,in the eastern IOD heat budget(Hong et al.,2008;Cai and Qiu, 2013).In addition,Zheng et al.(2010)noticed,using coupled model simulations,that SST-thermocline depth feedback is asymmetric.Furthermore,this asymmetry might induce the negative skewness of temperature anomalies in the equatorial southeastern Indian Ocean.Ogata et al.(2013)confirmed the above assumption through observations and an ocean general circulation model(OGCM)study.However, the rectified effect of IOD asymmetry on the climatological mean state has not been considered until now.

    Although many climate models,including standalone OGCMs and coupled general circulation models,have succeeded in reproducing the skewness of temperature anomalies,the mechanism of temperature skewness remains unclear.In this study,we conduct numerical experiments with an OGCM and then evaluate how well the model reproduces the observed asymmetry of temperature anomalies in two ocean basins.Many studies have shown a close relationship between temperature asymmetry and climatological mean state;however,it is difficult to identify a physical mechanism and causality.In our study,the OGCM willbe forced by the climatological mean wind and anomalous wind on an interannualtimescale.The simulated mean state is clearly different between the two experiments.Therefore,we can compare the differences between the two experiments and further explore how the physical mechanism influences the temperature asymmetry and mean state.Subsequently, the numerical experiments would help us comprehend the asymmetry and rectification resulting from the external forcing or ocean dynamics.

    2.Methodology

    2.1 Model

    The OGCM used in this study is a climate sys-tem ocean model,LICOM2.0(Liu et al.,2004a;Liu et al.,2012),which was developed at the State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics,Institute of Atmospheric Physics.This model has been widely applied in many studies(Liu et al.,2004b;Yu et al., 2011;Li et al.,2012).It is a global ocean model that covers the region from 75°S to 88°N,and the North Pole is modeled as an isolated island.In particular,the model’s dynamical framework is based on a latitudelongitude grid system with 1°×1°horizontal resolution;however,the meridional resolution increases to 0.5°between 10°S and 10°N,and there are 30 layers with 15 equal-depth levels in the upper 150 m. The model adopts some new physical parameterization schemes,such as a new turbulent mixing scheme (Canuto and Dubovikov,2005),solar radiation penetration,and an improved isopycnal-mixing scheme (Gent and Williams,1990;Large et al.,1997).It can well reproduce the basic temporal and spatial structures of the observed ENSO and IOD,especially the significant asymmetry between El Ni?no and La Ni?na events,as well as positive and negative IOD events (Hua et al.,2010).By using climatological monthlymean wind stress,LICOM2.0 was first integrated from a motionless state to an equilibrium state through a 500-yr spin up.The end state of year 500 was then used as the initial condition for the numerical experiments that will be described below.

    2.2 Experiment design

    In the two sensitivity runs,the surface wind stresses were climatological winds plus the interannual monthly anomalies,and the surface boundary conditions for temperature and salinity were the same as those for the control run.In the first sensitivity run(hereafter referred to as WSA),which was the same as the sensitivity run described by Sun et al.(2014),the interannual monthly anomalies during 1958-2001(obtained from the European Centre for Medium-Range Weather Forecasts 40-yr reanalysis(ERA40;Uppala et al.,2005))over the tropical Pacific(30°N-30°S,120°E-80°W)and Indian Ocean (30°N-30°S,30°-120°E)were superimposed on the climatological surface wind stress.The second sensitivity run(hereafter referred to as WSA?R)was the same as WSA except that the signs of the interannual monthly anomalies were reversed.Notably,because of the restoring boundary conditions,the simulated surface temperature and salinity were restored to the climatologicalmean value with a seasonalcycle.Then, as a negative feedback,the restoring boundary conditions would have a damping effect on any perturba-tions in SST and salinity that deviated from the climatological mean state.Consequently,the anomalies in temperature and other variables and the skewness of the anomalies could only result from the wind stress forcing or the ocean dynamical process in the numerical experiments.

    According to the experimental design,the simulated temperature asymmetry,as well as the ENSO/ IOD rectification,should originate only from the external forcing and internal nonlinear dynamics in the ocean.Comparison between the two sensitivity runs may help to quantify the relative contributions of the above two factors.The ocean reanalysis data from simple ocean data assimilation(SODA;Carton and Giese,2008)were compared with model simulations. The ocean model of SODA has a horizontal resolution of 0.5°by 0.5°and 40 vertical levels.The monthly dataset used in this study was the temperature data from SODA2.0.2.Climatology fields were defined as the averages for the period 1958-2001.

    2.3 Diagnosis method

    The diagnosis methods used in the paper follow:

    2)The time mean temperature equation can be expressed as

    3.Model verification

    Forced with the observed wind stress and heat flux,the ocean model well captured the basic characteristics of the observed mean state.Figure 1 presents the climatological mean temperature and its annual standard deviation in the equatorial upper ocean(2°N-2°S)as a function of longitude and depth in the two basins,in which both the simulation and reanalyses were used from 1958 to 2001.These plots demonstrate that the model could capture the basic structure of the observed temperature,including the depth and tilt of the equatorial thermocline.In the Pacific,the model succeeds in describing the spatial pattern of the thermocline,which is deeper in the west and shallower in the east,and the large-scale pattern of the simulated temperature coincides well with the reanalysis;however,the simulated temperature near the western boundary is approximately 1℃higher than that obtained from the reanalysis.Meanwhile,an excessive cold tongue and cold bias are shown in the equatorial eastern Pacific in WSA.In the Indian Ocean,the model also well reproduces the spatial structure of the observed thermocline,which is shallower in the western basin and deeper in the easternbasin.Both the observed and simulated tilts of the thermocline are opposite to those in the Pacific because the annual averaged wind stress prevailing in the equatorial Indian Ocean is westerly,but that in the Pacific is easterly.The simulated depth of the mean thermocline is almost 120 m,which agrees well with the reanalysis.However,there still exists some model bias;for instance,the simulated vertical temperature gradient is smaller than that obtained from the reanalysis.In general,the overall pattern sim-ulated by LICOM2.0 resembles the reanalysis in the two basins(implying that the model has enough simulation ability),and a few discrepancies,such as cold bias in the equatorial eastern Pacific,are attributed to the coarse resolution of the model.The simulated standard deviation in the tropical Pacific agrees well with the reanalysis except for the smaller magnitude along the thermocline.In the Indian Ocean,the spatial pattern of the standard deviation is basically similar to the reanalysis,although there are differences in magnitude.

    Fig.1.Temperature(℃;contour)and annual standard deviation(>0.3 with color)in the equatorial upper Pacific and Indian Oceans(2°N-2°S)in(a,b)SODA,(c,d)WSA,and(e,f)the corresponding temperature(contour)and standard deviation(color)differences between WSA and SODA.

    Fig.2.(a,b)Nino3 and IOD index,(c,d)temperature anomalies in subsurface of eastern Pacific(2°N-2°S,100°W; 100 m)and eastern Indian Ocean(2°N-2°S,95°E;60 m),and(e,f)temperature anomalies in subsurface of western Pacific(2°N-2°S,160°W;100 m)and western Indian Ocean(2°N-2°S,45°E;120 m)from SODA data(black line),WSA (red line),and WSA?R(green line).

    One of the principal causes of the simulated sea surface temperature anomalies as a response to the observed wind stress is ENSO in the tropical Pacific and IOD in the Indian Ocean.The observed ENSO/IOD variability was well reproduced by the OGCM,as shown in Fig.2,which compares the time series of the Nino3/IOD index from the reanalysis with thatfrom the WSA/WSA?R experiments.The simulated time series of the Nino3 from the WSA/WSA?R indices have correlations of 0.87/-0.79 with those from the reanalysis and 0.66/-0.50 for IOD index.The time series of the Nino3 from WSA has correlations of -0.92 and-0.72 with those from WSA?R and from the IOD index,respectively.The correlation between the time series of Nino3 and the IOD index is 0.36, 0.29,and 0.27 for the reanalysis,WSA,and WSA?R, respectively.The simulated time series from WSA has a similar temporal evolution to that obtained from the reanalysis,but the simulated amplitude is weaker,especially for the IOD index.This discrepancy is expected,given that the thermal boundary conditions are simplified and that the SST is restored to the observed climatological datasets(Levitus and Boyer,1994;Levitus et al.,1994).Moreover,because of the stronger magnitude associated with positive ENSO/IOD phases than with negative phases, the positive asymmetry can be found in both Nino3 and IOD indices(WSA experiment).More details can be seen in Table 1.Both the Nino3 and IOD indices show the temporal evolution at the surface,but the changes in the temperature anomalies in the subsurface also play an important role in ENSO/IOD dynamics.Therefore,the time series of the subsurface temperature anomalies in the eastern/western Pacific and Indian Ocean are presented for SODA,WSA,and WSA?R in Fig.2.In addition,the skewness and standard deviation of subsurface temperature anomalies can be seen in Tables 2 and 3.

    4.Asymmetry

    4.1 Subsurface

    4.1.1 ENSO/IOD asymmetry

    It is wellknown that ENSO is the most important coupled ocean-atmosphere phenomenon in the tropical Pacific(Zebiak and Cane,1987;Battisti,1988)and the first EOF mode of the equatorial subsurface temperature anomalies is ENSO(Cai et al.,2003;Rodger et al.,2004;Choi et al.,2013).Therefore,temperature asymmetry in the equatorial Pacific is mostly induced by ENSO asymmetry.Similarly,Cai et al.(2008) found,using temperature data from Scripps Institute of Oceanography,that the first EOF mode of the subsurface temperature anomalies is described by an eastwest dipole in the equatorial Indian Ocean.With SODA reanalysis data,we reached the same conclusion (figure omitted).In other words,the IOD mode is the first EOF mode for subsurface temperature anomalies,and the subsurface temperature asymmetry in the equatorial Indian Ocean also results from the IOD asymmetry.As suggested by An et al.(2005),skewness is an effective metric to quantify the asymmetryoftemperature anomalies.Figure 3 illustrates temperature skewness as a function of longitude and depth in the upper equatorial ocean(2°N-2°S).The positive (negative)skewness implies that the strength of warm (cold)anomalies is larger than that of cold(warm) anomalies.The skewness from WSA coincides well with that of SODA in the upper ocean.The maximum positive skewness is at 100 m at around 100°W,and the largest negative skewness is at the same depth at around 160°E.However,the positive skewness occurs in the region(170°E-160°W;0-80 m)that is different from SODA;and a small amount of negative skewness around 140°Win the upper 40 m in WSA is not found in SODA.Although there exist some discrepancies,the model could still represent the most significant skewness in the subsurface.In WSA?R,the location ofthe extreme center is not exactly the same as that in WSA,but there is a common phenomenon at a depth of 40-160 m.Meanwhile,the quantitative difference between the skewness from WSA and WSA?R is not negligible.From WSA minus WSA?R(figure omitted),we found that the largest positive value was 2.5 and was located at 100 m at around 100°W,while the largest negative value was-2 at 100 m at around 160°E.This result means that both the positive skewness in the equatorial eastern Pacific and the negative skewness in the western Pacific are stronger in WSA than in WSA?R.

    Table 1.Skewness and standard deviation of time series of Nino3 and IOD index from SODAdata,WSA, and WSA?R experiments

    Table 2.As in Table 1,but for time series of temperature anomalies in subsurface of eastern Pacific(2°N-2°S,100°W;100 m)and eastern Indian Ocean(2°N-2°S,95°E;60 m)

    Table 3.As in Table 1,but for time series of temperature anomalies in subsurface of western Pacific(2°N-2°S,160°W;100 m)and western Indian Ocean(2°N-2°S,45°E;120 m)

    Fig.3.Skewness of temperature in equatorialupper Pacific and Indian Oceans(2°N-2°S)in(a,b)SODA,(c,d)WSA, and(e,f)WSA?R.

    In the Indian Ocean,the skewness results from WSA and WSA?R and those from SODA are similar in the upper 140 m.There exists negative skewness in the upper 100 m of the equatorial eastern basin and positive skewness at 100-140 m.The largest negative skewness is between 90°-100°E at 20-40 m for both simulations and SODA,but the simulated values are larger.The positive skewness for WSA and WSA?R is also greater than that for SODA.The quantitative difference between the skewness results of WSA and WSA?R in the Indian Ocean(figure omitted)revealed that the largest negative skewness was nearly-7 and was located at 140 m between 80°-90°E,with a negative minimum of around-1 at 20-40 m between 90°and 100°E.Thus,considering the skewness difference between WSA and WSA?R,the negative skewness was slightly larger in the upper eastern Indian Ocean in WSA;more significant positive skewness was exhibited below 100 m in the eastern basin in WSA?R.The quantitative difference between WSA and WSA?R cannot be ignored in the two basins,which means that the asymmetry of the wind stress anomalies would cause changes in the magnitude of the temperature skewness.Specifically,the thermocline responses to the reversed wind stress anomalies were asymmetric. The thermocline became shallower when easterly wind stress anomalies occurred and deeper when westerly anomalies ofthe same magnitudes occurred.The more significant positive skewness in WSA?R is attributed to the stronger oceanic responses to the westerly wind stress anomalies.Furthermore,even if the climatological mean wind stress was applied in the control run,the simulated temperature still exhibited some skewness that was entirely different from that obtained from the reanalysis,WSA,and WSA?R(figure omitted).We found that the large-scale spatial pattern in the control run was obviously different from that in other plots,for example,differences in the location of the significant skewness and in the magnitude of the skewness.Considering the above findings,we conclude that the temperature skewness in the upper 40-160 m in the Pacific,and in the upper 140 m in the Indian Ocean,is not associated with the skewness of the wind stress.Our results resemble Figs.6 and 7c in Ogata et al.(2013);the skewness was zero in numerical experiments with wind stress anomalies as the external forcing(Ogata et al.,2013),but the SST response was strongly negatively skewed.The numerical experimental design in our study was not identical to that of Ogata et al.,but the conclusions are similar. Therefore,the temperature asymmetry in the subsurface to the first order may be attributed to the internal dynamical response to the external forcing.

    4.1.2 Role of NDH in ENSO/IOD asymmetry

    Besides the skewness of the temperature anomalies,another metric to quantify asymmetry is the residual of temperature anomalies associated with positive and negative phases of ENSO/IOD events(Rodgers et al.,2004;Su et al.,2010;Ogata et al.,2013;Sun et al.,2014).We performed heat budget analysis in the positive and negative phases of ENSO/IOD events (figure omitted).The results showed that the zonal and meridional terms of NDH enhanced the positive temperature anomalies in the warm events in the upper equatorial eastern Pacific but weakened the negative anomalies in the cold events.The vertical term of NDH contributed greatly to temperature decrease in the warm events in the equatorial eastern Indian Ocean,but mitigated temperature increase in the cold events.Our findings based on the numerical experiments are consistent with those from reanalysis data by Su et al.(2010).

    4.2 Heat budget in the mixed layer

    Figure 4 illustrates the temperature asymmetry in the mixed layer.In the Pacific,negative skewnessexists in the southeastern basin(3°N-5°S,140°-80°W)in WSA?R,but the opposite occurs in WSA. This difference indicates that the skewness of the wind stress anomalies would play a role in the skewness of the temperature anomalies in such a region.However, there is a common negative skewness in the western Pacific(120°-160°E)and in the southeastern tropical Indian Ocean,which confirms the more significant role of ocean dynamics in such regions.

    Fig.4.As in Fig.3,but for mixed layer condition.

    5.Time-mean effect

    5.1 Subsurface

    5.1.1 ENSO/IOD rectification

    Because ofasymmetry,the residualof the positive and negative phases of ENSO/IOD events is eithermuch larger or smaller than zero.The spatial pattern of the residual of ENSO/IOD events is similar to that of the above temperature asymmetry.Given the fact that warm and cold temperature anomalies cannot completely cancel each other out,the long-term mean ENSO/IOD events should influence the climatological mean state or decadaloscillation,as has also been suggested by some recent studies(Rodger et al.,2004;Yu and Kim,2011;Choiet al.,2012).As shown in Fig.5, although the climatological mean wind stress used in all numericalexperiments was identical,the simulated climatological mean temperature differs between the sensitivity experiments and control run.In the Pacific,due to interannual wind stress forcing,the simulated mean temperature in both WSA and WSA?R is 1℃colder than that of the control run in the equatorial western Pacific,while the warm anomalies show a maximum of around 1℃in the equatorial eastern Pacific.In the Indian Ocean,both WSA and WSA?R show temperature decreases in the upper 80 m in theequatorial eastern basin and temperature increases along the thermocline at 80-160 m in the Indian Ocean.In terms of the subsurface,the large-scale spatial pattern in Fig.5 is similar to that of Fig.3,but the locations of the maximum/minimum values in Fig. 5 are not exactly identical to those in Fig.3.In Fig. 5,the time-mean effect from WSA has a high spatialcorrelation with WSA?R,implying that the mean temperature change in the subsurface is due to the asymmetry of the oceanic response to the anomalous surface forcing.Besides the temperature change,the current and salinity differences between the sensitivity experiments and control run showed identical behaviors(figure omitted).Furthermore,we conducted an additionalexperiment to confirm the above results. The experiment was designed to force the ocean model using idealized interannual wind stress anomalies but without any asymmetry.The time-mean temperature differences between the simulation results from this experiment and those from the control run were consistent with Fig.5.

    Fig.5.Temperature differences between control run and(a,c)WSA and(b,d)WSA?R in equatorial upper Pacific and Indian Oceans(2°N-2°S).

    5.1.2 Role of NDH in ENSO/IOD rectification

    Because the observed wind stress anomalies revealed obvious interannual variations,there existed dramatic interannual variations in the simulated temperatures and currents of the sensitivity experiments; however,these were negligible in the controlrun.Consequently,the climatological mean nonlinear advection terms in the sensitivity runs were different from those in the control run.Liang et al.(2012)and Sun et al.(2014)noted that ENSO rectification of the mean state is induced by the nonlinear heating terms of the heat budget formula.Therefore,we present NDH differences between the sensitivity runs(WSA and WSA?R)and the control run(Fig.6),respectively.The results in the Figs.6a and 6c are similar to those in Figs.6b and 6d for the two basins,indicating that the situation is concurrent with Fig.5. Only a minor difference between WSA and WSA?R can be found,in which the maximum difference is limited to 0.1℃mon?1.In the Pacific,the largest positive NDH occurs at 20-60 m between 100°and 80°W, but the negative NDH mainly exists in the upper 80 m between 120°and 160°E.The warming center is located at 20-60 m between 100°and 80°W in Fig.5, while the cooling center at 120 m is around 160°E. Thus,the mean temperature change in the equatorial eastern Pacific(Fig.5)could be explained by NDH,but the temperature decrease in the equatorial western Pacific could not.In the Indian Ocean,the largest negative NDH is located at 20-80 m between 85°and 98°E,and the largest positive NDH is at 80-120 m between 70°and 80°E.The cooling center is located at 20-80 m between 85°and 98°E in Fig.5, and the warming center is at 80-120 m between 70°

    and 100°E.Therefore,NDH could explain the temperature decrease in the Indian Ocean,as well as in the partial warming region.To sum up,NDH favors the rectified effect of temperature asymmetry on the mean state.Further research(Hua et al.,2015)has shown that the upper ocean temperature and velocity anomalies(T?,u?,v?,w?)respond linearly to wind stress anomalies.Therefore,the significant correlation between the temperature anomaly(T?)and velocity anomalies(u?,v?,w?)(and thus the NDH)retains the same sign even though the wind stress anomalies are reversed.In this way,the impact of the rectification on the mean state remains almost unchanged. We further compared the three components of the NDH with the total NDH(figure omitted)and found that the most important components were the zonal and meridional in the Pacific Ocean,while the key component was the vertical in the Indian Ocean.

    5.1.3 Role of mean advection in ENSO/IOD rectification Given the fact that NDH could not account for the mean state change in some regions,such as the western Pacific,Fig.7 shows the differences in mean advection between the sensitivity runs(WSA and WSA?R)and control run.In the Pacific,the cold advection corresponding with the cold temperature anomalies covers the western region(140°E-160°W)below 100 m.The warm advection is located at 0-20 m between 120°and 80°W,and it is also in line with the warm anomalies near the surface of the Pacific in Fig.5.In the Indian Ocean,the warm advection center is located at 60-120 m between 90°and 100°E,which could well explain the warm anomalies near the eastern boundary of theIndian Ocean in Fig.5.The cold advection below 120 m covers the partial region with cold anomalies in Fig.5.Considering that neither the total NDH nor the mean advection could provide an adequate explanation of the mean state change,the non-linearity of the vertical mixing(R in the temperature formula) should also play a role.The largest temperature differences mainly occurred around the thermocline,where vertical mixing cannot be ignored.Moreover,because any verticalmixing scheme is nonlinear and dependent on the mean state,it is vital to diagnose the contribution of vertical mixing to the heat budget in future studies.Although temperature change is related,to some extent,to vertical mixing,understanding the effects of NDH and mean advection on the mean state is essential.

    5.2 Heat budget in the mixed layer

    Fig.6.Differences in NDH between control run and(a,c)WSA and(b,d)WSA?R in equatorial upper Pacific and Indian Oceans(2°N-2°S).

    The pattern of mean temperature change in the mixed layer from WSA is very similar to that fromWSA?R(Fig.8),illustrating the importance of nonlinear dynamics.As can be seen,the temperature increases in the eastern Pacific(15°S-10°N,120°-80°W) and decreases in the western Pacific(140°-160°E).In the Indian Ocean,there is a warming center located over 0°-10°S,40°-80°E,and there are also some cooling regions such as north of equator,south of 10°S.

    Figure 9 further shows the corresponding differences in the NDH between the sensitivity runs and control run.The NDH contributes much to the cooling region in the western Pacific(140°-160°E)as well as to some warming regions(0°-3°N,120°-80°W).In the Indian Ocean,the NDH shows a negative value in most regions that cover the cooling regions in Fig.8.

    The mean advection differences are shown in Fig. 10.These can explain the temperature increases in most regions of the eastern Pacific(15°S-15°N,120°-80°W)and also the warming center(0°-10°S,40°-80°E)in the Indian Ocean.

    Fig.7.Differences in mean advection between control run and(a,c)WSA and(b,d)WSA R in equatorial upper Pacific and Indian Oceans(2°N-2°S).

    6.Conclusions

    The temperature anomalies associated with ENSO/IOD in the tropical Pacific and Indian Oceans not only show significant asymmetry,but also have a rectification effect on the mean state or decadal climate variations.We conducted several numerical experiments with an OGCM in this study.These experiments were designed(1)to evaluate the model’s ability to reproduce the interannual variations,especially the temperature asymmetry in tropical basins;(2)to help quantify the relative contributions of wind stress forcing and ocean dynamical processes to the asymmetry of temperature anomalies;and(3)to understand the rectification effect of the interannual variations on the mean state and its physical mechanism.

    Fig.8.As in Fig.5,but for mixed layer condition.

    When forced by the observed wind stress with interannual anomalies,the model well reproduced thelarge-scale patterns of the mean states in these two basins.Meanwhile,the generalspatialpattern of temperature asymmetry simulated by the ocean model was similar to the reanalysis in the upper 40-160 m in Pacific and 0-140 m in Indian Ocean.The asymmetry in such subsurface regions mainly results from the nonlinear response to the surface wind stress anomalies through internal dynamic processes,rather than the skewness of winds.The climatological mean states simulated by the two numerical experiments were alike,but different from that of the control run, which indicates that the two sensitivity runs,in which the wind stress anomalies were opposite to each other, exerted similar time-mean effects on the mean state. Further analysis showed that NDH plays a crucial role in the warming of the equatorial eastern Pacific below20 m and in the cooling ofthe upper equatorialeastern Indian Ocean.Mean advection contributes much to the cooling in the equatorial western Pacific below 100 m as well as to the warming at the eastern boundary of the Indian Ocean.

    Fig.9.As in Fig.6,but for most important components,i.e.,(a,b)zonal NDH for Pacific and(c,d)vertical NDH for Indian Ocean,in mixed layer condition.Units are℃month?1.

    There are some interannual variation modes in SST anomalies in the tropical Atlantic(Servain,1991; M′elice and Servain,2003),but previous studies have rarely addressed the temperature asymmetry in this basin.Besides numerical experiments in the Pacific and Indian oceans,we conducted similar numericalexperiments in the tropical Atlantic(30°N-30°S,80°W-20°E)in this study.Comparison of these simulations indicated that the asymmetric responses of the temperature and velocity anomalies to the wind stress anomaly are highly universal in the tropical oceans; temperature anomalies are less asymmetric in the tropical Atlantic than they are in the tropical Pacificand Indian oceans,and the corresponding rectification of the mean state is also much weaker.

    Fig.10.As in Fig.7,but for most important components,i.e.,vertical advection for(a,b)Pacific and(c,d)Indian Oceans,in mixed layer condition.Units are℃month?1.

    REFERENCES

    An,S.I.,and Jin Feifei,2001:Collective role of thermocline and zonal advective feedbacks in the ENSO mode.J.Climate,14,3421-3432.

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

    An,S.I.,Y.G.Ham,J.S.Kug,et al.,2005:El Ni?no-La Ni?na asymmetry in the coupled model intercomparison project simulations.J.Climate,18,2617-2627. Battisti,D.S.,1988:Dynamics and thermodynamics of a warming event in a coupled tropical atmosphereocean model.J.Atmos.Sci.,45,2889-2919.

    Burgers,G.,and D.B.Stephenson,1999:The“normality”of El Ni?no.Geophys.Res.Lett.,26, 1027-1030,doi:10.1029/1999GL900161.

    Cai Wenju,H.H.Hendon,and G.Meyers,2005:Indian Ocean dipolelike variability in the CSIRO Mark 3 coupled climate model.J.Climate,18,1449-1468.

    Cai Wenju and Qiu Yun,2013:An observation-based assessment of nonlinear feedback processes associated with the Indian Ocean dipole.J.Climate,26, 2880-2890.

    Cai Yi,Wang Zhanggui,Yu Zhouwen,et al.,2003:The EOF analysis of temperature and zonal flow in the equatorial Pacific Ocean and the study of the El Ni?no forecasting.Acta Oceanol.Sinica,25,12-18. (in Chinese)

    Cai Yi,Li Hai,and Zhang Renhe,2008:A study on the relationship between ENSO and tropical Indian Ocean temperature.Acta Meteor.Sinica,66,120-124.(in Chinese)

    Canuto,V.M.,and M.S.Dubovikov,2005:Modeling mesoscale eddies.Ocean Modeling,8,1-30.

    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.

    Choi,J.,S.I.An,and S.W.Yeh,2012:Decadal amplitude modulation of two types of ENSO and its relationship with the mean state.Climate Dyn.,38, 2631-2644.

    Choi,J.,S.I.An,S.W.Yeh,et al.,2013:ENSO-like and ENSO-induced tropical Pacific decadal variability in CGCMs.J.Climate,26,1485-1501.

    Fedorov,A.V.,and S.G.Philander,2001:A stability analysis of tropical ocean-atmosphere interactions: Bridging measurements and theory for El Ni?no.J. Climate,14,3086-3101.

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

    Hua Lijuan,Yu Yongqiang,and Yin Baoshu,2010:Numerical modeling of the asymmetry of Indian Ocean dipole and its mechanism.Chinese J.Atmos.Sci.,34,1046-1058.(in Chinese)

    Hua Lijuan,Yu Yongqiang,and Sun Dezheng,2015:A further study of ENSO rectification:Results from an OGCM with a seasonal cycle.J.Climate,28, 1362-1382.

    Hong,C.C.,T.Li,L.Ho,et al.,2008:Asymmetry of the Indian Ocean dipole.Part I:Observational analysis. J.Climate,21,4834-4848.

    Large,W.G.,G.Danabasoglu,S.C.Doney,et al.,1997: Sensitivity to surface forcing and boundary layer mixing in a Global Ocean Model:Annual-mean climatology.J.Phys.Oceanogr.,27,2418-2447.

    Levitus,S.,and T.P.Boyer,1994:World Ocean Atlas 1994 Volume 4:Temperature.NOAA Atlas NESDIS 4.US Department of Commerce,Washington D.C.,1-117.

    Levitus,S.,R.Burgett,and T.P.Boyer,1994:World Ocean Atlas 1994 Volume 3:Salinity.NOAA Atlas NESDIS 3.US Department of Commerce,Washington D.C.,1-99.

    Li Yangchun,Xu Yongfu,Chu Min,et al.,2012:Influences of climate change on the uptake and storage of anthropogenic CO2in the global ocean.Acta Meteor.Sinica,26,304-317,doi:10.1007/s13351-012-0304-z.

    Liang Jin,Yang Xiuqun,and Sun Dezheng,2012:The effect of ENSO events on the tropical Pacific mean climate:Insights from an analytical model.J.Climate,25,7590-7606.

    Liu Hailong,Zhang Xuehong,Yu Yongqiang,et al., 2004a:Manual for LASG/IAP Climate System Ocean Model.Science Press,Beijing,1-108.(in Chinese)

    Liu Hailong,Zhang Xuehong,Li Wei,et al.,2004b:An eddy-permitting oceanic general circulation model and its preliminary evaluation.Adv.Atmos.Sci.,21,675-690.

    Liu Hailong,Lin Pengfei,Yu Yongqiang,et al.,2012: The baseline evaluation of LASG/IAP Climate System Ocean Model(LICOM)Version 2.Acta Meteor. Sinica,26,318-329,doi:10.1007/s13351-012-0305-y.

    M′elice,J.L.,and J.Servain,2003:The tropical Atlantic meridional SST gradient index and its relationships with the SOI,NAO,and southern Ocean.Climate Dyn.,20,447-464.

    Ogata,T.,Xie Shangping,Lan Jian,et al.,2013:Importance of ocean dynamics for the skewness of the Indian Ocean dipole mode.J.Climate,26,2145-2159.

    Ohba,M.,and H.Ueda,2009:Role of nonlinear atmospheric response to SST on the asymmetric transition process of ENSO.J.Climate,22,177-192.

    Okumura,Y.M.,and C.Deser,2010:Asymmetry in the duration of El Ni?no and La Ni?na.J.Climate,23, 5826-5843.

    Rodgers,K.B.,P.Friederichs,and M.Latif,2004:Tropical Pacific decadal variability and its relation to decadalmodulations of ENSO.J.Climate,17,3761-3774.

    Roeske,F.,2001:An Atlas of Surface Flues Based on the ECMWF Reanalysis-A Climatological Dataset to Force Global Ocean General Circulation Models. Report No.23,Hamburg:Max-Planck-Institut f¨ur Meteorologie,1-31.

    Rong Xinyao,Zhang Renhe,T.Li,et al.,2011:Upscale feedback of high-frequency winds to ENSO.Quart. J.R.Meteor.Soc.,137,894-907.

    Sa ji,N.H.,B.N.Goswami,P.N.Vinayachandran,et al., 1999:A dipole mode in the tropical Indian Ocean. Nature,401,360-363.

    Sa ji,N.H.,and T.Yamagata,2003:Structure of SST and surface wind variability during Indian Ocean dipole mode events:COADS observations.J.Climate,16,2735-2751.

    Servain,J.,1991:Simple climatic indices for the tropical Atlantic Ocean and some applications.J.Geophys. Res.,96,15137-15146.

    Su Jingzhi,Zhang Renhe,T.Li,et al.,2010:Causes of the El Ni?no and La Ni?na amplitude asymmetry in the equatorial eastern Pacific.J.Climate,23, 605-617.

    Sun Dezheng,2003:A possible effect of an increase in the warm-pool SST on the magnitude of El Ni?no warming.J.Climate,16,185-205.

    Sun Dezheng and Zhang Tao,2006:A regulatory effect of ENSO on the time-mean thermal stratification of the equatorial upper ocean.Geophys.Res.Lett.,33,L07710.

    Sun Dezheng,Zhang Tao,Sun Yan,et al.,2014:Rectification of El Ni?no-Southern Oscillation into climate anomalies of decadal and longer time scales:Results from forced ocean GCM experiments.J.Climate,27,2545-2561.

    Sun Yan,Sun Dezheng,Wu Lixin,et al.,2012:The western Pacific warm pool and ENSO asymmetry in CMIP3 models.Adv.Atmos.Sci.,30,940-953.

    Uppala,S.M.,P.W.K°allberg,A.J.Simmons,et al., 2005:The ERA-40 re-analysis.Quart.J.Roy.Meteor.Soc.,131,2961-3012.

    Yan Bangliang and Zhang Renhe,2002:The role of atmosphere climate basic state in the formation of westerly over the equatorial Pacific.Acta Oceanol. Sin.,24,39-50.(in Chinese)

    Yu Jinyi and S.T.Kim,2011:Reversed spatial asymmetries between El Ni?no and La Ni?na and their linkage to decadal ENSO modulation in CMIP3 models.J. Climate,24,5423-5434.

    Yu Yongqiang,Zheng Weipeng,Wang Bin,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.

    Zebiak,S.E.,and M.A.Cane,1987:A model El Ni?no-Southern oscillation.Mon.Wea.Rev.,115,2262-2278.

    Zhang Tao,Sun Dezheng,R.Neale,et al.,2009:An evaluation of ENSO asymmetry in the community climate system models:A view from the subsurface. J.Climate,22,5933-5961.

    Zheng Xiaotong,Xie Shangping,G.A.Vecchi,et al., 2010:Indian Ocean dipole response to global warming:Analysis of ocean-atmospheric feedbacks in a coupled model.J.Climate,23,1240-1253.

    Zhong Aihong,H.H.Hendon,and O.Alves,2005:Indian Ocean variability and its association with ENSO in a global coupled model.J.Climate,18,3634-3649.

    Hua Lijuan and Yu Yongqiang,2015:Nonlinear responses of oceanic temperature to wind stress anomalies in the tropical Pacific and Indian Oceans:A study based on numerical experiments with an OGCM.J.Meteor.Res.,29(4),608-626,

    10.1007/s13351-015-4115-x.

    Supported by the“Strategic Priority Research Program Climate Change:Carbon Budget and Relevant Issues”of the Chinese Academy of Sciences(XDA05110302),National(Key)Basic Research and Development(973)Program of China(2013CB956204), and Jiangsu Collaborative Innovation Center for Climate Change.

    ?yyq@lasg.iap.ac.cn.

    ?The Chinese Meteorological Society and Springer-Verlag Berlin Heidelberg 2015

    November 6,2014;in final form June 8,2015)

    国产成人av教育| 97碰自拍视频| 久久久水蜜桃国产精品网| a级毛片在线看网站| 精品国产一区二区久久| 免费观看精品视频网站| 精品一区二区三区视频在线观看免费| 日韩精品免费视频一区二区三区| 亚洲精品中文字幕一二三四区| 久久久精品国产亚洲av高清涩受| 热re99久久国产66热| 国产一区二区三区视频了| 少妇被粗大的猛进出69影院| 最近最新免费中文字幕在线| 国产片内射在线| 少妇熟女aⅴ在线视频| 亚洲全国av大片| 国产精品二区激情视频| 波多野结衣一区麻豆| 亚洲天堂国产精品一区在线| 一二三四在线观看免费中文在| 一级a爱片免费观看的视频| 亚洲中文av在线| 久久午夜亚洲精品久久| 黄片播放在线免费| 麻豆久久精品国产亚洲av| 亚洲视频免费观看视频| 禁无遮挡网站| 看免费av毛片| 亚洲全国av大片| 亚洲精品国产精品久久久不卡| 久久性视频一级片| 精品国内亚洲2022精品成人| 999久久久国产精品视频| 天堂动漫精品| 欧美日韩福利视频一区二区| 亚洲色图综合在线观看| 韩国精品一区二区三区| 久久中文字幕人妻熟女| 大码成人一级视频| 免费无遮挡裸体视频| 丁香欧美五月| 欧美在线一区亚洲| 色尼玛亚洲综合影院| 一边摸一边抽搐一进一出视频| 在线视频色国产色| 色播亚洲综合网| 自拍欧美九色日韩亚洲蝌蚪91| 香蕉久久夜色| 99国产极品粉嫩在线观看| 女性生殖器流出的白浆| 两人在一起打扑克的视频| 丝袜美足系列| 亚洲国产精品sss在线观看| 国内毛片毛片毛片毛片毛片| 久久久久国产一级毛片高清牌| 欧美乱码精品一区二区三区| 欧美在线黄色| 中国美女看黄片| 午夜久久久久精精品| 亚洲欧美精品综合一区二区三区| 日韩欧美三级三区| 精品熟女少妇八av免费久了| 国产亚洲精品av在线| 国产精品久久久久久亚洲av鲁大| 国产一区二区三区在线臀色熟女| 精品欧美一区二区三区在线| 超碰成人久久| 少妇被粗大的猛进出69影院| 一区在线观看完整版| 欧美一级a爱片免费观看看 | 国产熟女xx| 久久欧美精品欧美久久欧美| 黄频高清免费视频| 法律面前人人平等表现在哪些方面| 亚洲午夜精品一区,二区,三区| 亚洲第一av免费看| 麻豆av在线久日| 久久久久精品国产欧美久久久| 亚洲专区国产一区二区| 久久久水蜜桃国产精品网| 操出白浆在线播放| 不卡av一区二区三区| 日本欧美视频一区| 久久精品91无色码中文字幕| 看免费av毛片| 嫩草影视91久久| 国产精品爽爽va在线观看网站 | 色综合欧美亚洲国产小说| 在线十欧美十亚洲十日本专区| 久久国产精品男人的天堂亚洲| 两个人视频免费观看高清| 自拍欧美九色日韩亚洲蝌蚪91| 免费无遮挡裸体视频| 丰满的人妻完整版| 操出白浆在线播放| 757午夜福利合集在线观看| 国产真人三级小视频在线观看| 高清黄色对白视频在线免费看| 美女扒开内裤让男人捅视频| 亚洲精品在线美女| or卡值多少钱| 熟女少妇亚洲综合色aaa.| 亚洲专区国产一区二区| 女生性感内裤真人,穿戴方法视频| 国产国语露脸激情在线看| 精品久久久精品久久久| 午夜福利一区二区在线看| 久久精品成人免费网站| 欧美一级a爱片免费观看看 | 99久久久亚洲精品蜜臀av| 午夜激情av网站| 成人三级黄色视频| a在线观看视频网站| 一卡2卡三卡四卡精品乱码亚洲| 91大片在线观看| 十八禁人妻一区二区| 99国产精品一区二区蜜桃av| 欧美成人一区二区免费高清观看 | 中国美女看黄片| 国产麻豆69| 亚洲天堂国产精品一区在线| 国产亚洲精品综合一区在线观看 | 国产精品亚洲av一区麻豆| 人人澡人人妻人| 亚洲va日本ⅴa欧美va伊人久久| 97人妻天天添夜夜摸| 一区二区三区激情视频| 午夜a级毛片| or卡值多少钱| 十八禁人妻一区二区| 亚洲精品国产一区二区精华液| 亚洲精品一卡2卡三卡4卡5卡| 深夜精品福利| 亚洲免费av在线视频| 少妇熟女aⅴ在线视频| 精品国产美女av久久久久小说| 国产精华一区二区三区| 91在线观看av| 久久中文看片网| 欧美色视频一区免费| 99国产精品99久久久久| 婷婷六月久久综合丁香| 欧美人与性动交α欧美精品济南到| 亚洲av片天天在线观看| 两个人免费观看高清视频| 久久久久九九精品影院| 9191精品国产免费久久| 午夜福利欧美成人| 久9热在线精品视频| 久久精品国产综合久久久| 亚洲国产精品999在线| 老鸭窝网址在线观看| 人人澡人人妻人| 久久欧美精品欧美久久欧美| 在线观看午夜福利视频| 波多野结衣巨乳人妻| 婷婷精品国产亚洲av在线| 亚洲精品中文字幕一二三四区| 级片在线观看| 日本五十路高清| 国产成人精品无人区| 高清在线国产一区| 欧美日本亚洲视频在线播放| 国产色视频综合| 国产亚洲精品综合一区在线观看 | 男女午夜视频在线观看| 国产亚洲av嫩草精品影院| 看片在线看免费视频| 亚洲成人精品中文字幕电影| 青草久久国产| 搡老熟女国产l中国老女人| 精品国产一区二区三区四区第35| 69精品国产乱码久久久| 亚洲一区中文字幕在线| 久久香蕉国产精品| 高清在线国产一区| 久久久国产欧美日韩av| 丝袜人妻中文字幕| 精品高清国产在线一区| 精品国产美女av久久久久小说| 可以在线观看的亚洲视频| 欧美乱色亚洲激情| 一二三四社区在线视频社区8| 亚洲人成伊人成综合网2020| 亚洲人成网站在线播放欧美日韩| 久久人妻福利社区极品人妻图片| 亚洲在线自拍视频| 免费少妇av软件| 亚洲五月婷婷丁香| 老司机深夜福利视频在线观看| 国产精品,欧美在线| 国产麻豆成人av免费视频| 欧美国产精品va在线观看不卡| 亚洲成人国产一区在线观看| 色综合欧美亚洲国产小说| 999久久久精品免费观看国产| 久久中文字幕人妻熟女| 精品国内亚洲2022精品成人| 丰满的人妻完整版| www国产在线视频色| 久久香蕉精品热| 欧美日韩中文字幕国产精品一区二区三区 | 亚洲精品av麻豆狂野| 美女大奶头视频| 嫁个100分男人电影在线观看| 亚洲中文日韩欧美视频| 国产主播在线观看一区二区| x7x7x7水蜜桃| 欧美精品亚洲一区二区| 午夜亚洲福利在线播放| 中文字幕高清在线视频| 精品乱码久久久久久99久播| 亚洲欧美激情综合另类| 黄色 视频免费看| 国产伦人伦偷精品视频| 丁香六月欧美| 国产成年人精品一区二区| 首页视频小说图片口味搜索| 老司机靠b影院| 国产成人精品无人区| 国产激情欧美一区二区| 男人舔女人下体高潮全视频| 国产区一区二久久| 亚洲精品在线美女| 母亲3免费完整高清在线观看| 亚洲精品av麻豆狂野| 久久亚洲精品不卡| 欧美中文日本在线观看视频| 成人av一区二区三区在线看| 一本久久中文字幕| 色老头精品视频在线观看| 国产精品日韩av在线免费观看 | 国产欧美日韩精品亚洲av| 欧美中文日本在线观看视频| 12—13女人毛片做爰片一| av视频在线观看入口| 国内精品久久久久精免费| 中文字幕人成人乱码亚洲影| 极品人妻少妇av视频| 亚洲国产高清在线一区二区三 | 高清毛片免费观看视频网站| 欧美日韩精品网址| 91av网站免费观看| 国产亚洲精品综合一区在线观看 | 午夜a级毛片| 久久这里只有精品19| 精品少妇一区二区三区视频日本电影| 18禁美女被吸乳视频| 精品国产乱子伦一区二区三区| 男女下面进入的视频免费午夜 | 日韩欧美免费精品| 国产成人免费无遮挡视频| 一区福利在线观看| 亚洲成国产人片在线观看| 老汉色∧v一级毛片| 两个人视频免费观看高清| 成人精品一区二区免费| 国产精品二区激情视频| 91大片在线观看| 国产又爽黄色视频| 国产成人av教育| 夜夜看夜夜爽夜夜摸| 国产伦人伦偷精品视频| 好看av亚洲va欧美ⅴa在| 久久精品影院6| 成人18禁在线播放| 亚洲精品久久成人aⅴ小说| 国产麻豆69| 欧美日韩黄片免| 国产成人精品久久二区二区91| 在线国产一区二区在线| 日韩欧美三级三区| 日韩欧美国产在线观看| 黄色女人牲交| 精品欧美一区二区三区在线| 天天躁夜夜躁狠狠躁躁| 99精品久久久久人妻精品| 国语自产精品视频在线第100页| 精品久久久久久久毛片微露脸| 我的亚洲天堂| 美女高潮喷水抽搐中文字幕| 国产精品永久免费网站| 国产主播在线观看一区二区| 亚洲av电影不卡..在线观看| 欧美av亚洲av综合av国产av| 亚洲一区高清亚洲精品| 此物有八面人人有两片| 亚洲人成电影免费在线| 一区二区三区高清视频在线| 久久九九热精品免费| 国产精品国产高清国产av| 亚洲av电影在线进入| 国产在线精品亚洲第一网站| 亚洲成人国产一区在线观看| 日韩欧美国产在线观看| 精品国产美女av久久久久小说| 中出人妻视频一区二区| 亚洲国产精品999在线| 国产麻豆69| 女人被躁到高潮嗷嗷叫费观| 久久午夜亚洲精品久久| 妹子高潮喷水视频| 桃红色精品国产亚洲av| 精品少妇一区二区三区视频日本电影| 中亚洲国语对白在线视频| 亚洲最大成人中文| 两个人视频免费观看高清| 国产熟女xx| 国产精品日韩av在线免费观看 | 国产亚洲精品久久久久5区| 日日爽夜夜爽网站| 国产一区二区激情短视频| 久久午夜综合久久蜜桃| 国产蜜桃级精品一区二区三区| 首页视频小说图片口味搜索| 波多野结衣av一区二区av| 国产精品美女特级片免费视频播放器 | 狠狠狠狠99中文字幕| 丝袜美足系列| 亚洲狠狠婷婷综合久久图片| 热re99久久国产66热| 人人妻人人爽人人添夜夜欢视频| 波多野结衣一区麻豆| 亚洲在线自拍视频| 国产99白浆流出| 两个人免费观看高清视频| 久久久久久久久久久久大奶| 亚洲成国产人片在线观看| 亚洲专区字幕在线| 人人妻人人澡欧美一区二区 | 熟妇人妻久久中文字幕3abv| 一区二区日韩欧美中文字幕| 电影成人av| 脱女人内裤的视频| 久久精品91蜜桃| 久久久精品欧美日韩精品| 国产97色在线日韩免费| 久久精品亚洲精品国产色婷小说| 精品久久久久久,| 最近最新免费中文字幕在线| 99精品在免费线老司机午夜| 久久久国产欧美日韩av| 在线观看午夜福利视频| 精品久久久久久久人妻蜜臀av | 夜夜爽天天搞| 脱女人内裤的视频| 嫩草影视91久久| 丝袜美腿诱惑在线| 国产亚洲精品一区二区www| 欧美日韩中文字幕国产精品一区二区三区 | 精品国产超薄肉色丝袜足j| 午夜福利在线观看吧| 欧美精品亚洲一区二区| 在线观看午夜福利视频| 亚洲 欧美 日韩 在线 免费| 男人操女人黄网站| 12—13女人毛片做爰片一| 麻豆国产av国片精品| 亚洲 欧美 日韩 在线 免费| 国产欧美日韩一区二区精品| 亚洲专区中文字幕在线| 热99re8久久精品国产| 18禁观看日本| 无人区码免费观看不卡| 亚洲av成人不卡在线观看播放网| 国产麻豆69| 亚洲伊人色综图| 日韩一卡2卡3卡4卡2021年| 天天躁夜夜躁狠狠躁躁| 日本撒尿小便嘘嘘汇集6| 黑人巨大精品欧美一区二区蜜桃| 国内精品久久久久久久电影| 人人妻人人澡欧美一区二区 | 少妇被粗大的猛进出69影院| 久久国产精品人妻蜜桃| 99精品欧美一区二区三区四区| 亚洲在线自拍视频| 国产一区二区激情短视频| 久久 成人 亚洲| 成人亚洲精品一区在线观看| 国产伦一二天堂av在线观看| 亚洲一码二码三码区别大吗| 国产精品免费视频内射| 国产精华一区二区三区| 久久久久久人人人人人| 精品国产一区二区久久| 无遮挡黄片免费观看| 久久久久久亚洲精品国产蜜桃av| 国产成人系列免费观看| √禁漫天堂资源中文www| 在线永久观看黄色视频| 欧美一级a爱片免费观看看 | www.999成人在线观看| 国产成人av激情在线播放| 国产欧美日韩精品亚洲av| 可以在线观看毛片的网站| АⅤ资源中文在线天堂| 波多野结衣巨乳人妻| 国产精品自产拍在线观看55亚洲| 国产成+人综合+亚洲专区| 欧美日韩亚洲综合一区二区三区_| 又黄又粗又硬又大视频| 国产成人精品久久二区二区91| 国产乱人伦免费视频| 老鸭窝网址在线观看| 麻豆成人av在线观看| 一边摸一边抽搐一进一小说| 精品国产超薄肉色丝袜足j| 久9热在线精品视频| 精品电影一区二区在线| 99久久精品国产亚洲精品| 最好的美女福利视频网| 亚洲男人的天堂狠狠| 久热这里只有精品99| 国产aⅴ精品一区二区三区波| 欧美激情 高清一区二区三区| 久久国产精品人妻蜜桃| av有码第一页| 搞女人的毛片| 久久精品成人免费网站| 少妇粗大呻吟视频| 国内精品久久久久久久电影| 久久国产精品男人的天堂亚洲| 日韩欧美在线二视频| 国产精品,欧美在线| 中文字幕久久专区| 99国产综合亚洲精品| av片东京热男人的天堂| 国产成人系列免费观看| 国产av一区在线观看免费| 88av欧美| 波多野结衣巨乳人妻| 亚洲无线在线观看| 中文字幕精品免费在线观看视频| 国产精品久久久久久人妻精品电影| 久久精品亚洲熟妇少妇任你| 大型av网站在线播放| 亚洲欧美精品综合一区二区三区| 男人操女人黄网站| 日韩欧美国产一区二区入口| 人妻久久中文字幕网| 一区二区日韩欧美中文字幕| 精品国产乱子伦一区二区三区| 久久精品人人爽人人爽视色| 欧美国产精品va在线观看不卡| 中亚洲国语对白在线视频| 国产精品一区二区精品视频观看| 99riav亚洲国产免费| 在线观看舔阴道视频| e午夜精品久久久久久久| 少妇 在线观看| 亚洲狠狠婷婷综合久久图片| √禁漫天堂资源中文www| 天天躁夜夜躁狠狠躁躁| 热99re8久久精品国产| 亚洲国产高清在线一区二区三 | 一进一出好大好爽视频| 一卡2卡三卡四卡精品乱码亚洲| 亚洲伊人色综图| 久久午夜亚洲精品久久| 国产精品亚洲一级av第二区| 男女之事视频高清在线观看| 香蕉久久夜色| 男女之事视频高清在线观看| 天天一区二区日本电影三级 | 亚洲中文日韩欧美视频| 长腿黑丝高跟| 国产aⅴ精品一区二区三区波| 两个人免费观看高清视频| 成人三级做爰电影| 中文字幕人成人乱码亚洲影| 女人被狂操c到高潮| 人人妻,人人澡人人爽秒播| 99精品久久久久人妻精品| 午夜视频精品福利| 在线观看免费午夜福利视频| 久久精品91蜜桃| 免费人成视频x8x8入口观看| 丁香六月欧美| 成人免费观看视频高清| 亚洲午夜理论影院| 黄色视频,在线免费观看| 亚洲国产欧美网| 99国产综合亚洲精品| 欧美日韩瑟瑟在线播放| 丝袜美腿诱惑在线| 搞女人的毛片| 老鸭窝网址在线观看| 侵犯人妻中文字幕一二三四区| 男人的好看免费观看在线视频 | 性欧美人与动物交配| 美女大奶头视频| 成人亚洲精品av一区二区| 精品一品国产午夜福利视频| 午夜久久久久精精品| 亚洲av成人av| 精品国产乱子伦一区二区三区| 1024香蕉在线观看| 黄色a级毛片大全视频| av欧美777| 老司机午夜福利在线观看视频| 老熟妇仑乱视频hdxx| a在线观看视频网站| 777久久人妻少妇嫩草av网站| 午夜亚洲福利在线播放| 久久久久久久久中文| 露出奶头的视频| 国产99久久九九免费精品| 欧美午夜高清在线| 国产精品永久免费网站| 日韩欧美一区视频在线观看| 久久草成人影院| 久久狼人影院| 国产单亲对白刺激| 亚洲少妇的诱惑av| 在线观看免费午夜福利视频| 日本免费一区二区三区高清不卡 | 午夜视频精品福利| 在线观看日韩欧美| 99久久久亚洲精品蜜臀av| 波多野结衣一区麻豆| 成人特级黄色片久久久久久久| 99riav亚洲国产免费| 在线观看免费视频日本深夜| 最近最新免费中文字幕在线| 9热在线视频观看99| 亚洲精华国产精华精| 亚洲自偷自拍图片 自拍| 欧美成人性av电影在线观看| 亚洲色图av天堂| 亚洲人成电影观看| 9热在线视频观看99| 这个男人来自地球电影免费观看| 欧美激情久久久久久爽电影 | 亚洲国产精品合色在线| 51午夜福利影视在线观看| 亚洲av成人一区二区三| 久久久久久国产a免费观看| 国产99久久九九免费精品| 91精品三级在线观看| 国产1区2区3区精品| 亚洲精品国产精品久久久不卡| 麻豆一二三区av精品| 国产成年人精品一区二区| 中文字幕av电影在线播放| 亚洲中文av在线| 真人做人爱边吃奶动态| videosex国产| 亚洲av片天天在线观看| 国产精品爽爽va在线观看网站 | 亚洲av电影不卡..在线观看| 日韩成人在线观看一区二区三区| 成人亚洲精品av一区二区| 黄色女人牲交| 久久国产精品人妻蜜桃| 黄片小视频在线播放| 怎么达到女性高潮| 十八禁网站免费在线| 九色国产91popny在线| 19禁男女啪啪无遮挡网站| 亚洲 欧美一区二区三区| 嫩草影院精品99| 亚洲av片天天在线观看| 丝袜美腿诱惑在线| 九色亚洲精品在线播放| 嫁个100分男人电影在线观看| 精品国产亚洲在线| 国产欧美日韩综合在线一区二区| 好男人电影高清在线观看| 精品国产一区二区久久| 国产精品久久久久久亚洲av鲁大| 亚洲人成伊人成综合网2020| 999久久久国产精品视频| 性少妇av在线| 最好的美女福利视频网| 一卡2卡三卡四卡精品乱码亚洲| 91麻豆精品激情在线观看国产| 精品第一国产精品| 99精品欧美一区二区三区四区| 麻豆一二三区av精品| 亚洲性夜色夜夜综合| 无限看片的www在线观看| 狠狠狠狠99中文字幕| 亚洲欧美日韩另类电影网站| 国产精华一区二区三区| 神马国产精品三级电影在线观看 | 校园春色视频在线观看| 淫秽高清视频在线观看| 亚洲精品一卡2卡三卡4卡5卡| 久久国产精品男人的天堂亚洲| 免费一级毛片在线播放高清视频 | 日韩欧美国产一区二区入口| 高清在线国产一区| 欧美日韩乱码在线| 亚洲专区中文字幕在线| 久久影院123| 黄频高清免费视频| 亚洲 国产 在线| xxx96com| 啪啪无遮挡十八禁网站| 午夜影院日韩av| 成人欧美大片| 高清毛片免费观看视频网站| 波多野结衣av一区二区av| 日日摸夜夜添夜夜添小说| 18禁裸乳无遮挡免费网站照片 | 99热只有精品国产| 成人国产综合亚洲| 99在线视频只有这里精品首页| 亚洲国产精品999在线| 精品第一国产精品| 大陆偷拍与自拍| 国语自产精品视频在线第100页| 欧美乱妇无乱码|