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

    典型燃料點火延遲時間的一階和二階局部和全局敏感度分析

    2019-03-08 08:30:40席雙惠李象遠
    物理化學學報 2019年2期

    席雙惠,王 繁,*,李象遠

    1四川大學原子與分子物理研究所,成都 610065

    2四川大學化學工程學院,成都 610065

    1 Introduction

    Sensitivity analysis (SA)1–9is an important tool in model validation and evaluation and it provides quantitative information on importance of input parameters in a model on its outputs.The first-order sensitivity coefficient describes effect of one input parameter on model output,while the second-order sensitivity coefficient represents interactive effect between two parameters on the output1.SA has been employed extensively in analysis of chemical kinetic models for combustion processes.In those analysis,input parameters in the model are reaction rate parameters,thermodynamic data or transportation parameters and outputs of the model are combustion characters predicted by the model such as ignition delay time,laminar flame speeds,and concentration of critical species under various combustion conditions.Sensitivity analysis is also closely related to uncertainty analysis1,2,5,6.There is always some uncertainties in the input parameters of a chemical kinetic model.Effects of these uncertainties on the predicted combustion properties can be revealed with the help of sensitivity analysis.On the other hand,uncertainties of the input parameters should be considered when carrying out sensitivity analysis.Such uncertainty analysis has become rather popular recently10–16.In addition,SA has also been employed to develop reduced mechanisms by removing redundant species and reactions through local rate sensitivity analysis or directed relation graph aided sensitivity analysis17–23.

    Sensitivity analysis is divided into local sensitivity analysis and global sensitivity analysis.In local sensitivity analysis1–5,response of the model output due to a small change of the input parameters from their nominal values is evaluated.Local sensitivity coefficients for concentration of species,temperature or pressure in a chemical kinetic model satisfy certain differential equations in a spatially homogeneous,dynamic system and they can be determined by solving the corresponding equations numerically.This approach has been implemented in software packages such as Senkin24.On the other hand,when sensitivity of some global combustion characters such as ignition delay time is required,brute-force method is generally employed25.Sensitivity coefficient of an input parameter is determined from difference between the original model output and the output with the input parameter varied by a fixed percentage in this method.A total number of (m + 1) model runs need to be carried out to determine the first-order local sensitivity coefficients of m input parameters.In principle,both first- and second-order local sensitivity coefficients can be obtained.However,only first-order local sensitivity analysis has been reported previously for chemical kinetic model of combustion.In those first-order local sensitivity analyses,the input parameters are allowed to vary only in a small range so that linear sensitivity coefficients can be calculated.However,nonlinear effect could be important particularly when the input parameter can change in a relatively large range,or when uncertainty of an input parameter is large.In addition,coupled effects of two input parameters are not considered.

    To circumvent these problems,global sensitivity analysis can be carried out1–3,5,10,16,26–30.In global sensitivity analysis,the input parameters are allowed to change in a much larger space,or the whole uncertainty range and coupling between two input parameters is also evaluated.A large number of model runs are usually required in global sensitivity analysis where all the investigated input parameters vary simultaneously within the whole range of their possible values.Several methods have been developed to perform global sensitivity analysis8,10,11,16,26–28,31,32Global sensitivity analysis through an optimized high-dimensional model representation method using a random sampling of inputs(RS-HDMR) has been implemented in the GUI_HDMR program by Ziehn and Tomlin10,26,29–35.A relatively small number of samples is needed to carry out global sensitivity analysis using the GUI_HDMR program and it is presently one of the most popular global sensitivity analysis method in analysis of combustion models26,29,30,33,34,36–46.Recently,the artificial neural network method28has been adopted to further reduce computational effort in HDMR method.However,a relatively large number of samples are still required compared with local sensitivity analysis and global sensitivity analysis is rather expensive particularly when the number of input parameters is large.In practical calculations,global sensitivity coefficients with respect to the most important 20–30 input parameters in the chemical models for combustion processes are calculated with GUI_HDMR26.In a chemical kinetic model for combustion process,there usually exist several hundreds or even thousands of input parameters and the most important 20–30 input parameters are chosen based on physical consideration or firstorder local sensitivity coefficients.Third-order or even higherorder sensitivity coefficients can be determined in principle with HDMR,however,only first and second-order global sensitivity coefficients are calculated in practice26,29,30,33,34,36–46.

    The output of the model is written as a polynomial function of the input parameters in RS-HDMR and global sensitivity coefficients are calculated from this polynomial function26.A large number of samples is employed to obtain the parameters in this polynomial function through Monte-Carlo integration in RSHDMR.In fact it is possible to determine this polynomial function by only varying two input parameters simultaneously if the first- and second-order sensitivity coefficients are to be determined.When this polynomial function is available,firstand second-order sensitivity coefficients in both local and global sensitivity analysis can readily be achieved.The number of samples can be reduced dramatically with this approach and it is also possible to obtain sensitivity coefficients for a larger number of input parameters.Global sensitivity coefficients determined in this way are different from those using RS-HDMR,where some averaged effects are included.In this work,we propose to calculate first and second order sensitivity coefficients of ignition delay time of some hydrogen carbon fuels with respect to reaction rate constants in corresponding combustion kinetic mechanisms using this approach.These sensitivity coefficients can provide valuable information on effects of a single reaction as well as interactive effects of two reactions on ignition delay time.

    This paper is organized in the following manner:the method employed in this work as well as computational procedures to obtain sensitivity coefficients are provided in Method Section.The first and second-order sensitivity coefficients for ignition delay time of some important fuels:H2,methane,n-butane and n-heptane under different combustion conditions are provided in Results and Discussion Section.Sensitivity coefficients in local sensitivity analysis are compared with those in global sensitivity analysis in this section.Conclusion will be drawn in Conclusion Section.

    2 Method

    In the HDMR method for global sensitivity analysis,the output of the model f(x) is represented as a function of the input parameters based on the idea of analysis of variance(ANOVA)40decomposition using the following form:

    In Eq.(1),n is the number of input parameters,f0is the zeroth order component function which represents the mean effect of the input parameters on the output of the model,fi(xi) is the first order component function and represents independent effects of each input parameter xion the output f(x),fij(xi,xj) is the second order component function which describes correlated effects of both xiand xjon f(x).f1,2,...,n(x1,x2,...,xn) in Eq.(1) accounts for the n’th order correlated effects on f(x).For most systems,the HDMR expression up to the second-order can provide a satisfactory approximation for f(x)47:

    In previous works,each input parameter is always rescaled to be in the range of [0,1].We prefer to rescale the input parameters so that all xisatisfy -1 ≤ xi≤ 1 and the input parameter is zero when it takes its nominal value in this work.Advantage of rescaling the input parameters in this way will be manifested in the following.

    In the GUI-HDMR software,the first- and second-order component functions are further approximated from linear combinations of a set of orthonormal polynomials {φp(x)} over the range of the scaled input parameters:

    where u,v and w are the orders in the polynomial expansions,are the expansion coefficients.In the present work,we employ exactly the same relations,i.e.,Eqs.(2)–(4) to approximate relationship between the output of the model and the input parameters.Legendre polynomials which are orthogonal over the range of [-1,1] are chosen for the set of orthonormal polynomials {φp(x)} and they satisfy the following equations:

    The first and second-order global sensitivity coefficients are defined from the ratio between the partial variances and the total variance using the following equations:

    Siis the first order global sensitivity coefficient,or the main effect of the input variable xion the output;and Sijis the secondorder global sensitivity coefficient,or the interaction effect of xiand xjon the output.We only consider the first- and the secondorder sensitivity coefficients in calculations and all these global sensitivity coefficients add up to 1:

    One can see from the above equations that the global sensitivity coefficients are determined from the expansion coefficientsin Eqs.(3) and (4).According to their definition in Eqs.(7)–(10),global sensitivity coefficients are always positive.In the RS-HDMR approach,these coefficients40,41as well as the zeroth order term f0are determined from Monte Carlo integration over the input sample.

    In fact,f(x) expressed in Eqs.(2)–(4) is fully equivalent to an expression expanded with a set of ordinary polynomials:

    Considering that xiis zero when the corresponding input parameter takes its nominal value,it can readily be seen from Eq.(12) that B0is the value of the output when all the input parameters take their nominal values.In addition,can be determined from samples where only the i’th input parameter is different from its nominal value,whilecan be obtained from samples by only varying the i’th and j’th input parameters.When the expansion coefficientsandare determined,the expansion coefficientscan be calculated using the following equations:

    The first- and second-order global sensitivity coefficients can thus be achieved fromusing Eqs.(8)–(10) and Eqs.(12)–(14).To determine,we only need to vary one- or two-input parameters in the model runs.This means firstand second-order global sensitivity analysis can possibly be carried out from samples where only one or two input parameters change their values at the same time.The number of samples can be reduced significantly using this procedure.

    In local sensitivity analysis,the output of the model with respect to the input parameters is expressed using the following Taylor expansion3:

    We focus on first- and second-order sensitivity coefficients of ignition delay time with respect to reaction rate constants in this work.The brute-force approach is usually employed to obtain first-order local sensitivity coefficients in previous works using the following equation:

    In the present sensitivity analysis,the output of the model f(x)is chosen to be lnτignand the input parameters are the natural logarithm of the factor that scales the reaction rate constant.The input parameters are restricted in the range of [-1,1] and the factor that scales the reaction rate constant is thus in the range of[1/e,e] where e is the base of the natural logarithm.It is straightforward to scale the reaction rate constants in a larger range,or in its uncertainty range and still map the input parameters in the range of [-1,1].

    According to Eq.(12),B0is simplyin this work.To determine,we only need to run simulations on ignition delay times with scaled reaction rate constant for the i-th reaction.can be obtained from ignition delay times with scaled reaction rate constants for both the i-th and the j-th reaction.In the present work,the orders u,v and w in the polynomial expansions are set to be 4.The input parameters only vary in a small range and setting the orders to be 4 should provide reliable sensitivity coefficients.Higher order may be required if the input parameters are allowed to change in a much larger range.are obtained from ignition delay times with the i-th reaction rate con stant scaled by 1/2,1/2,2and 2.On the other hand,are determined from ignition delay times with both the ith and the j-th reaction scaled by these four factors.The expansion coefficientsare calculated using Eqs.(13) and(14) and the first- and second-order local and global sensitivity coefficients are achieved from these expansion coefficients.To account for the first- and second-order sensitivity coefficients of m input parameters,a total number of (8m(m - 1) + 4m + 1)model runs are required.One can also choose smaller numbers for u,v and w in Eq.(12) initially and increase them only when necessary.Number of model runs can possibly be reduced in this way.

    3 Results and discussion

    First and second-order global and local sensitivity coefficients for ignition delay time of H2,methane,n-butane and n-heptane with respect to reaction rate constants at various initial conditions are studied using the approaches described in Method Section.Global sensitivity coefficients calculated from Eqs.(8)–(14) will be termed as GSA,while the first-order global sensitivity coefficients with only theterm in Eq.(13) are represented with GSA-1.Global sensitivity coefficients with GSA-1 are those only considering the first order terms in global sensitivity analysis.LSA denotes traditional local sensitivity coefficients based on Eqs.(16) and (17) with the constant c taken to be e,while the local sensitivity coefficients determined fromare represented by LSA-1.Lists of the ten important reactions or reaction pairs on ignition delay times of these fuels under various conditions are given in tables in supporting information and they are discussed in detail in the following.

    3.1 Sensitivity analysis on ignition delay time of H2

    A small-size mechanism for combustion ofincluding 10 species and 33 reactions is employed for sensitivity analysis.Ignition delay times in a zero-dimensional closed homogeneous reactor27with constant pressure are simulated under the initial conditions of T = 1000 K,p = 1 Pa,φ = 1.0;T = 1000 K,p = 40 Pa,φ = 1.0;T = 1400 K,p = 1 Pa,φ = 1.0;and T = 1400 K,p = 40 Pa,φ = 1.0,respectively.First- and second order sensitivity analysis with respect to reaction rate constants of the 33 reactions were carried out for ignition delay times under these four conditions.

    3.1.1 First-order sensitivity coefficients

    Reactions with the largest absolute GSA and LSA for ignition delay time of H2under these four conditions are illustrated in Fig.1.It can be seen from this figure that the most important reactions on ignition delay times obtained from GSA are almost the same as those with LSA,except for reaction No.11 (H +O2(+ O2) = HO2(+ O2)).This reaction is rather important using LSA on ignition delay time at T = 1400 K and p = 40 Pa,while it is of minor importance with GSA.In fact,this reaction is shown to be the third most important reaction on ignition delay time at this condition using GSA-1,where the correlated effects of two input parameters are not considered.According to our results,difference between important reactions using GSA-1,LSA and LSA-1 are negligible for H2.This shows results of firstorder global sensitivity analysis using the present approach are affected to some extent by presence of higher order terms in expansion for f(x) in some cases.

    One can see from Fig.1 that the chain branching reaction:H +O2= OH + O is always one of the most important reaction on ignition delay time.On the other hand,HO2and H2O2play a more important role in combustion of H2at higher pressures.These results are consistent with previous findings on combustion processes of.

    3.1.2 Second-order sensitivity analysis

    Pairs of reactions that have the largest correlated effect on ignition delay times with GSA and LSA are illustrated in Fig.2.Results of GSA and LSA are again quite similar to each other.One can see from this figure that correlated effect of the two chain branching reaction:No.16 (H + O2= OH + O) and No.15(O+ H2= OH + H) are the most important reaction pair on ignition delay time at high temperature and normal pressure condition,i.e.,T = 1400 K,p = 1 Pa.On the other hand,reaction No.16 (H + O2= OH + O) and No.9 (H + O2(+M) = HO2(+M)) are the most important reaction pair on ignition delay time under the other conditions.First-order sensitivity coefficients of reactions No.9 and 16 are the largest under the conditions of T= 1000 K,p = 1 Pa,φ = 1.0 and T = 1400 K,p = 40 Pa,φ = 1.0,while they are the third and fourth most important reactions with T = 1000 K,p = 40 Pa,φ = 1.0,according to our first-order sensitivity analysis.On the other hand,reactions Nos.15 and 16 are the most important two reactions under the condition of T =1400 K,p = 1 Pa,φ = 1.0.This seems to indicate that correlated effect of two reactions is important when the two reactions also have large first-order sensitivity coefficient.However,reaction No.22 (H + HO2= H2+ O2) and reaction No.16 (H + O2= OH+ O) have sizable correlated effect on ignition delay time under the condition of T = 1000 K,p = 40 Pa,φ = 1.0,while reaction No.22 does not have a pronounced effect on ignition delay time under this condition.This shows it is possible that correlated effect of an important reaction and a minor reaction could still be important on ignition delay time.

    Fig.1 First-order global and local sensitivity coefficients of important reaction on ignition delay time of H2.

    3.2 Sensitivity analysis on ignition delay time of methane

    The GRI-Mech 3.0 mechanism is one of the most popular reaction mechanisms for combustion of methane and it contains 53 species and 325 reactions.Number of model runs is too large if first- and second-order sensitivity analysis are carried out for all the involved reactions even with the present approach.Skeletal mechanism reduction is performed before mechanism analysis.A skeletal mechanism is obtained based on intersection50of skeletal mechanisms obtained from direct relation graph(DRG)51–53,DRG with error propagation (DRGEP)54,revised DRG approach55,path flux analysis(PFA)56,and flux projection tree method (FPT)57approaches over range of temperature of 1000–1600 K,pressures of 1–20 Pa,and equivalence ratios of 0.5–1.5.The skeletal mechanism employed in sensitivity analysis contains 22 species and 109 reactions with a maximum error of 10.5% on ignition delay times compared with the detailed mechanism on these conditions.In addition,computational cost can be reduced when the skeletal mechanism is employed in simulation of ignition delay time compared with that using the detailed mechanism,particularly when a large number of model runs are required in sensitivity analysis.Firstand second order sensitivity coefficients of the 109 reactions in the skeletal mechanism are calculated for ignition delay times in a zero-dimensional closed homogeneous reactor with constant pressure under the initial condition of T = 1000 K,p = 1 Pa,φ =1.0;T = 1000 K,p = 20 Pa,φ = 1.0;T = 1400 K,p = 1 Pa,φ =1.0;and T = 1400 K,p = 20 Pa,φ = 1.0,respectively.

    Fig.2 Second-order global and local sensitivity coefficients of important pairs of reactions on ignition delay time of H2.

    3.2.1 First-order sensitivity coefficients

    Reactions that have the largest absolute sensitivity coefficients on ignition delay time under the four conditions in GSA and LSA are illustrated in Fig.3.It can be seen that important reactions identified with GSA are consistent with those using LSA in most cases.In fact,all the four methods,i.e.,GSA,GSA-1,LSA and LSA-1,provide quite similar information on important reactions in most cases except for the case of T = 1400 K and p = 20 Pa.Reaction No.93 CH3+ H2O2=HO2+ CH4,is the second most important reaction on ignition delay time at T = 1400 K and p = 20 Pa with GSA.However,effect of this reaction on ignition delay time at this condition is calculated to be unimportant with the other three approaches.Results of GSA depend on how the orthonormal polynomials are chosen,while LSA does not have such a problem.Difference between GSA and GSA-1 arises from the second-order term mathematically.In addition,the input parameters are allowed to vary in a relatively small range,which means generally LSA and GSA should not differ too much.The fact that the other three approaches give similar results indicates that GSA may overestimate importance of this reaction.

    Fig.3 First-order global and local sensitivity coefficients of important reaction on ignition delay time of CH4.

    According to results in Fig.3,the recombination reaction 2CH3(+ M) = C2H6(+ M) is always an important reaction on ignition delay times of CH4under all the considered conditions.One can see from LSA results that this recombination reaction suppresses ignition due to annihilation of radicals.On the other hand,global sensitivity coefficients are always positive and it is not straightforward to tell whether the reaction actually promotes or suppresses ignition only with global sensitivity coefficients.Reaction No.72:HO2+ CH3= OH + CH3O becomes more important to promote ignition at higher pressures.On the other hand,CH3+ O2= OH + CH2O is rather important on ignition at 1000 K and 1 Pa,while reaction CH3+ O2= O + CH3O becomes more important at 1400 K and 1 Pa.

    3.2.2 Second-order sensitivity coefficients

    Pairs of reactions that have important correlated effect on ignition delay time of methane using GSA and LSA are illustrated in Fig.4.It can be seen from this figure that second-order sensitivity coefficients of many pairs of reactions are quite similar.This is particularly the case for ignition delay time at T = 1400 K and p = 20 Pa.It is thus difficult to figure out certain pairs of reactions that have significantly more important correlated effects on ignition delay time than the other pairs.An interesting observation is that the pairs of reactions that have the largest second-order global sensitivity coefficients always contain reaction No.94 at T = 1400 K and p = 20 Pa in GSA,which also has the largest first-order global sensitivity coefficient at this condition.

    Fig.4 Second-order global and local sensitivity coefficients of important pairs of reactions on ignition delay time of CH4.

    In most cases the pairs of reactions that are important on ignition delay time using GSA are also identified as important reaction pairs with LSA,although the sequence of the absolute second-order sensitivity coefficients with GSA may be different from that using LSA.The pairs of reactions that are important on ignition delay times are usually composed of reactions that have large first-order sensitivity coefficients.However,reaction pair of No.71 HO2+ CH3= O2+ CH4and No.72 HO2+ CH3= OH +CH3O are shown to have large correlated effect on ignition delay time at T = 1000 K,although reaction No.71 has only a rather small first-order sensitivity coefficient.lnτ at T = 1000 K and p =1 Pa at different scaling factor for reaction No.71 is plotted in Fig.5a and it can be seen that ignition delay time does not change monotonically with respect to the scaling factor.To further investigate sensitivity of reaction No.71 on ignition delay time,we treat the forward and backward reactions for HO2+ CH3=O2+ CH4as two irreversible reactions and only the forward or the backward reaction is scaled in simulation of ignition delay time.The obtained ignition delay times at the temperature of 1000 K and pressure of 20 Pa with respect to the scaling factor for the forward and the backward reactions are plotted in Fig.5b,c,respectively.We can see that the forward reaction suppresses ignition,while the backward reaction actually promotes ignition.First-order sensitivity analysis is also carried out by treating the forward and backward reactions as two irreversible reactions.We found that these two reactions have sizeable first-order sensitivity coefficients.This can explain that this reaction and another important reaction could possibly have pronounced correlated effect on ignition delay time.This seems to imply that reversible reactions have better be treated as two irreversible reactions in sensitivity analysis.However,the forward reaction rate constant and the backward reaction rate constant are connected by the equilibrium constant.If uncertainty of the equilibrium constant is quite small or if we are not interested in sensitivity of thermodynamic parameters,it is reasonable to treat reversible reaction as a single reaction in sensitivity analysis.

    Fig.5 Ignition delay time versus scaling factor for rate constant of reaction No.71 in reaction mechanism of methane at T = 1000 K and p = 20 Pa.

    Fig.6 First-order global and local sensitivity coefficients of important reaction on ignition delay time of n-butane.

    3.3 Sensitivity analysis on ignition delay time of nbutane

    The detailed reaction mechanism for combustion of n-butane58containing 230 species and 1328 reactions is reduced before sensitivity analysis.A large number of model runs are required to perform sensitivity analysis and it will be too demanding to use the detailed mechanisms in the model runs.We thus choose to carry out sensitivity analysis using skeletal mechanisms derived from these detailed mechanisms.FPT is employed for mechanism reduction over range of temperature of 600–1600 K,pressures of 1 Pa,and equivalence ratios of 1.0.The resulting skeletal mechanism contains 80 species and 420 reactions and the largest relative error on ignition delay time with the skeletal mechanism is 8.8%,which is able to provide a reliable description on combustion of this fuel.It is still rather demanding to carry out sensitivity analysis on all the reactions in the skeletal mechanism.The 100 reactions with the largest absolute first-order local sensitivity coefficients from traditional local sensitivity analysis are chosen for global and local sensitivity analysis under the initial conditions of T = 650 K,p =1 Pa,φ = 1.0;T = 1000 K,p = 1 Pa,φ = 1.0;and T = 1400 K,p= 1 Pa,φ = 1.0.It should be noted that only about 20–30 reactions are chosen in global sensitivity analysis with RSHDMR.The 100 reactions under one initial condition is different from those under another initial condition.

    Fig.7 Second-order global and local sensitivity coefficients of important pairs of reactions on ignition delay time of n-butane.

    The reactions with the largest absolute first-order sensitivity coefficients on ignition delay time at T = 650 K,1000 K and 1400 K according to GSA and LSA are listed in Fig.6.According to this figure,important reactions with these four sensitivity analysis methods are quite similar to each other.At low temperature,decomposition of H2O2to form OH and H-abstraction reaction of the fuel by OH are the most important two reactions that promote ignition.Reactions involving C4H9OO and C4H8OOH also play an important role in ignition at this temperature.On the other hand,two reactions between CH3and HO2have the largest effect on ignition delay time at 1000 K.Recombination of the two radicals CH3and HO2to produce stable molecules CH4and O2suppress ignition,while reaction between these two radicals to produce CH3O and OH is the most important reaction that promote ignition at this temperature.One can see from this figure that a bunch of reactions have similar absolute sensitivity coefficients on ignition delay time at low temperature,while the chain branching reaction H + O2= OH +O has significantly larger absolute sensitivity coefficients than the other reactions at high temperature.

    The pairs of reactions that have important correlated effects on ignition delay times of n-butane at T = 650 K,1000 K and 1400 K,respectively,according to GSA and LSA are displayed in Fig.7.One can see from this figure that the four pairs of reactions that affect ignition delay time at T = 650 K the most according to GSA all contain reaction No.294:C4H10+ OH =pC4H9+ H2O.This reaction also has the largest absolute firstorder sensitivity on ignition delay time at this temperature.Furthermore,most of the important reaction pairs for ignition delay time at 650 K are composed of two reactions that have large first-order sensitivity coefficients.On the other hand,correlated effect between reaction No.16:H2O2(+M) =2OH(+M) and another reaction is minor although this reaction ranks No.2 and No.1 in first order GSA and LSA,respectively.The most important reaction pair on ignition delay time at 1000 K is reaction no.16 H2O2(+M) = 2OH(+M) and reaction No.294:C4H10+ OH = pC4H9+ H2O.First-order sensitivity coefficients of these two reactions rank Nos.5 and 6 in GSA,respectively.At high temperature,the chain branching reaction H + O2= OH + O together with another reaction have important correlated effect on ignition delay time.In fact,H + O2= OH +O and reaction No.128:C2H4+ H(+M) = C2H5(+M) are the most important reaction pair on ignition delay time according to our results,although reaction No.128 has a rather small first-order sensitivity coefficient.Its first-order sensitivity coefficient ranks No.21 in GSA and 23 in LSA.

    Fig.8 First-order global and local sensitivity coefficients of important reaction on ignition delay time of n-heptane.

    3.4 Sensitivity analysis on ignition delay time of nheptane

    The detailed mechanism for combustion of n-heptane59,60employed in this work contains 561 species and 2539 reactions.A skeletal mechanism based on intersection of skeletal mechanisms obtained from DRG,DRGEP and FPT over the range of temperature of 600–1800 K,pressures of 1–40 Pa,and equivalence ratios of 0.5–1.5 is adopted in sensitivity analysis.The skeletal mechanism is composed of 175 species and 799 reactions and the maximum error on ignition delay times of these conditions is 11%.We calculated first and second-order global and local sensitivity coefficients of the 80 reactions with the largest first-order local sensitivity coefficients on ignition delay times under the initial conditions of T = 650 K,p = 1 Pa,φ = 1.0;T = 1000 K,p = 1 Pa,φ = 1.0;and T = 1400 K,p = 1 Pa,φ = 1.0.

    Reactions with the largest absolute first-order sensitivity coefficients from GSA and LSA are illustrated in Fig.8.Reactions that have significant effect on ignition delay times at these conditions according to GSA are also shown to be important with LSA for combustion of n-heptane.Our results show that the chain termination reaction:HO2+ OH = H2O + O2has the largest effect to suppress ignition at low temperature.On the other hand,decomposition of H2O2to produces OH is the most important reaction that promotes ignition and this reaction also has large sensitivity coefficient on ignition delay time of nbutane at low temperature.However,sensitivity coefficients for reactions involving peroxyl radical of n -heptane are not very large compared with those of the listed reactions.Reactions involving ROO such as the C7H15O2–2= C7H14OOH2–4ranked 8th and 12th respectively in the global and local sensitivity analysis.On the other hand,sensitivity coefficient of reaction O2QOOH = keto-hydroperoxide + OH ranks around 20th.It is well-known that reactions involving ROO and O2QOOH play an important role in low temperature ignition of alkanes.C7H15OO has several isomers that undergo similar reactions.Sensitivity coefficient for any single reaction involving ROO or O2QOOH is not among the largest ones,but combinational effect of reactions involving isomers of ROO could still be significant.At intermediate temperature,decomposition reaction of H2O2still plays an important role on ignition.Similar to the case of ignition of n-butane at 1000 K,the two reactions between CH3and HO2are rather important on ignition delay time of n-heptane.At high temperature,absolute sensitivity coefficient of the chain branch reaction:H + O2= OH + O are much larger than that of the other reactions.

    Pairs of reactions that have the largest absolute second-order sensitivity coefficients on ignition delay times are listed in Fig.9.It can be seen from this figure that reactions involving nheptane are rather important in second-order sensitivity analysis for ignition at 650 K,while they are less important in first-order sensitivity analysis at this temperature.On the other hand,the most important pairs of reactions usually contain reactions related to CH3at intermediate temperature.At 1400 K,reaction No.8 H+O2= OH + O always appear in the pairs of reactions that have the largest absolute second-order sensitivity coefficients.In addition,first-order sensitivity coefficient of reaction No.564 n-C7H16+ H = C7H15-3+ H2is quite small,but this reaction together with reaction No.8 constitute the reaction pair that has the largest second-order sensitivity coefficient at 1400 K.

    Fig.9 Second-order global and local sensitivity coefficients of important pairs of reactions on ignition delay time of n-heptane.

    4 Conclusions

    We propose to calculate both first- and second-order global and local sensitivity coefficients of ignition delay time with respect to reaction rate constants of H2,methane,n-butane and n-heptane at various conditions in this work.Ignition delay time with respect to the factors that scale reaction rate constants is expanded using a normal polynomial function with at most coupling between two reactions.The polynomial function is determined from ignition delay times by scaling rate constants of one and two reactions and both global and local sensitivity coefficients are determined from this polynomial function.Number of samples required to determine global sensitivity coefficients can be reduced compared with global sensitivity analysis using RS-HDMR.In RS-HDMR,sensitivity coefficients are determined for rate constants of a limited number of reactions,while the present approach is able to afford sensitivity coefficients for a larger number of reactions.

    Reactions and reaction pairs with the largest sensitivity coefficients are listed for ignition delay time of these four fuels.The factor that scale reaction rate constants vary in a small range and important reactions or reaction pairs identified with global sensitivity analysis are usually the same as those from local sensitivity analysis.This is probably because nonlinear effect is still unimportant in this small range of parameters.It is possible to determine global sensitivity coefficients by varying the input parameters in a larger range using the present approach.Our results show that correlated effects between an important reaction and a minor reaction could also have significant secondorder sensitivity coefficient in some cases.On the other hand,first-order global sensitivity coefficients with the present approach will be affected by coupling between two reactions and occasionally some result with first-order global sensitivity analysis will be different from those with local sensitivity analysis or global sensitivity analysis by neglecting the correlated effects of two reactions.The present sensitivity analysis approach can provide valuable information on important reactions as well as correlated effects of two reactions reaction on combustion characters obtained from chemical kinetic mechanisms of combustion.Furthermore,it can also be employed to aid global sensitivity analysis using RS-HDMR where global sensitivity coefficients are determined more reliably.

    此物有八面人人有两片| 国产精品二区激情视频| 日韩欧美免费精品| 成年人黄色毛片网站| 在线观看免费午夜福利视频| 国内久久婷婷六月综合欲色啪| 一级毛片精品| 在线视频色国产色| 亚洲熟妇中文字幕五十中出| 黄频高清免费视频| 亚洲一卡2卡3卡4卡5卡精品中文| 亚洲av电影不卡..在线观看| 精品久久久久久久人妻蜜臀av| 激情在线观看视频在线高清| 美女大奶头视频| 在线天堂中文资源库| 免费女性裸体啪啪无遮挡网站| 久久性视频一级片| 午夜久久久在线观看| 亚洲 欧美一区二区三区| 超碰成人久久| 久久久久久人人人人人| 亚洲第一青青草原| 亚洲欧美一区二区三区黑人| 欧美日韩福利视频一区二区| 久久热在线av| 成人亚洲精品av一区二区| 欧美色欧美亚洲另类二区| 日韩有码中文字幕| 一本久久中文字幕| 国产v大片淫在线免费观看| 免费观看精品视频网站| 亚洲人成网站在线播放欧美日韩| 男女做爰动态图高潮gif福利片| 亚洲av电影在线进入| 日韩大码丰满熟妇| 91麻豆av在线| 国产激情欧美一区二区| 亚洲人成网站在线播放欧美日韩| 少妇被粗大的猛进出69影院| 亚洲aⅴ乱码一区二区在线播放 | 亚洲成人久久爱视频| 特大巨黑吊av在线直播 | 色哟哟哟哟哟哟| 在线观看免费视频日本深夜| 欧美日韩精品网址| 看黄色毛片网站| 窝窝影院91人妻| 日韩中文字幕欧美一区二区| 日本成人三级电影网站| 无限看片的www在线观看| www.自偷自拍.com| 一级a爱片免费观看的视频| 午夜福利在线在线| 欧美日韩中文字幕国产精品一区二区三区| 久久精品成人免费网站| 国产单亲对白刺激| 亚洲av片天天在线观看| 一级毛片女人18水好多| 成人国语在线视频| 欧美性长视频在线观看| 午夜两性在线视频| 波多野结衣高清无吗| 男人舔女人下体高潮全视频| 淫秽高清视频在线观看| 久久中文字幕一级| 精品国产乱码久久久久久男人| 久久久久久久精品吃奶| 99久久国产精品久久久| 国产精品永久免费网站| 国产成人精品久久二区二区免费| 亚洲精品一卡2卡三卡4卡5卡| 国产爱豆传媒在线观看 | 啦啦啦免费观看视频1| 国产91精品成人一区二区三区| 日本a在线网址| 这个男人来自地球电影免费观看| 久久久久久亚洲精品国产蜜桃av| 怎么达到女性高潮| 精品日产1卡2卡| 久久热在线av| 曰老女人黄片| 国内少妇人妻偷人精品xxx网站 | 一个人观看的视频www高清免费观看 | 少妇裸体淫交视频免费看高清 | 亚洲av电影在线进入| 国产成人精品久久二区二区91| 一边摸一边抽搐一进一小说| 久久香蕉激情| 久久久国产成人精品二区| 一本精品99久久精品77| 欧美日本视频| 日韩中文字幕欧美一区二区| 亚洲激情在线av| 久久久久国内视频| 日日摸夜夜添夜夜添小说| 97碰自拍视频| 香蕉丝袜av| 两人在一起打扑克的视频| 午夜日韩欧美国产| 精品国产一区二区三区四区第35| 久久精品成人免费网站| 国产99久久九九免费精品| 久久久久久免费高清国产稀缺| 久久天堂一区二区三区四区| 99久久国产精品久久久| 欧美不卡视频在线免费观看 | 免费看美女性在线毛片视频| 欧美日韩亚洲国产一区二区在线观看| 老汉色∧v一级毛片| 国产欧美日韩精品亚洲av| 欧美日韩福利视频一区二区| 又黄又粗又硬又大视频| 亚洲在线自拍视频| 久久久久久人人人人人| 极品教师在线免费播放| 高清在线国产一区| 侵犯人妻中文字幕一二三四区| 国产91精品成人一区二区三区| 亚洲色图 男人天堂 中文字幕| 看免费av毛片| 欧美+亚洲+日韩+国产| 欧美性猛交╳xxx乱大交人| 91麻豆av在线| 亚洲 欧美 日韩 在线 免费| 国产精品精品国产色婷婷| 久久久久国产一级毛片高清牌| 亚洲一区二区三区不卡视频| 国产av一区二区精品久久| 桃红色精品国产亚洲av| 亚洲av成人一区二区三| 1024香蕉在线观看| 国内精品久久久久久久电影| 一边摸一边做爽爽视频免费| www.999成人在线观看| 日本一区二区免费在线视频| 女性被躁到高潮视频| 国产免费av片在线观看野外av| 国产亚洲av嫩草精品影院| 国产精品电影一区二区三区| 亚洲国产毛片av蜜桃av| 成人手机av| 日日摸夜夜添夜夜添小说| 亚洲成人精品中文字幕电影| 国产精华一区二区三区| 九色国产91popny在线| 又黄又粗又硬又大视频| 亚洲第一电影网av| 日本在线视频免费播放| 亚洲av第一区精品v没综合| 日日干狠狠操夜夜爽| 国产精品美女特级片免费视频播放器 | www国产在线视频色| 韩国av一区二区三区四区| 色av中文字幕| 久久国产亚洲av麻豆专区| 激情在线观看视频在线高清| 日韩中文字幕欧美一区二区| 在线视频色国产色| 亚洲国产看品久久| 免费高清视频大片| 精品久久久久久久末码| 91麻豆av在线| 国产精品乱码一区二三区的特点| 欧美日韩黄片免| 精品国产一区二区三区四区第35| 欧美日韩中文字幕国产精品一区二区三区| 久久九九热精品免费| 一卡2卡三卡四卡精品乱码亚洲| 黑人操中国人逼视频| 黑人巨大精品欧美一区二区mp4| 国产精品久久久av美女十八| 大型av网站在线播放| 亚洲欧洲精品一区二区精品久久久| 白带黄色成豆腐渣| 韩国精品一区二区三区| 精品国产一区二区三区四区第35| 亚洲成av片中文字幕在线观看| 麻豆av在线久日| 午夜福利视频1000在线观看| 午夜福利欧美成人| 亚洲成人国产一区在线观看| 欧美激情 高清一区二区三区| 久久天堂一区二区三区四区| 免费人成视频x8x8入口观看| 狂野欧美激情性xxxx| 99精品在免费线老司机午夜| 亚洲成人免费电影在线观看| 国产蜜桃级精品一区二区三区| 欧美在线一区亚洲| 免费在线观看成人毛片| 日日干狠狠操夜夜爽| 国产精品,欧美在线| 国产野战对白在线观看| 国产成人欧美| 国产成人精品久久二区二区免费| 亚洲av成人一区二区三| av免费在线观看网站| 精品国产亚洲在线| 国产一区在线观看成人免费| 国产黄a三级三级三级人| 桃红色精品国产亚洲av| 久久久久国产一级毛片高清牌| 成人一区二区视频在线观看| 激情在线观看视频在线高清| 亚洲欧美日韩无卡精品| 黄色视频不卡| 国产亚洲av嫩草精品影院| 国产精品永久免费网站| 757午夜福利合集在线观看| 午夜福利在线在线| a级毛片a级免费在线| 黄色片一级片一级黄色片| 看片在线看免费视频| 亚洲精品久久成人aⅴ小说| 男人操女人黄网站| 激情在线观看视频在线高清| 欧美日韩亚洲综合一区二区三区_| 国产av一区在线观看免费| 大型黄色视频在线免费观看| 99riav亚洲国产免费| 亚洲av成人一区二区三| 哪里可以看免费的av片| 美女扒开内裤让男人捅视频| 99国产极品粉嫩在线观看| 最近最新免费中文字幕在线| 国产精品亚洲美女久久久| 久久中文看片网| av片东京热男人的天堂| 999久久久国产精品视频| 亚洲av熟女| 亚洲av中文字字幕乱码综合 | 免费看日本二区| 亚洲精品在线美女| 免费人成视频x8x8入口观看| 欧美乱色亚洲激情| 琪琪午夜伦伦电影理论片6080| 99热这里只有精品一区 | 欧美成人午夜精品| 亚洲av熟女| 一区二区三区国产精品乱码| 深夜精品福利| 十分钟在线观看高清视频www| 亚洲中文字幕日韩| 免费人成视频x8x8入口观看| 老司机深夜福利视频在线观看| 99在线人妻在线中文字幕| 国产亚洲欧美精品永久| 亚洲真实伦在线观看| 午夜免费激情av| 欧美国产精品va在线观看不卡| 亚洲av成人av| 久久伊人香网站| 亚洲精品在线观看二区| 国产精品 国内视频| 色综合婷婷激情| 亚洲精品久久成人aⅴ小说| 欧美在线一区亚洲| 两人在一起打扑克的视频| 精品国产乱子伦一区二区三区| 亚洲精品一区av在线观看| 变态另类丝袜制服| 男人操女人黄网站| 夜夜夜夜夜久久久久| 又大又爽又粗| 久久伊人香网站| 一级黄色大片毛片| 国产av一区二区精品久久| a在线观看视频网站| 亚洲国产精品久久男人天堂| 国产精品爽爽va在线观看网站 | 午夜免费激情av| 亚洲男人天堂网一区| 亚洲专区中文字幕在线| 国产精品av久久久久免费| 亚洲精品国产一区二区精华液| 精品国产亚洲在线| 色综合站精品国产| 狂野欧美激情性xxxx| 亚洲熟女毛片儿| 亚洲精品美女久久av网站| 别揉我奶头~嗯~啊~动态视频| 窝窝影院91人妻| 婷婷精品国产亚洲av在线| 免费观看人在逋| 亚洲男人的天堂狠狠| 亚洲 欧美 日韩 在线 免费| 亚洲中文av在线| 少妇的丰满在线观看| 婷婷精品国产亚洲av| 欧美黑人欧美精品刺激| 在线观看免费午夜福利视频| 亚洲一区二区三区不卡视频| avwww免费| 禁无遮挡网站| 热re99久久国产66热| 亚洲欧美日韩无卡精品| 久久久久九九精品影院| 色播在线永久视频| 99re在线观看精品视频| a在线观看视频网站| 婷婷亚洲欧美| 午夜亚洲福利在线播放| 久久精品国产综合久久久| 亚洲精品国产精品久久久不卡| 国产人伦9x9x在线观看| 丰满的人妻完整版| 亚洲国产精品久久男人天堂| 男人舔女人下体高潮全视频| 欧美另类亚洲清纯唯美| 一区福利在线观看| 一卡2卡三卡四卡精品乱码亚洲| 村上凉子中文字幕在线| xxxwww97欧美| 久久久久久九九精品二区国产 | 国产视频内射| 大型黄色视频在线免费观看| 女警被强在线播放| 久久精品人妻少妇| 国产熟女午夜一区二区三区| 12—13女人毛片做爰片一| 国产精华一区二区三区| 久久久久久免费高清国产稀缺| 在线观看午夜福利视频| 午夜两性在线视频| 禁无遮挡网站| 国产99白浆流出| 最近最新免费中文字幕在线| 在线观看66精品国产| 色播亚洲综合网| 亚洲国产精品合色在线| 日韩av在线大香蕉| 美女国产高潮福利片在线看| 极品教师在线免费播放| 中文字幕最新亚洲高清| 成人特级黄色片久久久久久久| 天堂√8在线中文| 色尼玛亚洲综合影院| 一级a爱片免费观看的视频| 在线av久久热| 美女大奶头视频| 性欧美人与动物交配| 亚洲精品一区av在线观看| 精品一区二区三区四区五区乱码| 国产av在哪里看| 久久这里只有精品19| 一区二区三区高清视频在线| 精品熟女少妇八av免费久了| 久久久久久九九精品二区国产 | 丁香六月欧美| 哪里可以看免费的av片| 免费观看精品视频网站| av有码第一页| 亚洲午夜理论影院| 午夜亚洲福利在线播放| 搞女人的毛片| 51午夜福利影视在线观看| 一本大道久久a久久精品| 亚洲av电影在线进入| 精品卡一卡二卡四卡免费| 久久久久久大精品| av电影中文网址| 丝袜美腿诱惑在线| 欧美国产精品va在线观看不卡| 搡老妇女老女人老熟妇| 少妇熟女aⅴ在线视频| 日日干狠狠操夜夜爽| www日本在线高清视频| 亚洲欧美一区二区三区黑人| 精品福利观看| 俺也久久电影网| 黑人巨大精品欧美一区二区mp4| 中国美女看黄片| 亚洲五月色婷婷综合| 日韩高清综合在线| 动漫黄色视频在线观看| 国产真实乱freesex| 国产成人欧美| 午夜老司机福利片| 午夜免费成人在线视频| 中出人妻视频一区二区| avwww免费| √禁漫天堂资源中文www| 男人舔女人下体高潮全视频| 伦理电影免费视频| 黄片小视频在线播放| 国产久久久一区二区三区| 国产蜜桃级精品一区二区三区| 久久久久精品国产欧美久久久| 亚洲激情在线av| 两个人看的免费小视频| 欧美激情极品国产一区二区三区| netflix在线观看网站| 国产一区二区三区在线臀色熟女| 亚洲一区二区三区不卡视频| 亚洲aⅴ乱码一区二区在线播放 | 88av欧美| 看片在线看免费视频| 国产精品爽爽va在线观看网站 | 久久久久久国产a免费观看| 人人妻人人澡欧美一区二区| 国产亚洲欧美98| 亚洲国产毛片av蜜桃av| 国产精品 欧美亚洲| 国产av一区二区精品久久| 中文在线观看免费www的网站 | 国产熟女午夜一区二区三区| 波多野结衣高清无吗| 亚洲一区高清亚洲精品| 国产在线观看jvid| 法律面前人人平等表现在哪些方面| 在线永久观看黄色视频| 精品久久久久久久久久免费视频| 色综合站精品国产| 又黄又粗又硬又大视频| 国内毛片毛片毛片毛片毛片| 啦啦啦韩国在线观看视频| 亚洲三区欧美一区| 嫩草影视91久久| a级毛片在线看网站| 精品卡一卡二卡四卡免费| 宅男免费午夜| 首页视频小说图片口味搜索| 国产爱豆传媒在线观看 | 亚洲成人精品中文字幕电影| 后天国语完整版免费观看| 久热爱精品视频在线9| 欧美激情高清一区二区三区| 久久青草综合色| 精品国产美女av久久久久小说| 91麻豆av在线| 国产三级在线视频| 久久香蕉激情| 动漫黄色视频在线观看| a级毛片在线看网站| 亚洲自偷自拍图片 自拍| 国产一区二区三区在线臀色熟女| 可以在线观看毛片的网站| 亚洲男人天堂网一区| 亚洲男人的天堂狠狠| 2021天堂中文幕一二区在线观 | 亚洲一区二区三区色噜噜| 亚洲五月天丁香| 国产伦在线观看视频一区| 国产亚洲欧美在线一区二区| 伦理电影免费视频| 亚洲自拍偷在线| 久99久视频精品免费| 婷婷精品国产亚洲av| 亚洲精品中文字幕一二三四区| 国产亚洲精品av在线| 黄色视频不卡| 黄片播放在线免费| 国产欧美日韩一区二区三| 亚洲男人的天堂狠狠| 在线视频色国产色| 在线看三级毛片| 一区二区三区高清视频在线| 嫁个100分男人电影在线观看| 搞女人的毛片| 男女之事视频高清在线观看| 国产成人精品久久二区二区免费| 视频区欧美日本亚洲| 丁香六月欧美| 久久久水蜜桃国产精品网| 久久久久精品国产欧美久久久| 久久精品夜夜夜夜夜久久蜜豆 | 十八禁人妻一区二区| 免费一级毛片在线播放高清视频| 亚洲熟妇中文字幕五十中出| 色综合欧美亚洲国产小说| 人人妻,人人澡人人爽秒播| 可以在线观看的亚洲视频| 12—13女人毛片做爰片一| 十八禁人妻一区二区| 熟女电影av网| 黄片大片在线免费观看| 亚洲熟女毛片儿| 成年人黄色毛片网站| 麻豆久久精品国产亚洲av| 大香蕉久久成人网| 亚洲精品中文字幕在线视频| 一区二区三区精品91| 在线观看免费视频日本深夜| 日韩大码丰满熟妇| 不卡一级毛片| 99久久无色码亚洲精品果冻| 91在线观看av| 一夜夜www| 在线十欧美十亚洲十日本专区| 精品免费久久久久久久清纯| 日韩国内少妇激情av| 国产亚洲精品久久久久5区| 9191精品国产免费久久| 成人亚洲精品一区在线观看| 精品国产乱子伦一区二区三区| 99精品欧美一区二区三区四区| 夜夜躁狠狠躁天天躁| 亚洲成人国产一区在线观看| 午夜久久久久精精品| 深夜精品福利| 免费观看精品视频网站| 国产欧美日韩精品亚洲av| 亚洲在线自拍视频| 1024香蕉在线观看| √禁漫天堂资源中文www| 男人舔女人下体高潮全视频| 国产v大片淫在线免费观看| 亚洲真实伦在线观看| e午夜精品久久久久久久| 午夜免费成人在线视频| 黄色片一级片一级黄色片| 亚洲全国av大片| 久久午夜亚洲精品久久| 亚洲电影在线观看av| 亚洲aⅴ乱码一区二区在线播放 | 亚洲人成77777在线视频| 国产av一区在线观看免费| 日韩欧美一区二区三区在线观看| 国产单亲对白刺激| 久久国产精品男人的天堂亚洲| 成熟少妇高潮喷水视频| 午夜a级毛片| www.熟女人妻精品国产| 香蕉av资源在线| 久久久久精品国产欧美久久久| 久久天堂一区二区三区四区| 久久国产精品人妻蜜桃| 欧美黑人欧美精品刺激| 久久青草综合色| 视频区欧美日本亚洲| 欧美乱码精品一区二区三区| 黄频高清免费视频| 757午夜福利合集在线观看| 69av精品久久久久久| 成人特级黄色片久久久久久久| 国产视频内射| 欧美黄色片欧美黄色片| 色播在线永久视频| 母亲3免费完整高清在线观看| 国产高清视频在线播放一区| 一卡2卡三卡四卡精品乱码亚洲| 亚洲精品色激情综合| 嫁个100分男人电影在线观看| 亚洲av成人不卡在线观看播放网| 国产精品永久免费网站| 国产精品日韩av在线免费观看| 757午夜福利合集在线观看| 国产色视频综合| 无遮挡黄片免费观看| 精品第一国产精品| 叶爱在线成人免费视频播放| 亚洲中文字幕日韩| 十分钟在线观看高清视频www| 免费av毛片视频| 国产亚洲精品综合一区在线观看 | 日韩高清综合在线| 男女午夜视频在线观看| 精品第一国产精品| 亚洲av第一区精品v没综合| 欧美亚洲日本最大视频资源| 成熟少妇高潮喷水视频| 久久久久久久午夜电影| 亚洲av五月六月丁香网| 别揉我奶头~嗯~啊~动态视频| 国产真人三级小视频在线观看| 啦啦啦观看免费观看视频高清| 久久中文字幕一级| 中文亚洲av片在线观看爽| 十分钟在线观看高清视频www| 欧美激情 高清一区二区三区| 美国免费a级毛片| 黄片大片在线免费观看| 91成人精品电影| 黄片播放在线免费| 成在线人永久免费视频| 一边摸一边做爽爽视频免费| 桃色一区二区三区在线观看| 1024香蕉在线观看| 亚洲专区字幕在线| 男人操女人黄网站| 男女下面进入的视频免费午夜 | 91麻豆精品激情在线观看国产| 精品卡一卡二卡四卡免费| 又大又爽又粗| 色在线成人网| 免费在线观看日本一区| 亚洲国产精品成人综合色| 男人舔女人的私密视频| 欧美国产日韩亚洲一区| 国产精品精品国产色婷婷| 亚洲熟女毛片儿| 精品一区二区三区四区五区乱码| 1024香蕉在线观看| 国产伦一二天堂av在线观看| 精品欧美国产一区二区三| 美女午夜性视频免费| 99精品在免费线老司机午夜| 一个人观看的视频www高清免费观看 | 最近最新中文字幕大全免费视频| 一进一出抽搐gif免费好疼| 国产亚洲欧美精品永久| 亚洲无线在线观看| 少妇 在线观看| 国产日本99.免费观看| 在线免费观看的www视频| 老司机深夜福利视频在线观看| 别揉我奶头~嗯~啊~动态视频| 日韩欧美一区视频在线观看| 色播在线永久视频| 禁无遮挡网站| 精品国产一区二区三区四区第35| 男男h啪啪无遮挡| 一本久久中文字幕| 国产真人三级小视频在线观看|