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

    The Influence of Arctic Sea Ice Concentration Perturbations on Subseasonal Predictions of North Atlantic Oscillation Events※

    2023-12-26 09:12:16GuokunDAIMuMUZheHANChunxiangLIZhinaJIANGMengbinZHUandXueyingMA
    Advances in Atmospheric Sciences 2023年12期

    Guokun DAI, Mu MU*, Zhe HAN, Chunxiang LI, Zhina JIANG,Mengbin ZHU, and Xueying MA

    1Department of Atmospheric and Oceanic Sciences and Institute of Atmospheric Sciences,Fudan University, Shanghai 200438, China

    2CMA-FDU Joint Laboratory of Marine Meteorology, Shanghai 200438, China

    3Shanghai Frontiers Science Center of Atmosphere-Ocean Interaction, Shanghai 200438, China

    4Innovation Center of Ocean and Atmosphere System, Zhuhai Fudan Innovation Research Institute, Zhuhai 509031, China

    5CAS Key Laboratory of Regional Climate-Environment for Temperate East Asia, Institute of Atmospheric Physics,Chinese Academy of Sciences, Beijing 100029, China

    6State Key Laboratory of Severe Weather (LaSW), Chinese Academy of Meteorological Sciences, Beijing 100081, China

    7Beijing Institute of Applied Meteorology, Beijing 100029, China

    8Yantai Vocational College, Yantai 264670, China

    ABSTRACT The influence of Arctic sea ice concentration (SIC) on the subseasonal prediction of the North Atlantic Oscillation(NAO) event is investigated by utilizing the Community Atmospheric Model version 4.The optimal Arctic SIC perturbations which exert the greatest influence on the onset of an NAO event from a lead of three pentads (15 days) are obtained with a conditional nonlinear optimal perturbation approach.Numerical results show that there are two types of optimal Arctic SIC perturbations for each NAO event, with one weakening event (marked as type-1) and another strengthening event (marked as type-2).For positive NAO events, type-1 optimal SIC perturbations mainly show positive SIC anomalies in the Greenland, Barents, and Okhotsk Seas, while type-2 perturbations mainly feature negative SIC anomalies in these regions.For negative NAO events, the optimal SIC perturbations have almost opposite patterns to those in positive events, although there are some differences among these SIC perturbations due to different atmospheric initial conditions.Further diagnosis reveals that the optimal Arctic SIC perturbations first modify the surface turbulent heat flux and the temperature in the lower troposphere via diabatic processes.Afterward, the temperature in the low troposphere is mainly affected by dynamic advection.Finally, potential vorticity advection plays a crucial role in the 500-hPa geopotential height prediction in the northern North Atlantic sector during pentad 4, which influences NAO event prediction.These results highlight the importance of Arctic SIC on NAO event prediction and the spatial characteristics of the SIC perturbations may provide scientific support for target observations of SIC in improving NAO subseasonal predictions.

    Key words: optimal Arctic SIC perturbation, NAO event, subseasonal prediction, CNOP approach

    1.Introduction

    The North Atlantic Oscillation is the seesaw-like phenomenon in sea-level pressure (SLP) between the Icelandic Low and Azores High in the North Atlantic sector (Walker and Bliss, 1932).It is one of the most dominant low-frequency atmospheric variabilities in the mid-to-high latitudes of the Northern Hemisphere during boreal winter.The NAO not only influences weather in the region but also greatly impacts the global climate (Hurrell, 1995; Diao et al., 2015).Both the Icelandic Low and the Azores High weaken during the negative phase of the NAO (NAO–), causing cold in Europe and North America, while both atmospheric centers strengthen and induce warm winters there during its positive phase (NAO+) (Thompson et al., 2000; Hurrell and Deser,2010).

    The NAO has trended toward its negative phase during past decades, which is considered to be related to Arctic sea ice reduction (Francis et al., 2009; Vihma, 2014; Cohen et al., 2020).It has been reported that the Arctic sea ice concentration (SIC) reduction in autumn could modulate the planetary waves of the Barents-Kara Seas (BKS) and favor an NAO– response in winter (Jaiser et al., 2012, 2013).However, Arctic SIC anomalies in different sectors could trigger various NAO responses (Screen, 2017; McKenna et al.,2018).By diagnosing the reanalysis data, it is found that the NAO tends to favor its negative phase in its response to SIC reduction in the BKS, while the positive phase often accompanies SIC decline in the Canadian Archipelago (Luo et al.,2016).Moreover, SIC variability in the Pacific can lead to various NAO responses by triggering Rossby waves there(Yamamoto et al., 2006; Mesquita et al., 2011).

    However, the intrinsic timescale is approximately 20 days for the NAO life cycle (Feldstein, 2003).On this timescale, NAO events are usually related to extreme weather events in Eurasia.For example, the Middle East was affected by an extreme heavy snow event in December 2013, whose primary cause was believed to be a decaying NAO+ event (Luo et al., 2015).The essence of such a relationship inspires us to investigate NAO events on subseasonal timescales.

    An NAO event can be treated as a nonlinear initial value problem on the synoptic timescale and the initial condition plays the dominant role in NAO evolution (Benedict et al., 2004; Franzke et al., 2004; Luo et al., 2007, 2008).Perturbations in both the troposphere and even the stratosphere are crucial for NAO evolution (Luo et al., 2008; Nie et al.,2019).With the help of a quasi-geostrophic model, Jiang et al.(2013) investigated the optimal precursors (OPR)which trigger the onset of NAO events 5–8 days in advance.Moreover, they investigated the optimally growing initial errors (OGIE) in NAO onset predictions and discussed the relationship between the OPR and OGIE (Dai et al., 2016).After that, they further discussed the target observations to improve the prediction skills of NAO event onsets (Dai et al., 2019).

    In addition to synoptic timescales, much attention has been dedicated to subseasonal forecasts over the last decade(Vitart et al., 2017).Their work reveals that skillful forecasts are available with lead times longer than two weeks (Buizza and Leutbecher, 2015).For example, with a lead time of 2–4 weeks, the ECMWF can skillfully forecast the extreme cold frequency in boreal winter (Xiang et al., 2019; 2020).For certain strong and long-lasting cold events in East Asia,skillful forecasts could be made two weeks in advance or earlier (Dai and Mu, 2020; Dai et al., 2021a).Apart from synoptic timescales, both the atmospheric initial condition and boundary conditions e.g., sea surface temperature (SST) and SIC, are crucial for atmospheric evolution on subseasonal timescales (Liu et al., 2017; Wang et al., 2018), especially for large-scale atmospheric circulations like the NAO (Mariotti et al., 2020).

    It is found that sea ice anomalies could trigger a baroclinic atmospheric response with surface heat flux modification and such atmospheric responses can reach an equivalent barotropic structure in two months (Deser et al., 2007).Moreover, utilizing the Integrated Forecasting System from ECMWF, Semmler et al.(2016) investigated the influence of sudden Arctic sea ice thinning on atmospheric circulations by introducing a 10°C warming in sea ice surface temperature.They found that the atmosphere has a baroclinic response in the first several days due to the boundary layer process and reaches a quasi-equilibrium state in two months.After that, Dai et al.(2021b) discussed the impact of sudden Arctic sea ice thinning on NAO events.Their results reveal that NAO events have different responses to sudden Arctic sea ice thinning due to various NAO circulation patterns.The above results show that Arctic sea ice is an important boundary condition for atmospheric evolution on subseasonal timescales.

    However, large uncertainties exist in Arctic sea ice datasets, which greatly impact subseasonal atmospheric circulations (Dai and Mu, 2020).Rather than integrating the numerical model with various boundary conditions, the conditional nonlinear optimal perturbation (CNOP) approach is efficient in revealing the largest influence of Arctic SIC on an NAO event.The CNOP approach was first proposed to solve for the initial uncertainty in ENSO prediction (Mu et al., 2003) and was further developed to solve for the uncertainty in model parameters (Mu et al., 2010), model tendencies (Duan and Zhou, 2013), and boundary conditions(Wang and Mu, 2015).The CNOP approach is widely utilized in investigating atmospheric and oceanic problems, as summarized in Wang et al.(2020).Recently, Ma et al.(2022) investigated the role of Arctic sea ice anomalies in the subseasonal forecasts of Ural blocking events with CNOP approaches.Their results indicate that SIC in the Greenland, Barents,and Okhotsk seas is crucial for the prediction of Ural blocking on subseasonal timescales.

    Due to the relationships between NAO events and Eurasian extreme weather, we attempt to investigate the following questions.(1) What kind of Arctic sea ice perturbations could exert the greatest influence on NAO event predictions on subseasonal timescales? (2) What is the physical mechanism by which Arctic sea ice perturbations affect NAO event prediction? Understanding their solutions may help us conceptualize the role of Arctic SIC in NAO formation.Additionally, these solutions can provide scientific support in the Arctic SIC target observations and the ensemble perturbation generations, which can further improve the NAO event predictions.

    The remainder of this paper is organized as follows.In section 2, the data and methods are described in detail.In section 3, the response time of NAO events to Arctic SIC perturbations is determined by sensitivity numerical experiments.The optimal Arctic SIC perturbations for NAO+ and NAO–events are explored in sections 4 and 5, respectively.Finally,the conclusion and discussion are presented in section 6.

    2.Data and methods

    2.1.Data

    Daily atmospheric and oceanic data with a 1.0° × 1.0° resolution from the ERA-Interim reanalysis are used in this investigation (Dee et al., 2011).The atmospheric variables include sea level pressure (SLP), surface air temperature,300-hPa geopotential height (Z300) and 500-hPa geopotential height (Z500), and the 300-hPa and 500-hPa zonal winds in winter (December, January, and February), spanning the period from 1981/82 to 2000/01.The daily oceanic variables include SST and SIC spanning January 1982 to December 2001.The SIC data are derived from the Special Sensor Microwave Imager Sounder based on a smooth combination of Bootstrap and Bristol algorithms (Comiso, 1986; Dee et al., 2011).These daily oceanic variables are used as the boundary conditions for the numerical model, noting that the period 1982–2001 is consistent with the default forcing settings in the numerical model.

    2.2.The numerical model

    The Community Atmospheric Model, version 4(CAM4) with prescribed SST and SIC is utilized in this investigation (Neale et al., 2013).The atmospheric model has 26 vertical levels with a top at approximately 3.5 hPa, and the horizontal resolution is approximately 0.9° × 1.25°.The model settings are the same as those in Dai et al.(2021b).Briefly, the calendar means of the ERA-Interim daily SST and SIC ranging from January 1982 to December 2001 are used as the prescribed model boundary condition, while the CO2concentration is set to 367 ppmv.These are consistent with the other default settings in CAM4.Moreover, the sea ice thickness is uniformly set to 2 meters in the Arctic and 1 meter in the Antarctic.The numerical model has been run for 31 years with the first year of model data being discarded as the model spin-up period.The remaining 30 years of model data are adopted for analysis.Compared with the reanalysis data, CAM4 has a good representation of the atmospheric circulation during wintertime, which has been evaluated in Dai et al.(2021b).

    2.3.The NAO index and NAO events

    Similar to that in Dai et al.(2021b), the NAO index(NAOI) is used to describe the NAO phase and intensity.The NAOI is defined as the scaled projection of the daily SLP anomaly on the typical positive NAO anomaly mode.The formula can be expressed as:

    where S LPNAOrepresents the typical positive NAO mode in the SLP field, which is derived from the first empirical orthogonal function of the daily SLP anomaly in the North Atlantic sector (approximately 20°–80°N, 90°W–40°E).SLPdis the daily SLP anomaly relative to the climatology.The angle brackets denote the Euclidean inner product.

    Similar to the NAO event definition in Dai et al.(2021b), an NAO+ (NAO–) event is identified when the daily NAOI is larger than 1.0 (smaller than –1.0) for five or more consecutive days.The first day with NAOI larger than 1.0 (smaller than –1.0) is defined as the NAO+ (NAO–)onset day (marked as lag 0 day).Similarly, 15, 10, and 5 days before NAO onset are marked as lag–15, lag–10, and lag–5 days, respectively, while the day corresponding to 5 days after onset is marked as lag+5 day.According to the above criteria, 22 NAO+ and 22 NAO– events are recognized from the 30-year CAM4 simulation, showing a similar frequency to the 13 NAO+ and 14 NAO– events in the 20-year ERA-Interim reanalysis.The NAO events in the CAM4 simulation are not expected to coincide with those in the reanalysis.However, the compositions show that CAM4 describes the evolution of such NAO events well and thus can be utilized for further investigation [see Figs.4 and 5 in (Dai et al.,2021b)].

    2.4.CNOP method

    To investigate the influence of Arctic SIC on NAO event predictions on subseasonal timescales, the CNOP for boundary uncertainty (CNOP-B, Wang and Mu, 2015) is utilized since the Arctic SIC is treated as the boundary condition in the CAM4 model.In this investigation, CNOP-B represents the specific Arctic SIC perturbation that causes the largest uncertainty of the NAO evolution under a given physical constraint.Since the NAOI can measure the phase and intensity of an NAO event, the corresponding objective function is defined as the NAOI difference induced by the Arctic SIC perturbations at the optimization time, which can be expressed as:

    To solve the nonlinear optimization problem [i.e., Eq.(2)], a particle swarm optimization (PSO) intelligent algorithm is utilized in this investigation due to the inaccessible adjoint model in CAM4.The PSO algorithm was inspired by the predation behavior of birds and is effective in solving the nonlinear optimization problem without adjoint models(Kennedy and Eberhart, 1995; Shi and Eberhart, 1999).The specific algorithm has been detailed and illustrated by Mu et al.(2015).However, the rotated empirical orthogonal function is utilized for the Arctic SIC dimension reduction due to the high dimensionality of Arctic SIC perturbations,which is similar to that of Ma et al.(2022).The first 20 leading patterns are used to represent the Arctic SIC variability,which explains more than 80% of the observed daily Arctic SIC variability in boreal winter.The PSO method is considered to have obtained a solution if the objective function value does not increase in the following 15 iterations.We adopted 20 particles to solve the nonlinear optimization problem and the optimal Arctic SIC perturbations are usually obtained in 30 iterations.

    3.Influence of SIC on NAO events

    This investigation explores the influences of Arctic SIC on NAO-event formations on subseasonal timescales.As mentioned in a previous study, the case of successive NAO events or phase transition events is more complicated, and their evolutionary mechanism may be distinct from that in traditional NAO formations (Luo et al., 2014).Therefore, we only consider the traditional NAO event in this investigation, which requires the NAOI to be between –1.0 and 1.0 within 15 days before its onset.According to this criterion,there are five NAO+ and five NAO– events in the 30-year CAM4 simulation.These NAO events also have strong intensities.

    The composite SLP and Z500 anomalies of these five NAO+ events are shown in Fig.1.It shows that the Rossby waves propagate from the North Pacific and evolve to a negative-over-positive structure in the North Atlantic.Regarding the NAO– events, their evolution shows an in situ developmental feature (Fig.2).The local development in NAO–events is consistent with that in observations (Benedict et al.,2004).The evolutionary behaviors of each NAO event are shown in Figs.S1 and S2 in the electronic supplementary materials (ESM).The evolutionary behaviors of the composite simulated NAO events are similar to those in observation,indicating that CAM4 has a good description of the NAO evolutions.Therefore, these cases are used to further explore the Arctic SIC influence on NAO event predictions.

    Similar to those in Dai et al.(2021b), some sensitivity experiments are conducted to determine the response time of NAO events to SIC perturbations.Specifically, 60 monthly Arctic SIC anomalies in winter (December, January, and February) ranging from 1981/82 to 2000/01 are obtained from ERA-Interim.These 60 SIC perturbations have been superimposed on the boundary conditions from 15 days (lag–15 day), 10 days (lag–10 day), and 5 days(lag–5 day) before the NAO onset respectively, and last to 15 days after the NAO onset.It should be noted that these observed Arctic SIC perturbations are added at the start of each integration, and unchanged during the whole integration period.The numerical model is run for 30 days to explore the Arctic SIC influence on NAO events.As shown in Fig.3 and Figs.S3-S4 in the ESM, the NAO shows obvious responses from the lag–5 day with the Arctic SIC perturbed 15 days in advance of the NAO onset.The Arctic SIC perturbations are likely to cause uncertainties in the NAO onset(lag 0 day), as well as the NAO development and decay periods (Figs.3a, b).With a start time 10 days in advance of the NAO onset (lag–10 day), there is little influence on the NAO onset (lag 0 day) but there are many impacts on the NAO development and decay periods (Figs.3c, d).Regarding the Arctic SIC perturbed from 5 days ahead of the NAO onset (lag–5 day), there is almost no obvious impact on the NAO onset or development periods, but there are some uncertainties during their decay periods (Figs.3e, f).That is, the entire NAO life cycle would be influenced if the Arctic SIC has been perturbed 15 days before its onset, indicating that Arctic SIC is a predictability source for subseasonal NAO event predictions.This result is similar to the effect of the sudden thinning of Arctic sea ice (Dai et al., 2021b).Therefore,the Arctic SIC perturbation 15 days in advance of the NAO onset is adopted to explore the Arctic SIC influence on NAO event predictions.

    4.Optimal SIC perturbations for NAO+event predictions

    Fig.1.Composite SLP anomaly (shading, units: hPa) and Z500 anomaly (contour, units: gpm) during the NAO+ events at(a) lag–4, (b) lag–2, (c) lag 0, (d) lag+2, (e) lag+4 and (f) lag+6 day in CAM4 simulation.Lag 0 denotes the NAO+ onset day.Solid and dashed lines represent positive and negative anomalies, respectively.The contour interval (CI) is 40 gpm and the zero line is omitted.The red box in each panel indicates the North Atlantic sector (20°–80°N, 90°W–40°E).

    Fig.2.Same as in Fig.1, but for the five NAO– events in the CAM4 simulation.

    Fig.3.NAOI responses to Arctic SIC perturbations with different lead times.The left panels correspond to Arctic SIC perturbations at (a) 15, (c) 10, and (e) 5 days before the first simulated NAO+ event onset where the x-axis is time, and the y-axis is the NAOI.In each panel, the thick black line indicates the NAOI with no Arctic SIC perturbation, while 60 gray lines correspond to the NAOI derived with different Arctic SIC perturbations.The red dots represent the lead time for superimposing Arctic SIC perturbations, while lag 0 corresponds to the NAO onset day.Panels (b), (d), and (f) are similar to (a), (c), and (e), but for the first simulated NAO– event.

    As mentioned previously, the SIC perturbations have been superimposed on the basic boundary conditions 15 days in advance of the NAO onset (lag–15) and last until the end of model integration.Since we are going to investigate the influence of Arctic SIC perturbations on NAO event predictions on subseasonal timescales, the pentad average has been adopted to remove the high frequency impacts.Consequently, lag–15 to lag–11, lag–10 to lag–6, lag–5 to lag–1,and lag 0 to lag+4 days correspond to pentads 1, 2, 3, and 4,respectively, and pentad 4 is the formation stage for NAO event, including the onset day of NAO event and following four days.Therefore, the corresponding objective function evolves to the difference between the 4th pentad-average NAOI, which can be expressed as:

    4.1.The optimal SIC perturbations

    Fig.4.(a) Type-1 and (b) type-2 optimal SIC perturbations for the first simulated NAO+ event prediction.Panel(c) shows the NAOI evolution derived in the reference state (black), type-1 optimal SIC perturbation (red), and type-2 optimal SIC perturbation (blue).

    Taking the five NAO+ events as the reference states,the corresponding nonlinear optimal problems have been solved.The first NAO+ event is used as an example to illustrate the influences of optimal Arctic SIC perturbations on NAO+ event predictions.For the first NAO+ event, two types of optimal SIC perturbations have been obtained by solving the corresponding nonlinear optimal problem.Type-1 optimal SIC perturbations exhibit positive SIC anomalies mainly in the Greenland Sea, Barents Sea, and Okhotsk Sea(Fig.4a), while type-2 optimal SIC perturbations show negative SIC perturbations in the Bering Sea, Greenland Sea,and Barents Sea (Fig.4b).Superimposing these two types of optimal SIC perturbations on the reference boundary conditions and allowing them to persist for the whole integration period, the corresponding NAOI evolutions are shown in Fig.4c.These two types of SIC perturbations have nearly the same evolutionary NAOI as the reference state in the first two pentads (lag–15 to lag–6 days).However, they have distinct influences on NAO events in pentads 3 and 4(lag–5 to lag 4 days).Type-1 optimal SIC perturbations lead to a smaller NAOI evolution in pentads 3 and 4 compared with the reference state (red line in Fig.4c), underestimating the NAO+ event.Moreover, the NAOI derived with type-1 optimal SIC perturbations is smaller than 1.0 during pentad 4.Therefore, the NAO+ event no longer exists due to type-1 SIC perturbations.However, for type-2 optimal SIC perturbations, the corresponding NAOI has a small difference from the reference state during pentad 3 but exhibits a larger NAOI evolution during pentad 4 (blue line in Fig.4c),which overpredicts the NAO+ event.The influence of type-1 optimal Arctic SIC perturbations is larger than that of type-2 perturbations.Moreover, both types of SIC perturbations exert a much greater influence on NAO event predictions compared with the random SIC perturbations under the same constraint (not shown).This suggests that only perturbations in specific patterns are likely to trigger great impacts on NAO predictions, while randomly structured perturbations have little influence on them.

    Figure 5 shows the SLP and Z500 responses induced by type-1 optimal SIC perturbations.During pentad 1, there is a positive SLP anomaly and negative geopotential height anomaly at 500 hPa in the Barents Sea (Fig.5a).Soon thereafter, the SLP and Z500 anomalies spread in the Arctic(north of 60°N) during pentad 2.However, the SLP anomalies are less than 1 hPa in amplitude and almost outside of the North Atlantic sector, which did not have an obvious influence on the NAO+ event (Fig.5b).For pentad 3, there are some positive SLP anomalies over Baffin Bay, as well as some negative SLP anomalies over the Canadian Arctic Archipelago (Fig.5c).Therefore, a somewhat negative NAOI response is shown in pentad 3.However, the 4th pentad shows a positive SLP response over Greenland and Barents Seas and a negative SLP response over the East Siberian Sea.It also shows a positive geopotential height response over the North Atlantic sector and a negative 500 hPa geopotential height response over East Siberia.The responses in SLP and geopotential height near 500 hPa exhibit a quasibarotropic structure (Fig.5d).Thus, the positive response in the North Atlantic sector weakens the north center of the NAO+ dipole and causes the NAO+ event to cease to exist.However, the SLP and Z500 responses induced by type-2 optimal SIC perturbations have nearly opposite evolutionary behaviors compared to those triggered by type-1 optimal SIC perturbations but with a much smaller amplitude(Fig.S5).

    Fig.5.SLP (shading, units: hPa) and Z500 response (contour, units: gpm) to type-1 optimal SIC perturbations during the (a) 1st pentad (CI = 0.2 gpm), (b) 2nd pentad (CI = 2.0 gpm), (c) 3rd pentad (CI = 20 gpm) and (d) 4th pentad(CI = 40 gpm) for the first simulated NAO+ event.The red box in each panel indicates the North Atlantic sector (20°–80°N, 90°W–40°E).

    4.2.Mechanism of SIC perturbations on NAO+predictions

    Both types of optimal SIC perturbations have notable influences on NAO+ event predictions in all four pentads,and type-1 SIC perturbations even result in the NAO+ event being too weak to exist.Therefore, type-1 SIC perturbations are taken as an example to illustrate the physical mechanism by which Arctic SIC perturbations influence NAO+ event prediction in the four pentads.

    As the Arctic SIC is perturbed, the corresponding heat fluxes between the sea ice and air change.Figure 6 shows the longwave radiation, latent heat flux, and sensible heat flux responses to the type-1 optimal SIC perturbations on the first day.Shortwave radiation has little influence on NAO events during boreal winter and is not discussed here.As shown in Figs.6b and c, there are some negative latent heat flux and sensible heat flux responses in Greenland and Barents Seas.These turbulent heat flux response patterns correspond to positive SIC anomalies, which prevent heat transport from the ocean to the air.The response of longwave radiation is smaller than that of turbulent heat flux (Fig.6a).

    As the heat flux between the sea ice and air changes,the air temperature in the low troposphere causes a consequential response.Figure 7 shows the 700-hPa temperature response for the four pentads.Triggered by type-1 optimal SIC perturbations, there is some cooling over the Barents and Greenland seas during the 1stpentad (Fig.7a).The cooling is consistent with the SIC perturbations and the heat flux responses, suggesting that the 700-hPa temperature responds via the diabatic process.For pentad 2, the temperature response in the North Atlantic sector is small, but there is a significant temperature response around the Canadian Archipelago (Fig.7b).However, for the 3rd pentad, the 700-hPa temperature shows a cooling response in northern Greenland but somewhat of a warming response over the Canadian Archipelago and Barents-Kara Seas.The 4th pentad shows a warming response in Greenland at 700 hPa, which corresponds to the nearly barotropic geopotential height response there.Moreover, there is a negative-over-positive temperature response in northern Siberia and the Barents-Kara Seas,which agrees with the geopotential height response there(Fig.5d).

    The heat budget equation in pressure coordinates is utilized to reveal the temperature evolution induced by the optimal SIC perturbations (Hoskins et al., 1989), which can be formulated as:Figure 8 shows that both the horizontal and vertical temperature advection play important roles in the 700-hPa temperature response in these two pentads.However, the diabatic process is less important in these two pentads than in dynamic terms.That is, the 700-hPa temperature response is mainly influenced by diabatic processes in the 1st pentad but controlled by dynamic processes (horizontal advection and vertical convection) later.Thus, it is not surprising that temperature responses cover the Northern Hemisphere.Furthermore, there is a baroclinic response in Z500 due to the geopotential height increases in the warming sectors and the geopotential height decline in the cooling sectors.The perturbed Z500 could further influence the atmospheric circulation via dynamic processes afterward.

    As mentioned previously, the response in SLP (Z500)is less than 1 hPa (10 gpm) during pentads 1-2, which has little influence on NAOI.However, for pentads 3-4, there are large responses in SLP and Z500, which obviously influence the NAOI evolution.Therefore, the response in SLP and Z500 during pentads 3-4 are further investigated.For the Z500 evolution, its behavior could be explained with the quasi-geostrophic potential vorticity equation (Hoskins,1997), which can be expressed as:

    Fig.6.(a) Longwave radiation, (b) latent heat flux, and (c) sensible heat flux responses (units: W m–2) to type-1 SIC perturbations on the first day.

    Fig.7.700-hPa temperature response (units: K) to type-1 SIC perturbations during the (a) 1st pentad, (b) 2nd pentad,(c) 3rd pentad, and (d) 4th pentad.The red box in each panel indicates the North Atlantic sector (20°–80°N,90°W–40°E).

    Therefore, the influence of type-1 optimal SIC perturbations on this NAO+ event prediction can be concluded as follows.Both the latent heat flux and sensible heat flux experienced local adjustments due to Arctic SIC perturbations on the first day.During the 1st pentad, the temperature in the low troposphere has a local response due to diabatic processes.After that, the dynamic processes (horizontal temperature advection and vertical temperature convection) play crucial roles in modifying the temperature response.The Z500 has a baroclinic response to the temperature with Z500 increasing in warming sectors while decreasing in cooling areas in pentads 1-2.In pentads 3-4, the response in Z500 is mainly influenced by dynamic potential vorticity advection and finally exhibits a positive response in the northern North Atlantic sector.Therefore, the SLP around Greenland increases and underestimates the NAO+ event as a consequence.

    Fig.8.Decomposition of the 700-hPa temperature response during pentads 2 and 3 (units: K d–1).Panel (a) is the horizontal advection term, (b) is the vertical convection term, (c) is for diabatic heating processes, and (d) is the dynamic term (i.e., the sum of horizontal and vertical advection) in pentad 2.Panels (e)–(h) are similar to (a)–(d) but for pentad 3.

    Fig.9.The response of the 500-hPa potential vorticity advection (in terms of geopotential height, units: gpm d–1) in(a) pentad 3 and (b) pentad 4.The red box in each panel indicates the North Atlantic sector (20°–80°N, 90°W–40°E).

    Regarding type-2 optimal SIC perturbations on NAO+event prediction, their influence has been investigated similarly.As shown in Fig.10a, on the first day, the turbulent heat flux increases in response to the Arctic SIC decline in the Greenland Sea, Barents Sea, and Bering Sea (Fig.4b).The turbulent heat flux response warms the temperature in the low troposphere via diabatic processes in pentad 1.After that, the temperature response is mainly affected by horizontal temperature advection and vertical convective processes.During pentads 1-2, Z500 has a baroclinic response to the temperature by raising Z500 in warm sectors and declining Z500 in cooling areas.For pentads 3-4, the atmospheric dynamic process (potential vorticity advection) dominates the Z500 response.Finally, the 500-hPa geopotential height decreases in Greenland, as does the SLP in pentad 4.This strengthens the northern center of the NAO+ dipole and overpredicts the NAO+ event (Fig.10f).event predictions are similar to those in the first NAO+event, which are not further stated here.

    5.Optimal SIC perturbations for NAO–events

    4.3.OptimalSIC perturbations for other NAO+ event predictions

    Similar processes have been applied to explore the influences of Arctic SIC on the remaining four NAO+ event predictions.There are two types of optimal Arctic SIC perturbations for each NAO+ event, mainly with positive SIC perturbations in the Barents Sea and Okhotsk Sea for type-1 and negative SIC perturbations in these regions for type-2(Fig.11).Similarly, type-1 optimal SIC perturbations would underestimate the NAO+ event and even make the event disappear, while type-2 optimal SIC perturbations overpredict the NAO+ events.However, it should be noted that there are some differences among the optimal SIC perturbations for these NAO+ event predictions.For example, there are obvious positive SIC perturbations in the Bering Sea in the type-1 optimal SIC perturbations for the 3rd NAO+ event,while there are obvious negative SIC perturbations in the Bering Sea in the type-2 optimal SIC perturbations for the 2nd and 4th NAO+ events.These differences are mainly due to the various initial atmospheric conditions while the optimal SIC perturbations have to match with the corresponding atmospheric states.Further diagnoses reveal that the mechanisms of optimal Arctic SIC perturbations on these NAO+

    In addition to the NAO+ events, the optimal Arctic SIC perturbations for NAO– event predictions are investigated.Similar to the investigations for NAO+ event predictions,the optimal SIC perturbations for the five NAO– event predictions are obtained by solving the corresponding nonlinear optimal problems.

    Taking the first NAO– event as an example, two types of optimal Arctic SIC perturbations are obtained and shown in Figs.12a and b.They mainly show negative SIC anomalies in type-1 optimal SIC perturbations while positive SIC anomalies emerge in type-2 optimal SIC perturbations.Specifically, type-1 optimal SIC perturbations mainly show negative SIC perturbations in the BKS and west Okhotsk Sea, with some positive SIC perturbations in the Bering Sea.However,it mainly shows positive SIC perturbations in the Greenland Sea, BKS, and the east Okhotsk Sea in type-2 optimal SIC perturbations.Regarding their influences on NAO– event predictions, there is little influence on NAOI during the first three pentads (Fig.12c).However, for the 4th pentad, the type-1 optimal SIC perturbation results in a higher NAOI than that in the reference state, which underestimates the NAO– event.However, type-2 optimal SIC perturbations result in a lower NAOI than that in the reference state leading to an overprediction of the NAO– event.This is also obvious in the SLP and Z500 response (Fig.S6).

    Fig.10.Physical mechanisms associated with type-2 optimal Arctic SIC perturbations on the first simulated NAO+ event prediction.Panel (a) is the turbulent heat flux response on day 1 (units: W m–2), while panels (b) and (c) illustrate the 700-hPa temperature response (shading, units: K) and Z500 response (CI = 0.2 gpm for (b) and CI = 4.0 gpm for (c)) during pentads 1 and 2, respectively.Panels (d) and (e) are the 500-hPa potential vorticity advection (in terms of geopotential height, units: gpm d–1)during pentads 3 and 4, respectively.Panel (f) is the SLP response (shading, units: hPa)and Z500 response (CI = 40 gpm) in pentad 4.The red box in panels indicates the North Atlantic sector (20°–80°N, 90°W–40°E).

    Fig.11.Similar to those in Fig.4, but for the (a) (e) (i) 2nd, (b) (f) (j) 3rd, (c) (g) (k) 4th, and (d) (h) (l) 5th simulated NAO+ event predictions.

    To illustrate the mechanisms of these Arctic SIC perturbations on this NAO– event prediction, similar processes have been investigated.As shown in Fig.13a, there is increasing turbulent heat flux from ocean to air due to the negative SIC perturbations in the Greenland Sea, BKS, and Okhotsk Sea.The turbulent heat flux further warms the air in the lower troposphere due to the diabatic process during pentad 1.After that, the temperature in the lower troposphere is modified by horizontal advection and vertical convective processes,thus further influencing the 500-hPa geopotential height in pentad 2.During pentads 3-4, the Z500 response is mainly affected by the potential vorticity advection.Finally, type-1 optimal SIC perturbations trigger a negative-over-positive 500-hPa geopotential height dipole as well as the SLP in the North Atlantic sector, which underestimates the NAO–event in pentad 4 (Fig.13d).Similarly, type-2 optimal Arctic SIC perturbations overpredict the NAO– event in pentad 4 but with smaller amplitudes compared with type-1(Figs.13e–h).

    Moreover, the optimal Arctic SIC perturbations for the remaining four NAO– event predictions have also been investigated similarly.The two types of optimal Arctic SIC perturbations and their influences on NAO– event predictions are shown in Fig.14.It is obvious that type-1 optimal SIC perturbations for these NAO– event predictions mainly show negative SIC anomalies in the Greenland Sea and Barents Sea,and they would like to underestimate the NAO– events and even disappear.However, the type-2 optimal SIC perturbations mainly show positive SIC anomalies in the Barents and Okhotsk Seas, and they would likely lead to a lower NAOI in pentad 4 and overpredict the NAO– events.However, there are still some differences among these optimal SIC perturbations.For example, there are some positive SIC anomalies in the Bering Sea for type-1 optimal SIC perturbations in the 2nd, 3rd, and 5th NAO– events, while type-2 optimal SIC perturbations for the 2nd NAO– event show some negative SIC anomalies in the Okhotsk and Bering Seas.These differences may be due to the different atmospheric initial conditions among these NAO– events.Further diagnoses show that the mechanisms of the optimal SIC perturbations on the corresponding NAO– event predictions are similar to those in the first NAO– event.This mechanism could be interpreted as the turbulent heat flux, acting in response to the SIC perturbations, serving to modify the low tropospheric temperature via diabatic processes in the first pentad.Horizontal and vertical temperature advection dominate the temperature response in pentads 2-3, allowing Z500 to have a consequent baroclinic response to the temperature perturbations.Subsequently, the dynamic potential vorticity advection becomes dominant in pentads 3-4.Finally, the Z500 demonstrates a positive (negative) geopotential height response in the northern North Atlantic sector, as does the SLP, which underestimates (overpredicts) the NAO– event.

    Fig.12.(a) Type-1 and (b) type-2 optimal SIC perturbations for the first simulated NAO– event predictions.Panel (c)is the NAOI evolution derived in the reference state (black), type-1 optimal SIC perturbations (red), and type-2 optimal SIC perturbations (blue).

    6.Conclusion and discussion

    In this investigation, we focused on the influence of Arctic SIC perturbations on NAO event predictions on subseasonal timescales.Utilizing the CAM4 model with climatological SST and SIC, we found that the model well describes NAO event evolution on subseasonal timescales in the control simulations.Five NAO+ and five NAO– events were selected from the control simulations and were further investigated to explore the influence of Arctic SIC on their predictions.

    First, the response time for NAO events to Arctic SIC perturbations was determined with a series of sensitivity numerical experiments.Numerical results show that NAO demonstrates notable responses to Arctic SIC perturbations in 10 days and SIC perturbations have large influences on NAO formations when the SIC perturbations are 15 days ahead of the onset, indicating that Arctic SIC is a predictability source for NAO events on subseasonal timescales.After that, the CNOP-B approach was utilized to find the optimal Arctic SIC perturbations, which have the greatest impacts on NAO event prediction.The optimal SIC perturbations were obtained by solving the corresponding nonlinear optimal problem for each NAO event.Numerical results show that there are two types of optimal Arctic SIC perturbations for each NAO event, with type-1 optimal SIC perturbations tending to underestimate the NAO event and type-2 optimal SIC perturbations tending to overpredict the event in 4 pentads.The optimal SIC perturbations have greater influences on NAO event predictions compared with the random SIC perturbations under the same physical constraint.For NAO+events, type-1 optimal SIC perturbations mainly show positive SIC anomalies in the Greenland Sea, Barents Sea, and Okhotsk Sea, while type-2 mainly exhibit negative SIC anomalies in these regions.However, for NAO– events,type-1 optimal SIC perturbations mainly show negative SIC anomalies in Greenland and the Barents Seas, while type-2 optimal SIC perturbations exhibit positive SIC anomalies in the Greenland Sea, BKS, and Okhotsk Sea.The type-1 optimal SIC perturbations almost share an opposite pattern to the type-2 optimal SIC perturbations for both NAO+ and NAO– events, although there are some slight differences among these optimal SIC perturbations due to the different atmospheric initial conditions.These optimal SIC perturbations also correspond to the ice-edge regions with strong SIC variations and large uncertainties in boreal winter, highlighting that sea-ice anomalies in both the Atlantic (Luo et al., 2016) and Pacific (Mesquita et al., 2011) sectors are important for the NAO.

    Fig.13.Physical mechanism for optimal Arctic SIC perturbations on NAO– event predictions.Panel (a) is the turbulent heat flux response on day 1 (units: W m–2) and (b) is the 700-hPa temperature response (shading, units: K) and 500-hPa geopotential height response (CI = 1.0 gpm) in pentad 2.Panel (c) is the potential vorticity advection (in terms of geopotential height, units: gpm d-1) at 500 hPa during pentad 4.Panel (d) is the SLP response (shading, units: hPa) and Z500 response (CI = 40 gpm) in pentad 4.Panels (e)–(h)are similar to (a)–(d) but for type-2 optimal SIC perturbations.The red box in each panel indicates the North Atlantic sector (20°–80°N,90°W–40°E).

    Fig.14.Similar to those in Fig.12, but for (a) (e) (i) 2nd, (b) (f) (j) 3rd, (c) (g) (k) 4th, and (d) (h) (l) 5th simulated NAO–event predictions.

    Further diagnosis reveals that the optimal Arctic SIC perturbations modify the turbulent heat flux on the first day and affect the temperature in the low troposphere via diabatic processes.After that, the temperature in the lower troposphere is mainly modified by the dynamic temperature advection in pentad 2.During pentads 3-4, dynamic processes such as potential vorticity advection play a crucial role in modifying the Z500 response and finally form positive or negative 500-hPa geopotential height anomalies in pentad 4.The SLP demonstrates an equivalent barotropic response to that at 500 hPa and further influence the NAO intensity in pentad 4.These results are dependent upon the reference state of NAO events.As for the reference states without any NAO events, the solution to the optimization problem refers to the optimal precursor if it can trigger an NAO event onset.This will be explored in our future work.

    Similar to previous studies, our results highlight the role of Arctic SIC on NAO event predictions.From a climatological perspective, many previous studies have revealed that Arctic SIC can affect NAO activities by modifying the temperature meridional gradients of temperature and potential vorticity as well as planetary wave activities in the North Atlantic sectors (Overland et al., 2011; Jaiser et al.,2012; Luo et al., 2019).Their results also indicate that Arctic SIC loss would result in a negative or negligible NAO response compared with atmospheric internal variability(Blackport and Screen, 2021; Siew et al., 2021; Warner et al., 2020).However, our results suggest that on subseasonal timescales, the NAO has positive responses to negative SIC anomalies in the Greenland Sea, Barents Sea, and Okhotsk Sea, and vice versa, while local dynamic and thermodynamic processes are crucial for NAO event prediction.The influence of thermal forcing on extratropical atmospheric circulations is complicated and includes multiple eddy feedback processes and flow dependencies.This indicates that both Arctic SIC and atmospheric internal variability are important for subseasonal NAO event prediction.We also highlight the matching mechanism between SIC perturbations and the atmospheric initial state for NAO event evolution.Moreover, our result indicates that the Greenland Sea, Barents Sea, and Okhotsk Sea may be sensitive areas for subseasonal NAO event predictions, which require further elaborate investigations in the future (Dai et al., 2019; Zhang et al., 2019, 2020).

    It should be further noted that atmospheric initial perturbations are also crucial for subseasonal NAO predictions since subseasonal ensemble forecast systems usually utilize specific atmospheric initial perturbations generated by various approaches to improve forecast skills (Vitart et al.,2017).However, the influence of atmospheric initial perturbations on subseasonal NAO predictions is not discussed here.Nevertheless, our result indicates that the specific Arctic SIC perturbations may further improve the subseasonal NAO prediction skill.Additionally, as pointed out by many investigations, the stratosphere is a key process in linking the Arctic and mid-latitudes (Sun et al., 2015, 2022; Zhang et al., 2018).Their results indicate that the elaborate stratospheric process in numerical models plays a crucial role in describing the seasonal or longer linkage between the Arctic and mid-latitudes.But for timescales of 4 pentads, the lowtop and coarse stratosphere in CAM4 may only have a negligible influence on the Arctic and mid-latitude linkage.We may further confirm our results with high-top models in the future.However, the feedback from atmospheric circulations to SIC variabilities cannot be estimated, since CAM4 is an atmospheric-only circulation model and the Arctic SIC is treated as a boundary condition.As revealed by Strommen et al.(2022), air–ice coupling might be important for simulating sea ice–NAO teleconnections.Therefore, the relationship between the NAO and sea ice should be further investigated with a fully coupled numerical model.Investigations of this nature may help to deepen the understanding of the mechanism of NAO event formation and provide scientific support for subseasonal NAO event predictions.

    AcknowledgementsThe authors are grateful to the editor and two anonymous reviewers for their helpful comments in improving this work.This study was supported by the National Natural Science Foundation of China (Grant Nos.42288101, 41790475, 42005046,and 41775001).The calculations within this work were performed on TianHe-2 and the Atmospheric-Oceanic Numerical Simulation Platform.The authors appreciate the support of the National Supercomputer Center in Guangzhou (NSCC-GZ) and the High Performance Computing Center in the Department of Atmospheric and Oceanic Sciences, Fudan University.

    Electronic supplementary material:Supplementary material is available in the online version of this article at https://doi.org/10.1007/s00376-023-2371-8.

    Data Availability StatementThe ERA-Interim reanalysis data are accessible from the European Centre for Medium-range Weather Forecasts (https://apps.ecmwf.int/datasets/), and the CAM4 model is accessible from the University Corporation for Atmospheric Research (https://www.cesm.ucar.edu/models/?ref=hp).

    欧美激情在线99| 欧美激情在线99| 精品一区二区三区视频在线观看免费| 在线免费观看不下载黄p国产 | 久久久久久久亚洲中文字幕 | 亚洲天堂国产精品一区在线| 国产黄色小视频在线观看| 日本免费一区二区三区高清不卡| av在线老鸭窝| 欧美xxxx黑人xx丫x性爽| 国产蜜桃级精品一区二区三区| 亚洲 国产 在线| 在线观看66精品国产| 赤兔流量卡办理| 老司机午夜十八禁免费视频| 国产精品自产拍在线观看55亚洲| 搡老岳熟女国产| 国产欧美日韩精品一区二区| 亚洲av电影在线进入| 欧美黑人欧美精品刺激| 美女高潮喷水抽搐中文字幕| av中文乱码字幕在线| 欧美日韩综合久久久久久 | 久久婷婷人人爽人人干人人爱| 亚洲电影在线观看av| 高清在线国产一区| www.www免费av| 亚洲五月天丁香| 亚洲欧美清纯卡通| 九九热线精品视视频播放| 国产一级毛片七仙女欲春2| 国产成人av教育| 久久中文看片网| 日本黄色视频三级网站网址| 色综合欧美亚洲国产小说| 搡女人真爽免费视频火全软件 | 91狼人影院| 亚洲人成电影免费在线| 精品乱码久久久久久99久播| 亚洲中文字幕一区二区三区有码在线看| 可以在线观看毛片的网站| 亚洲美女搞黄在线观看 | 俺也久久电影网| 丰满的人妻完整版| 亚洲经典国产精华液单 | 亚洲人成网站在线播放欧美日韩| 亚洲不卡免费看| 三级毛片av免费| 久久午夜亚洲精品久久| 性插视频无遮挡在线免费观看| 国产视频内射| 俺也久久电影网| 欧美xxxx黑人xx丫x性爽| 欧美成人a在线观看| 嫩草影院精品99| 俄罗斯特黄特色一大片| 国产一区二区激情短视频| 国产精品免费一区二区三区在线| 亚洲熟妇中文字幕五十中出| 精品久久国产蜜桃| 久久久久久久久久黄片| 久久中文看片网| 国产精品一区二区三区四区久久| 美女黄网站色视频| 亚洲国产欧洲综合997久久,| 欧美性猛交╳xxx乱大交人| 69av精品久久久久久| 国模一区二区三区四区视频| 亚洲黑人精品在线| 国产伦人伦偷精品视频| 亚洲最大成人av| 亚洲内射少妇av| 亚洲av成人av| 国产免费av片在线观看野外av| 国产真实伦视频高清在线观看 | av专区在线播放| 久久久精品大字幕| 久久人妻av系列| 国产成人a区在线观看| 又爽又黄a免费视频| 国产精品久久久久久久久免 | 婷婷亚洲欧美| 久久热精品热| 精品人妻视频免费看| 欧美精品国产亚洲| or卡值多少钱| 在线观看午夜福利视频| 久久中文看片网| 亚洲精品456在线播放app | 天天躁日日操中文字幕| 日韩欧美精品免费久久 | 伊人久久精品亚洲午夜| 亚洲色图av天堂| 在线观看66精品国产| 欧美最新免费一区二区三区 | 色在线成人网| 97人妻精品一区二区三区麻豆| 欧美又色又爽又黄视频| 精品人妻视频免费看| 男人舔女人下体高潮全视频| 99在线人妻在线中文字幕| 丰满乱子伦码专区| 亚洲熟妇熟女久久| 精品久久久久久久久av| a级毛片a级免费在线| 精品久久久久久久末码| 亚洲欧美日韩高清专用| 亚洲av.av天堂| 成人高潮视频无遮挡免费网站| 国产欧美日韩精品一区二区| 午夜日韩欧美国产| 免费一级毛片在线播放高清视频| 看黄色毛片网站| 日韩精品中文字幕看吧| 成人永久免费在线观看视频| 欧美区成人在线视频| 国产欧美日韩精品一区二区| 好看av亚洲va欧美ⅴa在| 国产熟女xx| 能在线免费观看的黄片| 美女xxoo啪啪120秒动态图 | 少妇裸体淫交视频免费看高清| 亚洲欧美日韩无卡精品| 亚洲国产精品合色在线| 直男gayav资源| 亚洲在线观看片| 欧美国产日韩亚洲一区| 性欧美人与动物交配| 精品欧美国产一区二区三| 亚洲内射少妇av| 丰满人妻熟妇乱又伦精品不卡| 国产精品久久久久久人妻精品电影| 免费高清视频大片| 一进一出好大好爽视频| 日韩中字成人| 国产爱豆传媒在线观看| 成人精品一区二区免费| 国产大屁股一区二区在线视频| 亚洲第一欧美日韩一区二区三区| 国产精品一区二区三区四区久久| 欧美三级亚洲精品| 亚洲欧美日韩东京热| 特级一级黄色大片| netflix在线观看网站| 成熟少妇高潮喷水视频| 欧美黄色淫秽网站| 夜夜爽天天搞| 国产精品爽爽va在线观看网站| 精品久久久久久久久av| 我的女老师完整版在线观看| 无遮挡黄片免费观看| 日韩欧美 国产精品| 国产欧美日韩精品一区二区| 免费黄网站久久成人精品 | 精品无人区乱码1区二区| 日韩欧美精品免费久久 | 直男gayav资源| 国产精品亚洲av一区麻豆| 真人一进一出gif抽搐免费| 午夜激情欧美在线| 天美传媒精品一区二区| 免费看美女性在线毛片视频| 90打野战视频偷拍视频| 成年女人永久免费观看视频| 久久九九热精品免费| 色5月婷婷丁香| 高潮久久久久久久久久久不卡| 婷婷丁香在线五月| 亚洲黑人精品在线| 亚洲aⅴ乱码一区二区在线播放| 长腿黑丝高跟| 亚洲激情在线av| 亚洲精品日韩av片在线观看| 97碰自拍视频| 亚洲中文字幕一区二区三区有码在线看| 欧美精品啪啪一区二区三区| 国产精品乱码一区二三区的特点| 国产欧美日韩精品一区二区| 成人一区二区视频在线观看| 别揉我奶头 嗯啊视频| 夜夜夜夜夜久久久久| 亚洲国产色片| 亚洲综合色惰| 看黄色毛片网站| 乱码一卡2卡4卡精品| 欧洲精品卡2卡3卡4卡5卡区| 国产在视频线在精品| 亚洲精品456在线播放app | 在线观看一区二区三区| 波野结衣二区三区在线| 午夜免费男女啪啪视频观看 | 12—13女人毛片做爰片一| 99视频精品全部免费 在线| 丁香六月欧美| 国产精品嫩草影院av在线观看 | 性欧美人与动物交配| 性色avwww在线观看| 极品教师在线免费播放| 午夜福利高清视频| 国产又黄又爽又无遮挡在线| 国产中年淑女户外野战色| 亚洲精品影视一区二区三区av| 中文字幕熟女人妻在线| 变态另类丝袜制服| 国产精品免费一区二区三区在线| 热99在线观看视频| 亚洲国产精品999在线| 精品国产亚洲在线| 亚洲经典国产精华液单 | 婷婷六月久久综合丁香| 国产国拍精品亚洲av在线观看| 天堂动漫精品| 长腿黑丝高跟| 国产精品爽爽va在线观看网站| 有码 亚洲区| 色哟哟·www| 国产一区二区三区在线臀色熟女| 亚洲熟妇中文字幕五十中出| 黄片小视频在线播放| 国产伦精品一区二区三区视频9| 亚洲激情在线av| 亚洲av免费在线观看| 搡老熟女国产l中国老女人| 丁香六月欧美| 夜夜看夜夜爽夜夜摸| 久久久成人免费电影| 少妇的逼好多水| 欧美xxxx黑人xx丫x性爽| 精品日产1卡2卡| 欧美+日韩+精品| 国产真实伦视频高清在线观看 | 日日夜夜操网爽| 人妻夜夜爽99麻豆av| www.999成人在线观看| 美女cb高潮喷水在线观看| 精品人妻一区二区三区麻豆 | 欧美性猛交╳xxx乱大交人| 性色avwww在线观看| 直男gayav资源| 午夜两性在线视频| 午夜a级毛片| 亚洲aⅴ乱码一区二区在线播放| 亚洲一区高清亚洲精品| 午夜精品在线福利| 免费看日本二区| 国产精品电影一区二区三区| 久久久久久久久久成人| 国产亚洲精品久久久久久毛片| 国内揄拍国产精品人妻在线| 69av精品久久久久久| 欧洲精品卡2卡3卡4卡5卡区| 99视频精品全部免费 在线| 亚洲精品在线观看二区| 国产三级黄色录像| 国产美女午夜福利| 在线国产一区二区在线| 日韩欧美国产一区二区入口| 国产三级黄色录像| 日韩有码中文字幕| av专区在线播放| 日韩人妻高清精品专区| 成年女人永久免费观看视频| 动漫黄色视频在线观看| 又黄又爽又刺激的免费视频.| 在线观看66精品国产| 香蕉av资源在线| 久久6这里有精品| 欧美国产日韩亚洲一区| 亚洲精品一卡2卡三卡4卡5卡| 亚洲国产精品sss在线观看| 国产精品三级大全| 真人一进一出gif抽搐免费| 国产色爽女视频免费观看| 日本成人三级电影网站| 欧美黑人欧美精品刺激| a在线观看视频网站| 久久精品国产亚洲av天美| 特级一级黄色大片| 午夜福利欧美成人| 欧美黑人巨大hd| 一本精品99久久精品77| 在线观看免费视频日本深夜| 国产老妇女一区| 国产精品野战在线观看| 国产精品乱码一区二三区的特点| 午夜免费成人在线视频| 亚洲av成人av| 狂野欧美白嫩少妇大欣赏| 成人精品一区二区免费| 欧美+亚洲+日韩+国产| 男女床上黄色一级片免费看| 亚洲最大成人av| 女人十人毛片免费观看3o分钟| 国产一区二区三区在线臀色熟女| 亚洲人成电影免费在线| 国产成人a区在线观看| 99久久久亚洲精品蜜臀av| 亚洲综合色惰| 欧美+亚洲+日韩+国产| 男女做爰动态图高潮gif福利片| 国产精品女同一区二区软件 | 男女视频在线观看网站免费| 色尼玛亚洲综合影院| 三级男女做爰猛烈吃奶摸视频| 免费电影在线观看免费观看| 国产真实乱freesex| 亚洲av中文字字幕乱码综合| 久久人人精品亚洲av| 国内毛片毛片毛片毛片毛片| 免费看日本二区| 国产精品不卡视频一区二区 | 999久久久精品免费观看国产| 五月玫瑰六月丁香| 免费人成视频x8x8入口观看| 一本一本综合久久| 男女那种视频在线观看| 国产精品98久久久久久宅男小说| 丰满的人妻完整版| 毛片一级片免费看久久久久 | 青草久久国产| 18禁黄网站禁片午夜丰满| 丁香六月欧美| 2021天堂中文幕一二区在线观| 91在线精品国自产拍蜜月| 欧美3d第一页| 白带黄色成豆腐渣| 亚洲片人在线观看| а√天堂www在线а√下载| 日本免费a在线| 久久99热6这里只有精品| 亚洲欧美日韩东京热| 一进一出好大好爽视频| 久9热在线精品视频| 小蜜桃在线观看免费完整版高清| 午夜福利在线在线| 精品人妻视频免费看| 国产黄a三级三级三级人| 欧美成狂野欧美在线观看| 一级黄片播放器| 黄色丝袜av网址大全| 99视频精品全部免费 在线| 久久久久久久久久黄片| 日本免费a在线| 午夜激情欧美在线| 亚洲中文日韩欧美视频| 国产爱豆传媒在线观看| 免费高清视频大片| 美女cb高潮喷水在线观看| 在线国产一区二区在线| 无遮挡黄片免费观看| 日韩av在线大香蕉| 国产av一区在线观看免费| 精品一区二区三区人妻视频| 国产精品av视频在线免费观看| 精品久久久久久久久亚洲 | 91麻豆精品激情在线观看国产| 给我免费播放毛片高清在线观看| 亚洲中文日韩欧美视频| 波多野结衣高清无吗| 成人特级av手机在线观看| 国产精品美女特级片免费视频播放器| 看免费av毛片| 精品一区二区三区人妻视频| 十八禁网站免费在线| 麻豆国产av国片精品| 国内揄拍国产精品人妻在线| 国产黄片美女视频| 一a级毛片在线观看| 亚洲性夜色夜夜综合| 在线看三级毛片| 少妇高潮的动态图| 九九久久精品国产亚洲av麻豆| 在线观看午夜福利视频| 真实男女啪啪啪动态图| 欧美日韩黄片免| 九九久久精品国产亚洲av麻豆| 久久6这里有精品| 欧美3d第一页| 亚洲色图av天堂| 国产三级黄色录像| 日韩欧美国产在线观看| 一进一出抽搐gif免费好疼| 男女之事视频高清在线观看| 精品一区二区三区av网在线观看| 日韩欧美精品v在线| 夜夜爽天天搞| 在线播放国产精品三级| 可以在线观看的亚洲视频| 精品人妻一区二区三区麻豆 | 精品久久国产蜜桃| 99精品久久久久人妻精品| 91在线观看av| av国产免费在线观看| 性色avwww在线观看| 国产亚洲精品久久久久久毛片| 精品久久久久久,| 亚洲精品影视一区二区三区av| 亚洲av成人精品一区久久| 亚洲专区中文字幕在线| 欧美日韩乱码在线| 日韩欧美在线二视频| 欧美日韩黄片免| 女人十人毛片免费观看3o分钟| 亚洲中文日韩欧美视频| 久久婷婷人人爽人人干人人爱| 成人性生交大片免费视频hd| 内射极品少妇av片p| 色综合亚洲欧美另类图片| 国产精品嫩草影院av在线观看 | 欧美性猛交黑人性爽| 超碰av人人做人人爽久久| 校园春色视频在线观看| 夜夜看夜夜爽夜夜摸| 综合色av麻豆| 日韩欧美在线二视频| 亚洲第一区二区三区不卡| 国产一区二区三区在线臀色熟女| 国产aⅴ精品一区二区三区波| 搡老岳熟女国产| 精品熟女少妇八av免费久了| 变态另类成人亚洲欧美熟女| 国产精品女同一区二区软件 | 又爽又黄无遮挡网站| 欧美绝顶高潮抽搐喷水| 九九久久精品国产亚洲av麻豆| 狂野欧美白嫩少妇大欣赏| 麻豆成人av在线观看| 嫩草影院新地址| 国产激情偷乱视频一区二区| 日韩欧美一区二区三区在线观看| 男人狂女人下面高潮的视频| 日韩中文字幕欧美一区二区| 亚洲最大成人av| 久久久久精品国产欧美久久久| av天堂在线播放| 亚洲欧美清纯卡通| 欧美日韩国产亚洲二区| 亚洲av免费在线观看| 1000部很黄的大片| 免费一级毛片在线播放高清视频| 91字幕亚洲| 亚洲欧美日韩卡通动漫| 亚洲人成网站高清观看| 亚洲精品在线观看二区| 两个人视频免费观看高清| 中文在线观看免费www的网站| av国产免费在线观看| 我的女老师完整版在线观看| 一级黄片播放器| 在线a可以看的网站| 久久久久久国产a免费观看| 欧美日韩亚洲国产一区二区在线观看| 三级国产精品欧美在线观看| 在线观看美女被高潮喷水网站 | 欧美xxxx性猛交bbbb| 国产毛片a区久久久久| 宅男免费午夜| 国产成人啪精品午夜网站| 亚洲精品一区av在线观看| 十八禁人妻一区二区| 欧美激情国产日韩精品一区| 国产精品一区二区三区四区久久| 99热6这里只有精品| 久久久久国产精品人妻aⅴ院| 天堂动漫精品| 成人亚洲精品av一区二区| 91午夜精品亚洲一区二区三区 | 91在线精品国自产拍蜜月| 九色国产91popny在线| 五月伊人婷婷丁香| 很黄的视频免费| 亚洲精品在线美女| 欧美在线一区亚洲| 在线播放无遮挡| 亚洲精品粉嫩美女一区| 国产色爽女视频免费观看| 美女高潮喷水抽搐中文字幕| 国产精品影院久久| av在线蜜桃| 午夜老司机福利剧场| 成人国产综合亚洲| 色哟哟·www| 真实男女啪啪啪动态图| 一级a爱片免费观看的视频| 国产精品一区二区三区四区久久| 欧美国产日韩亚洲一区| 精品人妻1区二区| 国产精品嫩草影院av在线观看 | 国产欧美日韩一区二区精品| 免费av不卡在线播放| 俄罗斯特黄特色一大片| 怎么达到女性高潮| av欧美777| 日本撒尿小便嘘嘘汇集6| 国产一区二区三区视频了| 最近中文字幕高清免费大全6 | 最近最新免费中文字幕在线| 看黄色毛片网站| 丰满人妻一区二区三区视频av| 99久国产av精品| 自拍偷自拍亚洲精品老妇| 国产欧美日韩一区二区精品| 白带黄色成豆腐渣| 国产一区二区在线av高清观看| 欧美日韩亚洲国产一区二区在线观看| 久久久久久九九精品二区国产| 一边摸一边抽搐一进一小说| 亚洲一区二区三区不卡视频| 我的女老师完整版在线观看| 特大巨黑吊av在线直播| 中文字幕av成人在线电影| 桃红色精品国产亚洲av| 三级男女做爰猛烈吃奶摸视频| 日本a在线网址| 一本综合久久免费| 欧美一区二区国产精品久久精品| 国产大屁股一区二区在线视频| eeuss影院久久| 久久99热6这里只有精品| 欧美高清性xxxxhd video| 嫩草影视91久久| 国产淫片久久久久久久久 | 成人国产一区最新在线观看| 真人一进一出gif抽搐免费| 久久久久国产精品人妻aⅴ院| 久久久精品欧美日韩精品| 窝窝影院91人妻| 哪里可以看免费的av片| 国产av一区在线观看免费| 3wmmmm亚洲av在线观看| 国产精品嫩草影院av在线观看 | 俺也久久电影网| 免费无遮挡裸体视频| 一区二区三区四区激情视频 | 男人狂女人下面高潮的视频| 国产午夜精品久久久久久一区二区三区 | 88av欧美| 色播亚洲综合网| 成年版毛片免费区| 亚洲精品久久国产高清桃花| 精品福利观看| 亚洲av一区综合| 欧美潮喷喷水| 国产成年人精品一区二区| 我的女老师完整版在线观看| 黄色女人牲交| 久久伊人香网站| 精品一区二区免费观看| 久久久久精品国产欧美久久久| 国产免费av片在线观看野外av| 欧美国产日韩亚洲一区| 精品一区二区免费观看| 久久伊人香网站| 久久亚洲精品不卡| 久久精品国产99精品国产亚洲性色| 大型黄色视频在线免费观看| 青草久久国产| 国产视频一区二区在线看| 亚洲国产精品成人综合色| 欧美+日韩+精品| 欧美性猛交黑人性爽| 国产在线男女| 国产精品自产拍在线观看55亚洲| 日日摸夜夜添夜夜添小说| 床上黄色一级片| 国产精品野战在线观看| 日韩欧美免费精品| 国产精品99久久久久久久久| 神马国产精品三级电影在线观看| 成年人黄色毛片网站| 午夜福利在线观看免费完整高清在 | 丁香六月欧美| 亚洲自偷自拍三级| 免费在线观看日本一区| 少妇的逼好多水| 色综合婷婷激情| 国产黄色小视频在线观看| 久久精品国产亚洲av涩爱 | 男人的好看免费观看在线视频| 免费av不卡在线播放| 性插视频无遮挡在线免费观看| 国产白丝娇喘喷水9色精品| 天天躁日日操中文字幕| 夜夜看夜夜爽夜夜摸| av女优亚洲男人天堂| 国产精品综合久久久久久久免费| 51国产日韩欧美| 午夜福利在线在线| 两人在一起打扑克的视频| www日本黄色视频网| 99久久成人亚洲精品观看| 欧美精品国产亚洲| 欧美潮喷喷水| 亚洲av一区综合| 久久国产精品影院| 久久久国产成人免费| 97碰自拍视频| 午夜福利成人在线免费观看| 欧美黑人欧美精品刺激| 精品久久久久久久人妻蜜臀av| 免费看美女性在线毛片视频| 国产精品国产高清国产av| 亚洲国产精品成人综合色| 国产老妇女一区| 亚洲精品久久国产高清桃花| 老熟妇仑乱视频hdxx| 一进一出好大好爽视频| 老师上课跳d突然被开到最大视频 久久午夜综合久久蜜桃 | 一a级毛片在线观看| 亚洲精品乱码久久久v下载方式| 国产精品免费一区二区三区在线| 最近最新中文字幕大全电影3| 国产成年人精品一区二区| 天堂网av新在线| 亚洲av电影在线进入|