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

    A comprehensive fractal char combustion model☆

    2016-06-12 03:48:24YutingLiuRongHe
    Chinese Journal of Chemical Engineering 2016年12期

    Yuting Liu *,Rong He

    1 Department of Thermal Engineering,Tsinghua University,Beijing 100084,China

    2 Institute of Materials,China Academy of Engineering Physics,Mianyang 621908,China

    1.Introduction

    The char particles are typical porous media with complex pore structures[1–3],so the combustion reactions can occur at both the external and the internal surfaces of the char particle.What's more,the pore structure in the char particles changes greatly during combustion[3].Many factors affect char combustion including the gas diffusion in the char pores,catalysis or inhibition effects of the ash,graphitization of the char and so on,among which the pore diffusion is the most important[4].The pore diffusion and the chemical reactions inside the pores are coupled with each other[5]which increases the complexities of the char combustion modeling.Since the char particle pore structure determines the pore diffusion resistance[2]and affects the collisions between the gas molecules and the pore surfaces,the pore structure significantly affects the char combustion kinetics[4,6].

    However,the char pore structures are too complex to describe with the classical geometric method[7],so the fractal theory has been used instead[1,7–13].The fractal dimension of the char pore structures,Df,based on mercury intrusion measurements[7]is the most representative one:

    whererporeis the char pore radius(m)andSa(rpore)is the accumulated pore surface area(m2)corresponding to rporein the mercury intrusion measurements.Dfdescribes the complexity of the char pore structures with a largerDfindicating a more complex structure and a higher pore diffusion resistance.The gas diffusion in a bulk space can be described by Fick's diffusion law.However,in small char pores with pore sizes close to the gas molecular mean free path,the gas molecular motion is restricted and the collision frequency between the gas molecules and the pore surface is much higher than the intermolecular collision frequency.In this situation,Fick's diffusion law is not suitable[14–16].In either case,fractal theory can be used to describe the pore diffusion effects.

    The pore diffusion effects have been extensively included in various char combustion models among which the traditional intrinsic model[6,17,18]is the most famous one.The tortuosity factor,τ,is used to describe the complexity of the pore structures in the intrinsic model.However,τ has no scientific definition and cannot be measured or calculated which makes the intrinsic model difficult to use.In the intrinsic model,the pore diffusion is described by the Fick's and Knudsen's diffusion law both of which cannot be conveniently used for the fractal char pores[14–16].In addition,the intrinsic model does not include the evolution of the char pore structure during combustion.

    In this article,a comprehensive fractal char combustion model was developed.The apparent reaction order,n,during combustion was allowed to vary.Four parameters(porosity,θ,specific surface area,Sg,fractal dimension Dfand particle diameter,dp)were used to describe the char pore structures with these parameters also allowed to vary during combustion.The model gives reasonable predictions of the CO/CO2ratio in the combustion products.A correction factor,α,was introduced to describe the effect of the high CO2concentration on the total char combustion rate to expand the scope of the char combustion model.Eleven different chars were combusted in two drop tube furnaces with the conversions of the partly burned char samples measured using thermogravimetric analysis.The combustion processes of these chars were then simulated by the model with the predicted char conversions compared with the measured data to determine the model parameters and verify the model.Finally,the model was used to study some phenomena in the char combustion.

    2.Development of the Comprehensive Char Combustion Model

    The comprehensive char particle combustion model in this article consists of five parts with the char combustion rate expression,the char particle pore structure evolution,the variation of the apparent reaction order,the CO/CO2ratio in the combustion products and the correction factor for oxy-char combustion.The first four parts have already been studied by our research group[3–5,19]and they are combined with the fifth part in this paper to form a complete model.

    2.1.Char combustion rate expression

    Acharcombustion rate expression was developed from a theoretical derivation in[4].The coupling effects between the chemical reactions and the gas diffusion inside the char pores were considered.The final form of the char combustion rate expression is:

    wherekais the apparent reaction rate constant,kiis the intrinsic reaction rate constant,Qis an intermediate variable,Aiis the preexponential factor,Eiis the intrinsic activation energy(J·mol-1),Tpis the particle temperature(K),Ris the universal gas constant,ζ is a scale factor,θ is the porosity,Sgis the specific surface area(m2·kg-1),ρiis the char true density(kg·m-3),dpis the particle diameter(m),Dfis the fractal dimension defined in Eq.(1),χ is a dimensionless modulus which is a combination of pore structure parameters to indicate the pore diffusion resistance with a largerχindicating a lower pore diffusion resistance[3–5].In Eq.(2),Ai,Ei,andζare model parameters to be determined from experimental data.

    2.2.Char particle pore structure evolution during combustion

    Fractal pore models[20]generated by a random walk algorithm can be used to simulate real char particles.Fig.1(a)shows a schematic diagram of a fractal pore model.The pore shapes in the model are irregular as in natural char particles(Fig.1(b))and the model also has good fractal properties.Thus,the fractal pore model is also called the numerical char particle model.The random walk process parameters can be adjusted to generate numericalchar particle models with various porosities and fractal dimensions.

    The complete char combustion processes were simulated[3]with this numerical char particle model with four chemical reactions:

    They were described by simple collision theory[21]with reaction rates expressed as:

    Reaction rate=(Collision frequency)

    Since the Fick's and Knudsen's diffusion law cannot be conveniently used for the fractal char pores[14–16],the diffusion law[22]based on the kinetic theory of gases[23]was used.

    21 different numerical char particle models were used in the char combustion simulations with 756 cases in total for various operating conditions[3]to calculate the variations of the pore structure parameters(θ,Sg,Dfanddp)with the char conversion,X,during combustion.Commonly used empirical expressions with some modifications were then used to associate these pore structure parameters withX.Then,all the coefficients in each expression were related to the operating conditions and the initial pore structure parameters.Finally,the variations of the pore structure parameters during combustion were expressed as:

    Fig.1.Schematic diagram of the fractal pore model(a)and SEM photo of a bituminous char particle(b).

    whereTis the temperature(K),PO2is the O2partial pressure at the particle external surface(Pa),ashis the ash content,Xis the char conversion,a1,a2,b1,b2,b3,c1,c2,d1andd2are constants for a given combustion process and the subscript 0 represents the initial value.

    2.3.Variation of the apparent reaction order during combustion

    The char combustion apparent reaction order,n,has been extensively studied in experiments,but only the values of n for specific chars for specific reaction conditions have been given and no clear universal relationships has been proposed[5,24].This may be the reason why the assumption n=1 is still commonly used in many simulations.

    Recently,an expression fornwas developed for different chars and different reaction conditions from the numerical simulations with the expression validated by experimental data[5].The simulations used the numerical char particle model in Fig.1(a)and the chemical reactions and pore diffusion model described in Section 2.2.

    10 different numerical char particle models were used in the char combustion simulations with 600 cases in total for various operating conditions[5].The values ofnin each case were obtained.First,the relationships betweennand the operating conditions were analyzed for each numerical char particle model to form 10 preliminary expressions.Then,the coefficients in the 10 expressions were correlated by introducing the effect of the particle category.Finally,the expression fornwas written as:

    where Tpis the char particle temperature(K),PO2is the O2partial pressure at the particle external surface(Pa)and χ is the modulus which is the same as in Eq.(2).

    2.4.Apparent CO/CO2 ratio in the combustion products

    The apparent CO/CO2ratio in the combustion products is an important factor affecting the char combustion[19].Since the heat released when CO2is the combustion product is 3.5 times that when CO is the combustion product,a small difference in the apparent CO/CO2ratio will lead to differences of hundreds of degrees in the combustion temperature[25].Usually,the apparent CO/CO2ratio is thought to be a function of the char particle temperature,Tp,and the O2partial pressure at the particle external surface,PO2,based on experimental studies[26].

    The simulations in Section 2.3 also gave the values of the apparent CO/CO2ratio by counting the numbers of CO and CO2molecules diffusing out of the numerical char particle[5].A correlation of the apparent CO/CO2ratio calculated in the simulations gave:

    which is very close to the experimental results in[25].

    2.5.Correction factor for oxy-char combustion

    In oxy-char combustion,CO2is used as the background gas instead of N2.Since the physical properties of CO2and N2are very different and CO2is not an inert gas like N2,the char combustion characteristics in oxy-mode will be a little different.The existing studies[27,28]of oxy-char combustion focused on the effects of the high-concentration CO2in an external environment on the combustion temperature,the flame stability and so on,but the effects of the high-concentration CO2inside the char pores were usually ignored.In the complex char pores,the four chemical reactions in Eqs.(3)–(6)are coupled with each other.If the CO2concentration in the pores is excessive,the original balance will be broken.

    Since the char combustion model in this paper included all the reactions within the char particle pores,a correction factor,α,was defined as:

    whereRoxyis the charcombustion rate(mol·s-1)for the CO2background gas andRairis the char combustion rate(mol·s-1)for the N2background gas.Other operating conditions like the temperature,the O2concentration,and the charparticle model were the same when Roxyand Rairwere solved.In this way,the effect of using CO2instead of N2is reflected in α.

    Then,the char combustion simulations were conducted to calculate α for the different operating conditions.10 different numerical char particle models(as shown in Fig.1(a))were used with the particle diameters and the pore structure parameters listed in Table 1.The chemical reactions and the pore diffusion were described in the same way as in Section 2.2.The simulation conditions are summarized in Table 2.The combustion simulations were conducted for each char particle model for each temperature and each atmosphere listed.Thus,there were 600 cases in total,with one half in air-mode(in an O2/N2atmosphere)and the other half in oxy-mode(in an O2/CO2atmosphere).The total pressure for all cases was 101.325 kPa.The simulations gave the total char combustion rate including the effects of the chemical kinetics and the pore diffusion for the corresponding operating conditions.Then,300 values of α were calculated using Eq.(24).

    The C+CO2reaction rate is generally much lower than the C+O2reaction rate with the measurements giving a C+O2reaction rate that was 5 orders of magnitude higher than for the C+CO2reaction at1073 K and 10 kPa pressure[29].At higher temperatures,the difference was not so extreme,but the C+CO2reaction rate was still much lower than the C+O2reaction rate as long as comparable O2concentrations were present.So,the α obtained in the simulations were very close to 1.

    Table 1Particle diameters and pore structure parameters of the numerical char particle models

    Table 2Simulation conditions for char combustion in air-mode and oxy-mode

    At lower temperatures,α <1,while at higher temperatures,α >1.Also,the demarcation temperatures are different for different char particles.With lower O2partial pressures and higher CO2partial pressures,α differed more from 1.Thus,α was correlated as:

    where PO2and PCO2are the partial pressures of O2and CO2at the char particle external surface(Pa),Tpis the particle temperature(K),andm1andm2are constants for a given char particle.The α calculated in the simulations were used correlated to obtain them1andm2listed in Table 3.

    All the correlation coefficients in Table 3 are higher than 0.98,which shows the validity ofEq.(25).The values ofm1for all the chars are in the range of(6.02±0.017)×10-4with small standard deviations,so m1is set to 6.02×10-4for all the chars.m1seems to be a scale factor.However,m2changes greatly for the different char particles,so it cannot be taken as a constant.The major difference among the char particles is the pore structure which affects the diffusion and reactions in thepores.The modulus χ defined in Eq.(2)is a combination of the pore structure parameters andm2correlates well with the modulus χ with Fig.2 showing the relationship betweenm2and χ which is correlated as:

    Table 3m1 and m2 for each numerical char particle model

    Fig.2.Exponential relationship between m2 and χ.

    Thus,the correction factor,α,has the following form:

    The α calculated by Eq.(27)are compared with those obtained in the simulations in Fig.3.All the points are distributed along the diagonal which shows that the data correlation works well.

    The physical properties of CO2are very different from those of N2.For example,the molar specific heat of CO2is significantly higher than that of N2and the O2diffusion coefficient in CO2is much lower than that in N2.Both factors restrict the char combustion in an O2/CO2atmosphere[27,28].Thus,the char combustion rate in oxy-mode is usually lower than that in air-mode.However,CO2is not inert,especially at high temperatures.As the temperature increases,the C+CO2reaction rate increases,so the char combustion rate in oxy-mode may exceed that in air-mode.

    Fig.3.Data correlation results for α.

    The transition process can be quantitatively described by defining a critical temperature based on Eq.(27)asTcritical=m2=258.45 exp(-59,051.09χ)+1123.3.WhenTp>Tcritical,the char combustion in oxy-mode is enhanced(α>1).WhenTp<Tcritical,the charcombustion in oxy-mode is reduced(α<1).Tcriticalis a characteristic of the char particle which only depends on the particle structure.When χ increases,Tcriticalwill decrease which indicates that the decreasing pore diffusion resistance tends to compensate for the slower oxy-char combustion caused by the lower O2diffusion coefficient in CO2.However,Tcriticalhas a limit of 1123.3 K.Thus,if the particle temperature is below 1123.3 K,the char combustion rate in oxy-mode will be always lower than that in air-mode.

    These conclusions are only for the combustion process inside the char particle,without considering the effect of the bulk diffusion and the reactions in the bulk space.At the same time,only a very short combustion time was simulated since only the char combustion rate for the given operating conditions was needed.Thus,the environmental conditions did not change much in such a short time.As a result,the results do not show the effect of the lower combustion temperature caused by the higher CO2molar specific heat in oxy-mode.However,these factors will not affect the use of Eq.(27)in oxy-char combustion simulations as α only depends on the current conditions and α can be used in conjunction with some process simulation programs.α makes the char combustion rates in oxy-mode more accurate and makes the reaction process closer to the real situation.

    2.6.Model summary

    Sections 2.1–2.5 introduced the five parts of the char particle combustion model.Thus,the combustion rate(mol·s-1)of a char particle,Rs,can be expressed as:

    where α is the correction factor for the oxy-char combustion given by Eq.(27),kais the apparent reaction rate constant described by Eq.(2),mis the char particle mass(kg),Sgis the specific surface area(m2·kg-1),ash is the ash content in the charparticle,Csis the O2concentration at the particle external surface(mol·m-3)andnis the apparent reaction order of the char combustion given by Eq.(22).The variations of all the pore structure parameters in Eqs.(2),(22),(27)and(28)during combustion are described by Eqs.(8)–(21).The apparent CO/CO2ratio in the combustion products is given by Eq.(23).All these model equations together can give a better description of the complex pore diffusion effects on the char combustion.

    3.Drop Tube Furnace Experiments

    Experiments were conducted[4]to determine the unknown model parameters(Ai,Eiand ζ)and to verify the model.

    11 different chars were combusted in two drop tube furnaces(DTFs)with the proximate and ultimate analyses of the 11 parent coals listed in Table 4.These chars were divided into two groups with Nos.1–3 in Group 1 tested in furnace 1 and Nos.4–11 in Group 2 tested in furnace 2.The two furnaces were similar with the same schematic diagram shown as in Fig.4.Both furnaces had 50 mm inner diameters.The length of furnace 1 was 750 mm and the length of furnace 2 was 1000 mm.In furnace 1,the char particle injector was 310 mm from the inlet and the sampling points were distributed along the furnace tube to obtain char samples with different conversions.In furnace 2,the char particle injector was 200 mm from the inlet and the chars were only sampled at the furnace tube outlet.

    Fig.4.Schematic diagram of the drop tube furnace.

    The chars were produced by pyrolysis of the parent coals in a N2atmosphere at around 1473 K in the drop tube furnaces.The experimental conditions for the char combustion are listed in Table 5.The conversions of the char samples,X,were calculated as:

    Table 4Proximate and ultimate analyses of the parent coals

    Table 5Drop tube furnace experimental conditions

    whereash0is the initial ash content of the char andashis the ash content of the partly burned char sample both of which were measured by a Q500 thermogravimetric analyzer produced by TA Instruments,US.The measurement error of this analyzer was only±0.1%of the full scale.

    The average charparticle diameters were measured by a Mastersizer 2000 Laser Particle Size Analyzer.This analyzer could measure particle sizes from 0.02 to 2000 μm and the measurement errors were less than 1%.Repeated measurements show differences of less than 0.5%.The char particle pore structure parameters(θ,SgandDf)were measured by an Autos can 33 mercury porosimeter with an operating pressure range of 0–232 MPa.The smallest pore which this mercury porosimeter could measure was a diameter of 3.23 nm and the nonlinearity was±0.05%of the full scale.The measured average particle diameters and the pore structure parameters of the 11 chars listed in Table 6 were used in the following simulations.The measured pore structure parameters were only for the large pores,while the small pores which nearly had no effect on the char combustion were excluded[1,7].

    Table 6Measured average particle diameters and pore structure parameters of the 11 chars

    4.Model Validation

    Since the flow in the drop tube furnace was a typical gas–solid two phase flow,the Euler–Lagrange approach was used to simulate the drop tube furnace experiments.The air was treated as a continuum described by the Navier–Stokes equations,while the char particles were tracked through the calculated flow field.The char particles exchange mass,momentum and energy with the air in the simulations.Since both the average particle diameters in Table 6 and the feed rates in Table 5 are very small,the interactions among particles and the effects of the particles on the flow field were assumed to be negligible.The detailed equations were given in[4].

    A three-dimensional grid was generated for the cylindrical furnace tube with a total of 222,400 nodes as shown in Fig.5.The air entered the furnace at the gas inlet.The particle injector was inserted into the furnace in the experiments but its influence on the flow field was neglected.Thus,the particles were assumed to be injected from a particle source location in the furnace.

    Fig.5.Grid for the drop tube furnace.

    The combustion rates of each char particle were described by the model introduced in Section 2.The model parameters(Ai,Eiand ζ)given in[4]were initially used as the initial values.The predicted char conversions were then compared with the measured values(only used the experimental results of chars Nos.1–3)with the model parameters(Ai,Eiand ζ)adjusted by trial and error to minimize the differences.The final model parameters for all the chars were determined to be:which completes the char particle combustion model.The value ofEiis the same with that in[4,5],which reconfirms its physical meaning as the intrinsic activation energy.

    Figs.6–8 compare the predicted char conversions with the experimental values for chars Nos.1–3 at all temperatures.The dots represent the measured char conversions and the lines represent the predicted char conversions.The good agreement shows that the comprehensive char particle combustion model is very accurate.Table 7 lists the root mean square errors of the predictions from the experimental results for both the present model and the model of He et al.[4]which shows that the present model gives a much better representation of the data.

    Fig.6.Measured and predicted conversions for char No.1.

    Fig.7.Measured and predicted conversions for char No.2.

    Fig.8.Measured and predicted conversions for char No.3.

    Only the final conversions were measured for chars Nos.4–11 in the experiments since furnace 2 only had one sampling point at the furnace tube outlet.The measured final conversions are listed in Table 8 together with the ash contents of the char samples.Then,the final char conversions were predicted by both the intrinsic model[6,17,18]and the present model with the results also listed in Table 8.All the predicted final conversions are plotted in Fig.9 versus the measured final conversions to show that the present model is more accurate than the intrinsic model.

    5.Discussion

    The combustion of a single char particle in a bulk space was simulated with the present model to get more details about the charcombustion process.The char particle was assumed to be spherical and isothermal.The equations describing the char particle were the same with those in the drop tube furnace experiment simulations[4].The one-dimensional conservation equations for multicomponent reacting systems in spherical coordinates[30]were used to describe the heat and mass transfer processes and the chemical reactions in the flow field around the char particle.The effects of the char particle combustion on the flow field were used as the boundary conditions.The No.2 char in Tables 4 and 6 was used with the detailed simulation conditions listed in Table 9.In all the cases,the total pressure is 101.325 kPa and the initial char particle temperature is 300 K.

    Table 7Root mean square errors for the model of He et al.[4]and the present model relative to the experimental data

    Table 8Comparisons of the predicted final conversions with the measured data for chars in Group 2

    Fig.9.Comparison of the present model with the traditional intrinsic model.

    Figs.10–15 show the simulation results in which “Air”means that the case used the air-mode(in an O2/N2atmosphere)and “Oxy”means the case used the oxy-mode(in an O2/CO2atmosphere).The results in Figs.10–14 were obtained with an initial O2mole fraction at the char particle external surface of 21%while the results in Fig.15 were obtained with an initial temperature of 1200 K.The other simulation conditions are indicated in the figures.

    Table 9Simulation conditions for a single char particle combustion

    Fig.10.Variations of the char particle temperature during combustion(21%O2).

    Fig.11.Variations of the total char combustion rate for different temperatures during combustion(21%O2).

    Fig.10 shows the variations of the char particle temperature,Tp,during combustion.The differences inTpbetween the air-mode and the oxy-mode are very small andTpin air-mode is slightly higher than that in oxy-mode.Fig.11 shows the variations of the total char combustion rate during combustion.The combustion rates in air-mode are higher than those in oxy-mode for all temperatures.These phenomena are mainly due to the higher molarspecific heat of CO2and the lower O2diffusion coefficient in CO2which slow the char combustion in oxymode[27,28].Even though the char can react with CO2,the contribution is insignificant since the reaction rate of char with CO2is much lower than that of char with O2as long as comparable O2concentrations are present[29].Since the complete char combustion processes were simulated with the bulk diffusion effects included,the reaction rate differences between the air-mode and the oxy-mode in Fig.11 are higher than those in Fig.3.

    Fig.12 shows the variations of the gas mole fractions at the char particle external surface during combustion.Since O2is a reactant in the combustion reactions,the O2mole fractions decrease as the combustion proceeds.During the final stage of the combustion,the combustion rates decrease significantly and the amount of O2diffusing from the ambient to the charparticle surface is greater than the amount ofO2consumed in the combustion reactions.Thus,the O2mole fractions gradually increase.At lower temperatures,the O2mole fractions decrease more in the air-mode than in the oxy-mode since the combustion rates in the air-mode are faster.However,at higher temperatures,although the combustion rates in air-mode are still faster,the O2mole fractions decrease more in the oxy-mode.This is because more CO is generated in the oxy-mode at higher temperatures as can be seen in Fig.12 and the combustion of CO will consume more O2.The char-CO2gasification reaction rate is faster at higher temperatures,especially in oxy-mode with a large amount of CO2in the atmosphere.What's more,the apparent CO/CO2ratio is also higher at higher temperatures.

    The variations of the CO mole fractions in Fig.12 are the opposite of the O2mole fractions with the reasons being analogous to those for the variations of the O2mole fractions.In air-mode,there is no CO2in the atmosphere before combustion,so the CO2mole fractions first increase and then decrease during combustion,which is similar to the CO mole fractions.The CO mole fractions are higher than the CO2mole fractions especially at higher temperatures,which shows that CO is the main product of the char combustion.In oxy-mode,the CO2is used as the background gas with relatively high mole fractions before combustion,so the CO2mole fractions first decrease as the reactant.Then,the CO2mole fractions increase significantly and even exceed the original level,because CO2is also one of the combustion products.

    Fig.13 shows the variations of the char particle specific surface area,Sg,during combustion.Since the char particle pore structure evolution model in Section 2.2 has the same form for both the air-mode and the oxy-mode,only the relationships betweenSgand the char conversion,X,in air-mode are shown.The results show that the temporal variations ofSgdiffer for the two modes.Sgfirst increases to a maximum and then decreases with increasingXwhich is the same as in some experiments[31–38].The increase ofSgis caused by the pore growth and the opening of closed pores and the decrease ofSgis caused by the collapse of the pores and the coalescence of the pore walls[3].The variations ofSgare similar to the variations of the total char combustion rates in Fig.11.The porosity,θ,increases(not linearly)while the fractal dimension,Df,and the particle diameter,dp,decrease(not linearly)during combustion,which is the same as in[3].The decrease ofDfmeans the char pore structures become less complex after combustion[3].

    Fig.14 shows the variations of the apparent reaction order,n,during combustion.During the initial stage of combustion,ndecreases rapidly from 1.Then,nbecomes stable.During the final stage of combustion,nincreases gradually.As the O2concentration at the particle surface changes,the pore diffusion resistance and the time constant for the flow and diffusion into the pores postpone the variations of the O2concentration at the reactive surface inside the char pores andnis not always 1[5].In the simulations,the initial char particle temperature is very low(300 K)and the pore diffusion is not the key factor restricting the chemical reactions,sonis approximately equal to 1.As the char particle heats and burns,the particle temperature rises quickly and the effect of pore diffusion is significantly enhanced,sondecreases rapidly.As the combustion reaches steady state,nbecomes stable.During the final stage of combustion,most of the carbon has been consumed and the pore structure becomes less complex[3]which means that the pore structure provides less pore diffusion resistance,songradually increases.Liu and He[5]also studied the variations ofnduring combustion and found that n decreased first and then became stable towards the end of the combustion process without increasing because the char pore structure evolution was not included in their simulations.

    Fig.15 shows the variations of the total char combustion rate for different O2concentrations at a higher temperature of 1200 K.Overall,the total char combustion rates in the air-mode are still higher than those in the oxy-mode and the differences become larger at higher O2concentrations.However,with an initial O2mole fraction of 5%,the total char combustion rates in the two modes are nearly the same.The total char combustion rates are actually a little higher in oxy-mode at this time but the difference cannot be seen in Fig.15.This shows that the effect of the char-CO2gasification reaction is important only when the temperature is relatively high and the O2concentration is relatively low.

    Fig.12.Variations of the gas mole fractions at the char particle external surface during combustion(21%O2).

    6.Conclusions

    The char combustion mechanisms were analyzed in terms of(1)the pore diffusion effects on the char combustion,(2)the char particle pore structure evolution during combustion,(3)the variation of the apparent reaction order during combustion,(4)the apparent CO/CO2ratio in the combustion products,and(5)the effects of high CO2concentration on the char combustion.A comprehensive char particle combustion model was then developed with all these mechanisms.Nearly all the complex factors affecting the char combustion are considered except for the catalysis of the trace elements in the ash and the graphitization of the chars which will be studied in a future work.

    Fig.13.Variations of the char particle specific surface area during combustion(21%O2).

    Eleven different chars were used in combustion experiments in two drop tube furnaces.The conversions of the partly burned char samples were measured using thermogravimetric analysis.The predicted char conversions matched very well with the measured values which showed the reliability of the present model.What's more,this model is significantly more accurate than the model of Heet al.[4]and the traditional intrinsic model.

    Fig.14.Variations of the apparent reaction order during combustion(21%O2).

    Fig.15.Variations of the total char combustion rate for different O2 concentrations during combustion(T0=1200 K).

    The apparent reaction order first decreases,then stabilizes and finally increases during the combustion process.The specific surface area of the char particle first increases and then decreases with a maximum during combustion.The combustion rates in oxy-mode are generally slower than those in air-mode due to the higher CO2molar specific heat and the lower O2diffusion coefficient in CO2which slow the char combustion in oxy-mode.The effect of the char-CO2gasification reaction becomes important only when the temperature is relatively high and the O2concentration is relatively low.

    Nomenclature

    Aipre-exponential factor

    a1,a2constants in Eq.(8)

    ash char particle ash content

    b1,b2,b3constants in Eq.(11)

    CsO2molar concentration at the char particle external surface,mol·m-3

    c1,c2constants in Eq.(15)

    Dffractal dimension of the char pores

    d1,d2constants in Eq.(18)

    dpchar particle diameter,m

    Eiintrinsic activation energy,J·mol-1

    kaapparent reaction rate constant

    kiintrinsic reaction rate constant

    mchar particle mass,kg

    m1,m2constants in Eq.(25)

    napparent reaction order of char combustion

    PCO2CO2partial pressure at the char particle external surface,Pa

    PO2O2partial pressure at the char particle external surface,Pa

    Qintermediate variable

    R universal gas constant(R=8.314 J·mol-1·K-1)

    Rairchar particle combustion rate in the O2/N2atmosphere,mol·s-1

    Roxychar particle combustion rate in the O2/CO2atmosphere,mol·s-1

    Rschar particle combustion rate,mol·s-1

    rporechar pore radius,m

    Saaccumulated pore surface area in mercury intrusion measurements,m2

    Sgspecific surface area of the char particle,m2·kg-1

    Ttemperature,K

    Tcriticalcritical temperature defined based on Eq.(27),K

    Tpchar particle temperature,K

    X char conversion

    α correction factor for oxy-char combustion defined in Eq.(24)

    ΔHreaction heat,kJ·mol-1

    ζ scale factor

    θ char particle porosity

    ρitrue density of the char particle,kg·m-3

    τ tortuosity factor of the char pores

    χ dimensionless modulus defined in Eq.(2)

    Subscript

    0 initial value of the variable

    [1]R.He,J.I.Sato,C.H.Chen,Modeling char combustion with fractal pore effects,Combust.Sci.Technol.174(2002)19–37.

    [2]W.He,R.He,L.Cao,T.Ito,T.Suda,J.I.Sato,Numerical study of the relationships between pore structures and reaction parameters for coal char particles,Combust.Sci.Technol.184(2012)2084–2099.

    [3]Y.Liu,R.He,Modeling of the pore structure evolution in porous char particles during combustion,Combust.Sci.Technol.188(2016)207–232.

    [4]W.He,Y.Liu,R.He,T.Ito,T.Suda,T.Fujimori,H.Ikeda,J.I.Sato,Combustion rate for char with fractal pore characteristics,Combust.Sci.Technol.185(2013)1624–1643.

    [5]Y.Liu,R.He,Variation of apparent reaction order in char combustion and its effect on a fractal char combustion model,Combust.Sci.Technol.187(2015)1638–1660.

    [6]I.W.Smith,The combustion rates of coal chars:A review,Symp.Combust.19(1982)1045–1065.

    [7]R.He,X.Xu,C.Chen,H.Fan,B.Zhang,Evolution of pore fractal dimensions for burning porous chars,Fuel77(1998)1291–1295.

    [8]D.Avnir,D.Farin,P.Pfeifer,Surface geometric irregularity of particulate materials:The fractal approach,J.Colloid Interface Sci.103(1985)112–123.

    [9]C.Fairbridge,S.H.Ng,A.D.Palmer,Fractal analysis of gas-adsorption on syncrude coke,Fuel65(1986)1759–1762.

    [10]W.I.Friesen,O.I.Ogunsola,Mercury porosimetry of upgraded western Canadian coals,Fuel74(1995)604–609.

    [11]P.J.Mcmahon,I.K.Snook,S.D.Moss,P.R.Johnston,Influence of fractal pores on the oxidation behavior of brown coal,Energy Fuel13(1999)965–968.

    [12]P.Pfeifer,D.Avnir,Chemistry in noninteger dimensions between two and three:I.Fractal theory of heterogeneous surfaces,J.Chem.Phys.79(1983)3558–3565.

    [13]P.Salatino,F.Zimbardi,S.Masi,A fractal approach to the analysis of lowtemperature combustion-rate of a coal char:I.Experimental results,Carbon31(1993)501–508.

    [14]M.Costa,A.D.Araujo,H.F.Da Silva,J.S.Andrade,Scaling behavior of diffusion and reaction processes in percolating porous media,Phys.Rev.E67(2003)061408–061411.

    [15]Y.Gefen,A.Aharony,S.Alexander,Anomalous diffusion on percolating clusters,Phys.Rev.Lett.50(1983)77–80.

    [16]P.Levitz,From Knudsen diffusion to Levy walks,Europhys.Lett.39(1997)593–598.

    [17]D.Fortsch,R.H.Essenhigh,U.Schnell,K.R.G.Hein,On the application of the Thiele/Zeldovich analysis to porous carbon combustion,Energy Fuel17(2003)901–906.

    [18]L.Ma,R.Mitchell,Modeling char oxidation behavior under Zone II burning conditions at elevated pressures,Combust.Flame156(2009)37–50.

    [19]W.He,R.He,T.Ito,T.Suda,J.I.Sato,Numerical investigations of CO/CO2ratio in char combustion,Combust.Sci.Technol.183(2011)868–882.

    [20]Z.Liang,R.He,Q.Chen,X.Xu,J.I.Sato,Fractal generation of char pores through random walk,Combust.Sci.Technol.179(2007)637–661.

    [21]J.E.House,Principles of chemical kinetics,second ed.Academic Press,Pittsburgh,2007.

    [22]L.Cao,R.He,Gas diffusion in fractal porous media,Combust.Sci.Technol.182(2010)822–841.

    [23]T.I.Gombosi,Gas kinetic theory,Cambridge University Press,Cambridge,1994.

    [24]R.H.Hurt,J.M.Calo,Semi-global intrinsic kinetics for char combustion modeling,Combust.Flame125(2001)1138–1149.

    [25]L.Tognotti,J.P.Longwell,A.F.Sarofim,The products of the high temperature oxidation of a single char particle in an electrodynamic balance,Symp.Combust.23(1991)1207–1213.

    [26]A.N.Hayhurst,M.S.Parmar,Does solid carbon burn in oxygen to give the gaseous intermediate CO or produce CO2directly?Some experiments in a hot bed of sand fluidized by air,Chem.Eng.Sci.53(1998)427–438.

    [27]L.Chen,S.Z.Yong,A.F.Ghoniem,Oxy-fuel combustion of pulverized coal:characterization,fundamentals,stabilization and CFD modeling,Prog.Energy Combust.Sci.38(2012)156–214.

    [28]T.Wall,Y.Liu,C.Spero,L.Elliott,S.Khare,R.Rathnam,F.Zeenathal,B.Moghtaderi,B.Buhre,C.Sheng,R.Gupta,T.Yamada,K.Makino,J.Yu,An overview on oxyfuel coal combustion—State of the art research and technology development,Chem.Eng.Res.Des.87(2009)1003–1016.

    [29]P.L.Walker Jr.,F.Rusinko Jr.,L.G.Austin,Gas reactions of carbon,Adv.Catal.11(1959)133–221.

    [30]K.K.Kuo,Principles of combustion,second ed.John Wiley and Sons,New Jersey,2005.

    [31]I.Aarna,E.M.Suuberg,Changes in reactive surface area and porosity during char oxidation,Symp.Combust.27(1998)2933–2939.

    [32]K.E.Adams,D.R.Glasson,S.A.A.Jayaweera,Development of porosity during the combustion of coals and cokes,Carbon27(1989)95–101.

    [33]B.Feng,S.K.Bhatia,Variation of the pore structure of coal chars during gasification,Carbon41(2003)507–523.

    [34]H.Lorenz,E.Carrea,M.Tamura,J.Haas,The role of char surface structure development in pulverized fuel combustion,Fuel79(2000)1161–1172.

    [35]L.M.Lu,V.Sahajwalla,D.Harris,Coal char reactivity and structural evolution during combustion—Factors influencing blast furnace pulverized coal injection operation,Metall.Mater.Trans.B Process Metall.Mater.Process.Sci.32(2001)811–820.

    [36]A.K.Sadhukhan,P.Gupta,R.K.Saha,Characterization of porous structure of coalchar from a single devolatilized coal particle:coal combustion in a fluidized bed,Fuel Process.Technol.90(2009)692–700.

    [37]I.Sircar,A.Sane,W.Wang,J.P.Gore,Experimental and modeling study of pinewood char gasification with CO2,Fuel119(2014)38–46.

    [38]J.L.Su,D.D.Perlmutter,Effect of pore structure on char oxidation kinetics,AIChE J.31(1985)973–981.

    女人被狂操c到高潮| 成年版毛片免费区| 亚洲欧美日韩高清在线视频| 国产精品免费大片| 亚洲成国产人片在线观看| 精品无人区乱码1区二区| 欧美日韩一级在线毛片| 波多野结衣一区麻豆| 久久天躁狠狠躁夜夜2o2o| 黑丝袜美女国产一区| 午夜福利,免费看| 久久久精品免费免费高清| 午夜福利乱码中文字幕| 在线观看一区二区三区激情| 欧美+亚洲+日韩+国产| 婷婷成人精品国产| 亚洲精品一卡2卡三卡4卡5卡| 97人妻天天添夜夜摸| 亚洲一区二区三区不卡视频| 国产免费男女视频| 一区二区三区国产精品乱码| 一区二区三区国产精品乱码| 国产1区2区3区精品| 高潮久久久久久久久久久不卡| 久久人妻av系列| 国产人伦9x9x在线观看| 精品第一国产精品| 精品久久久久久久久久免费视频 | 大陆偷拍与自拍| 免费在线观看完整版高清| 亚洲精品美女久久av网站| 男女床上黄色一级片免费看| 欧美午夜高清在线| 国产精品久久久av美女十八| 色在线成人网| 伊人久久大香线蕉亚洲五| 99国产极品粉嫩在线观看| 下体分泌物呈黄色| 啦啦啦在线免费观看视频4| 男人操女人黄网站| 亚洲色图 男人天堂 中文字幕| 欧美在线一区亚洲| 亚洲av日韩在线播放| 啪啪无遮挡十八禁网站| 这个男人来自地球电影免费观看| 一级片免费观看大全| 中国美女看黄片| 十八禁人妻一区二区| 日韩大码丰满熟妇| 大香蕉久久成人网| 久9热在线精品视频| 免费观看a级毛片全部| x7x7x7水蜜桃| 久久久久视频综合| 国产精品亚洲av一区麻豆| 成人国语在线视频| 亚洲成国产人片在线观看| 在线av久久热| 十分钟在线观看高清视频www| 老司机午夜福利在线观看视频| 人人妻人人澡人人看| 日韩一卡2卡3卡4卡2021年| 麻豆国产av国片精品| 欧美精品一区二区免费开放| 成人黄色视频免费在线看| 免费观看a级毛片全部| 国产国语露脸激情在线看| 99久久国产精品久久久| 久久久精品国产亚洲av高清涩受| 夜夜夜夜夜久久久久| 欧美国产精品一级二级三级| 我的亚洲天堂| 国产男女超爽视频在线观看| 国产一区在线观看成人免费| 国产真人三级小视频在线观看| av在线播放免费不卡| 国产欧美日韩一区二区三| 91麻豆精品激情在线观看国产 | 99精品欧美一区二区三区四区| 欧美激情极品国产一区二区三区| 亚洲成av片中文字幕在线观看| 国产主播在线观看一区二区| 国产91精品成人一区二区三区| 日本黄色日本黄色录像| 久久中文字幕人妻熟女| 天天添夜夜摸| 一夜夜www| cao死你这个sao货| 人成视频在线观看免费观看| 中文欧美无线码| 热99re8久久精品国产| 自拍欧美九色日韩亚洲蝌蚪91| 国产高清videossex| 成人影院久久| 国产99久久九九免费精品| 美女福利国产在线| 亚洲一区中文字幕在线| 久久国产亚洲av麻豆专区| 精品久久久精品久久久| 欧美日韩精品网址| 欧美日韩黄片免| 999精品在线视频| 国产欧美日韩一区二区三区在线| 十八禁高潮呻吟视频| 国产一区在线观看成人免费| 亚洲第一青青草原| 成在线人永久免费视频| 啪啪无遮挡十八禁网站| 亚洲五月婷婷丁香| 无人区码免费观看不卡| 中文字幕色久视频| av一本久久久久| 免费在线观看日本一区| 天天躁日日躁夜夜躁夜夜| 97人妻天天添夜夜摸| 999久久久国产精品视频| 天天躁狠狠躁夜夜躁狠狠躁| 精品乱码久久久久久99久播| 久久精品人人爽人人爽视色| 性色av乱码一区二区三区2| 成年人黄色毛片网站| 亚洲第一欧美日韩一区二区三区| 久久热在线av| 91九色精品人成在线观看| 国产精品一区二区精品视频观看| 国产精品亚洲一级av第二区| 精品电影一区二区在线| 日韩成人在线观看一区二区三区| 成年动漫av网址| 免费在线观看完整版高清| 美女高潮到喷水免费观看| 男女免费视频国产| 香蕉国产在线看| 成人三级做爰电影| 日日夜夜操网爽| 在线十欧美十亚洲十日本专区| 黄网站色视频无遮挡免费观看| 精品电影一区二区在线| 久久狼人影院| 欧美日韩国产mv在线观看视频| 欧美日韩福利视频一区二区| 一进一出抽搐动态| 久久精品亚洲av国产电影网| 99久久综合精品五月天人人| 少妇猛男粗大的猛烈进出视频| 国产乱人伦免费视频| 久久久久久人人人人人| 80岁老熟妇乱子伦牲交| 国产免费现黄频在线看| 1024视频免费在线观看| 久久久国产成人精品二区 | 无人区码免费观看不卡| svipshipincom国产片| 亚洲精品av麻豆狂野| 在线十欧美十亚洲十日本专区| 欧美亚洲日本最大视频资源| 国产真人三级小视频在线观看| 欧美激情极品国产一区二区三区| 久久天躁狠狠躁夜夜2o2o| 夜夜夜夜夜久久久久| 国产精品电影一区二区三区 | 91av网站免费观看| 久久午夜亚洲精品久久| 国产男女内射视频| 色综合婷婷激情| 在线观看www视频免费| 日本欧美视频一区| 亚洲久久久国产精品| 久久天躁狠狠躁夜夜2o2o| 成人国语在线视频| 国产成人免费无遮挡视频| 在线观看www视频免费| 久久久国产欧美日韩av| av欧美777| 在线天堂中文资源库| 高清黄色对白视频在线免费看| 成人黄色视频免费在线看| 美女扒开内裤让男人捅视频| 啦啦啦在线免费观看视频4| 下体分泌物呈黄色| 日韩成人在线观看一区二区三区| 1024香蕉在线观看| 91成年电影在线观看| 亚洲av日韩精品久久久久久密| 18在线观看网站| 成人国语在线视频| 日日爽夜夜爽网站| 麻豆乱淫一区二区| 精品亚洲成a人片在线观看| 精品乱码久久久久久99久播| 久久人妻熟女aⅴ| 免费在线观看黄色视频的| 国内久久婷婷六月综合欲色啪| 精品卡一卡二卡四卡免费| 女同久久另类99精品国产91| 婷婷精品国产亚洲av在线 | 成人18禁在线播放| 亚洲欧美一区二区三区久久| 在线观看www视频免费| 久久亚洲精品不卡| 国产精品一区二区在线观看99| 久久久久国内视频| 中文字幕色久视频| 国产成人免费无遮挡视频| 亚洲精品成人av观看孕妇| 村上凉子中文字幕在线| 国产精品二区激情视频| 欧美精品人与动牲交sv欧美| 99久久精品国产亚洲精品| 国产在线精品亚洲第一网站| 亚洲va日本ⅴa欧美va伊人久久| 老汉色∧v一级毛片| 他把我摸到了高潮在线观看| 一夜夜www| av片东京热男人的天堂| 国产99久久九九免费精品| 欧美激情久久久久久爽电影 | av有码第一页| 日韩大码丰满熟妇| av网站在线播放免费| 亚洲av成人一区二区三| 国产亚洲av高清不卡| 99re在线观看精品视频| 国产av一区二区精品久久| 日韩欧美一区二区三区在线观看 | 久久性视频一级片| 人妻久久中文字幕网| 欧美 亚洲 国产 日韩一| 欧美亚洲日本最大视频资源| 一区福利在线观看| videosex国产| 欧美日韩视频精品一区| 免费av中文字幕在线| 国产在视频线精品| 大型黄色视频在线免费观看| 黄片小视频在线播放| 日本a在线网址| 精品一区二区三区四区五区乱码| 看免费av毛片| 亚洲欧美激情在线| 两个人免费观看高清视频| 中文字幕另类日韩欧美亚洲嫩草| 视频区图区小说| 精品国内亚洲2022精品成人 | 亚洲视频免费观看视频| 亚洲七黄色美女视频| 18禁裸乳无遮挡动漫免费视频| 看片在线看免费视频| 精品一区二区三区视频在线观看免费 | 久久精品91无色码中文字幕| 欧美黄色淫秽网站| 高潮久久久久久久久久久不卡| 午夜两性在线视频| 久久精品国产亚洲av香蕉五月 | 女人爽到高潮嗷嗷叫在线视频| 午夜精品国产一区二区电影| 欧美黑人欧美精品刺激| 亚洲精品自拍成人| 色综合欧美亚洲国产小说| 一本大道久久a久久精品| 欧美av亚洲av综合av国产av| 免费高清在线观看日韩| 人人妻人人澡人人看| 91精品国产国语对白视频| 国产99白浆流出| 亚洲,欧美精品.| 欧美日韩成人在线一区二区| 欧美精品av麻豆av| 亚洲欧美日韩另类电影网站| 一区二区日韩欧美中文字幕| 69av精品久久久久久| 大片电影免费在线观看免费| 成人永久免费在线观看视频| 久久久国产欧美日韩av| 成在线人永久免费视频| 十八禁人妻一区二区| 亚洲精品国产色婷婷电影| 少妇 在线观看| 久久精品国产亚洲av香蕉五月 | 久久精品亚洲熟妇少妇任你| 久久久国产欧美日韩av| 欧美激情高清一区二区三区| 免费在线观看亚洲国产| 热re99久久国产66热| 精品人妻1区二区| 大陆偷拍与自拍| 欧美丝袜亚洲另类 | 亚洲精品av麻豆狂野| 不卡一级毛片| 一二三四社区在线视频社区8| 中文字幕人妻熟女乱码| 免费在线观看日本一区| 999久久久国产精品视频| 精品国产乱码久久久久久男人| 欧美成人免费av一区二区三区 | 少妇被粗大的猛进出69影院| 1024视频免费在线观看| 丰满迷人的少妇在线观看| 欧美日韩黄片免| 成人免费观看视频高清| 亚洲成a人片在线一区二区| 久久精品亚洲精品国产色婷小说| 一级a爱片免费观看的视频| 日本a在线网址| 亚洲精品国产色婷婷电影| 不卡一级毛片| 亚洲av日韩精品久久久久久密| 99热网站在线观看| 欧美日韩av久久| 午夜福利,免费看| 99精品久久久久人妻精品| 亚洲欧美一区二区三区久久| 美女午夜性视频免费| 日韩大码丰满熟妇| 精品少妇久久久久久888优播| 亚洲专区中文字幕在线| 日韩视频一区二区在线观看| 亚洲精品乱久久久久久| 亚洲国产精品合色在线| 一进一出好大好爽视频| 欧美精品一区二区免费开放| 男男h啪啪无遮挡| 国产日韩欧美亚洲二区| 高清黄色对白视频在线免费看| av天堂久久9| 欧美 日韩 精品 国产| 亚洲性夜色夜夜综合| 国产精品永久免费网站| 黄色 视频免费看| 久久人妻福利社区极品人妻图片| 久久久久国内视频| 脱女人内裤的视频| 国产精品久久久av美女十八| 亚洲五月婷婷丁香| 亚洲国产欧美网| 午夜福利一区二区在线看| 精品熟女少妇八av免费久了| 亚洲自偷自拍图片 自拍| 涩涩av久久男人的天堂| 久久精品熟女亚洲av麻豆精品| www.自偷自拍.com| 亚洲欧美激情在线| 老熟妇仑乱视频hdxx| 老汉色∧v一级毛片| 在线观看免费视频日本深夜| a级毛片黄视频| 亚洲va日本ⅴa欧美va伊人久久| 亚洲全国av大片| 精品国产乱码久久久久久男人| 国产精品美女特级片免费视频播放器 | 黑人猛操日本美女一级片| 亚洲成人免费电影在线观看| 高清视频免费观看一区二区| 免费不卡黄色视频| 欧美精品啪啪一区二区三区| 午夜免费观看网址| 亚洲人成电影免费在线| 国产免费男女视频| 少妇被粗大的猛进出69影院| 1024视频免费在线观看| 丁香六月欧美| 精品国产一区二区三区四区第35| 欧美日韩成人在线一区二区| 在线国产一区二区在线| 在线观看一区二区三区激情| 久久久久国内视频| 窝窝影院91人妻| 99re在线观看精品视频| 久久午夜亚洲精品久久| 在线观看舔阴道视频| 欧美黑人精品巨大| 久久九九热精品免费| 又黄又粗又硬又大视频| 他把我摸到了高潮在线观看| 国产欧美日韩综合在线一区二区| 亚洲va日本ⅴa欧美va伊人久久| 亚洲av美国av| 亚洲免费av在线视频| 亚洲精品成人av观看孕妇| 美女高潮到喷水免费观看| 亚洲成国产人片在线观看| 久久国产精品影院| 大香蕉久久网| 国产成人免费观看mmmm| 欧美日韩精品网址| 国产欧美亚洲国产| 99热网站在线观看| 一二三四在线观看免费中文在| 99国产精品免费福利视频| 性少妇av在线| 麻豆成人av在线观看| 国产一区二区三区综合在线观看| 国产单亲对白刺激| 黄色成人免费大全| 大香蕉久久网| 亚洲精品国产一区二区精华液| 一进一出抽搐gif免费好疼 | 亚洲人成77777在线视频| 日本vs欧美在线观看视频| 水蜜桃什么品种好| 在线天堂中文资源库| 美女扒开内裤让男人捅视频| 999精品在线视频| 狠狠狠狠99中文字幕| 国产aⅴ精品一区二区三区波| 韩国精品一区二区三区| 欧美精品人与动牲交sv欧美| 丁香六月欧美| 国产免费av片在线观看野外av| 精品熟女少妇八av免费久了| a级毛片黄视频| 欧美日韩中文字幕国产精品一区二区三区 | 国产高清国产精品国产三级| 欧美+亚洲+日韩+国产| 久久天躁狠狠躁夜夜2o2o| 国产成人精品久久二区二区免费| 国产免费av片在线观看野外av| 久久久精品国产亚洲av高清涩受| a级毛片黄视频| 亚洲熟妇熟女久久| 成熟少妇高潮喷水视频| 亚洲全国av大片| 亚洲av日韩在线播放| 啦啦啦视频在线资源免费观看| 亚洲第一青青草原| 国产区一区二久久| 变态另类成人亚洲欧美熟女 | 王馨瑶露胸无遮挡在线观看| 日韩欧美一区二区三区在线观看 | 777米奇影视久久| 久久久久久免费高清国产稀缺| 狠狠婷婷综合久久久久久88av| 丰满迷人的少妇在线观看| 亚洲精品中文字幕在线视频| 男女免费视频国产| 日日摸夜夜添夜夜添小说| 一区二区日韩欧美中文字幕| 欧美人与性动交α欧美精品济南到| 久久久国产一区二区| 国产在线一区二区三区精| 制服诱惑二区| 国产精品1区2区在线观看. | 极品人妻少妇av视频| 日本黄色视频三级网站网址 | 国产三级黄色录像| av一本久久久久| 在线观看一区二区三区激情| 精品熟女少妇八av免费久了| 97人妻天天添夜夜摸| 成人亚洲精品一区在线观看| 成人18禁高潮啪啪吃奶动态图| 亚洲国产精品sss在线观看 | 黑丝袜美女国产一区| 搡老乐熟女国产| 男女午夜视频在线观看| 制服诱惑二区| 成年版毛片免费区| 日韩视频一区二区在线观看| 久久人人爽av亚洲精品天堂| 成年动漫av网址| 无人区码免费观看不卡| 丰满饥渴人妻一区二区三| 亚洲国产精品一区二区三区在线| 黄色毛片三级朝国网站| 久久精品国产清高在天天线| 99久久综合精品五月天人人| 亚洲全国av大片| 亚洲av日韩在线播放| 男男h啪啪无遮挡| 欧美黄色片欧美黄色片| av网站在线播放免费| av天堂在线播放| 亚洲色图 男人天堂 中文字幕| 午夜激情av网站| 国产有黄有色有爽视频| 一级毛片精品| 久久精品国产清高在天天线| 国产成人免费观看mmmm| 亚洲精品国产区一区二| 久久天躁狠狠躁夜夜2o2o| 欧美 亚洲 国产 日韩一| 视频区欧美日本亚洲| 91老司机精品| 日韩欧美一区二区三区在线观看 | 无限看片的www在线观看| 人人妻人人澡人人爽人人夜夜| 久久国产精品人妻蜜桃| 久久狼人影院| av片东京热男人的天堂| 18禁观看日本| a在线观看视频网站| 欧美日韩av久久| 少妇裸体淫交视频免费看高清 | 成人亚洲精品一区在线观看| 怎么达到女性高潮| 在线免费观看的www视频| 国产高清videossex| 国产又色又爽无遮挡免费看| 国产亚洲欧美98| 99国产极品粉嫩在线观看| tocl精华| 国产无遮挡羞羞视频在线观看| 嫁个100分男人电影在线观看| 少妇 在线观看| 在线观看午夜福利视频| 欧美日韩av久久| 99国产精品一区二区蜜桃av | 精品免费久久久久久久清纯 | 国产亚洲精品久久久久5区| 韩国av一区二区三区四区| 一个人免费在线观看的高清视频| 久久国产亚洲av麻豆专区| 亚洲精品美女久久久久99蜜臀| 久久国产精品人妻蜜桃| 美女高潮喷水抽搐中文字幕| www.熟女人妻精品国产| 国产精品香港三级国产av潘金莲| 老司机午夜福利在线观看视频| 国产精品久久久av美女十八| 首页视频小说图片口味搜索| 免费黄频网站在线观看国产| 国产区一区二久久| 色在线成人网| 久久久久视频综合| 一级a爱片免费观看的视频| 久久草成人影院| 美女午夜性视频免费| 欧美另类亚洲清纯唯美| 亚洲精品久久成人aⅴ小说| 国产精品欧美亚洲77777| 成年动漫av网址| 日韩欧美在线二视频 | 国产男女内射视频| 午夜福利欧美成人| 成熟少妇高潮喷水视频| 午夜福利乱码中文字幕| 黑人猛操日本美女一级片| 亚洲视频免费观看视频| 国产精品一区二区免费欧美| 高清视频免费观看一区二区| 午夜免费鲁丝| 亚洲精品国产精品久久久不卡| 啦啦啦视频在线资源免费观看| 欧美日韩中文字幕国产精品一区二区三区 | 手机成人av网站| 国产精品 国内视频| 亚洲五月色婷婷综合| 国产成人影院久久av| 精品国产超薄肉色丝袜足j| 欧美乱妇无乱码| 老熟妇仑乱视频hdxx| 久久这里只有精品19| 色在线成人网| 欧美日韩亚洲国产一区二区在线观看 | 精品一品国产午夜福利视频| 成人18禁在线播放| 亚洲第一青青草原| 国产精品电影一区二区三区 | 在线国产一区二区在线| 国产有黄有色有爽视频| www日本在线高清视频| 最近最新免费中文字幕在线| 亚洲一区中文字幕在线| 一级,二级,三级黄色视频| 亚洲精品成人av观看孕妇| 如日韩欧美国产精品一区二区三区| 欧美 日韩 精品 国产| 久久这里只有精品19| 欧美精品高潮呻吟av久久| 欧美另类亚洲清纯唯美| 80岁老熟妇乱子伦牲交| 两个人看的免费小视频| 亚洲精品粉嫩美女一区| 9191精品国产免费久久| 男女之事视频高清在线观看| 国产乱人伦免费视频| 桃红色精品国产亚洲av| 国产高清视频在线播放一区| 国产精品秋霞免费鲁丝片| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲精品国产一区二区精华液| 国产精品免费一区二区三区在线 | 不卡一级毛片| 欧美日韩一级在线毛片| 人人妻,人人澡人人爽秒播| 欧美一级毛片孕妇| 精品一区二区三卡| 国产人伦9x9x在线观看| 亚洲国产看品久久| 亚洲精品国产色婷婷电影| 亚洲成av片中文字幕在线观看| 成人18禁高潮啪啪吃奶动态图| 欧美激情 高清一区二区三区| av天堂久久9| 在线天堂中文资源库| 最近最新免费中文字幕在线| 黄色成人免费大全| 露出奶头的视频| 视频区图区小说| 女人高潮潮喷娇喘18禁视频| 免费久久久久久久精品成人欧美视频| 一区二区日韩欧美中文字幕| 日韩一卡2卡3卡4卡2021年| 成人av一区二区三区在线看| 91精品三级在线观看| 亚洲欧美一区二区三区黑人| 少妇猛男粗大的猛烈进出视频| 免费女性裸体啪啪无遮挡网站| 黑人巨大精品欧美一区二区mp4| 一级作爱视频免费观看| 精品亚洲成a人片在线观看| 亚洲一卡2卡3卡4卡5卡精品中文| 国产亚洲精品第一综合不卡| 国产亚洲欧美精品永久| 国产成人系列免费观看|