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

    Simulation of a Reverse Flow Reactor for the Catalytic Combustion of Lean Methane Emissions☆

    2014-07-12 08:33:00JiajinZhangZhigangLeiJianweiLiBiaohuaChen

    Jiajin Zhang,Zhigang Lei,Jianwei Li*,Biaohua Chen

    State Key Laboratory of Chemical Resource Engineering,Beijing University of Chemical Technology,Beijing 100029,China

    Simulation of a Reverse Flow Reactor for the Catalytic Combustion of Lean Methane Emissions☆

    Jiajin Zhang,Zhigang Lei,Jianwei Li*,Biaohua Chen

    State Key Laboratory of Chemical Resource Engineering,Beijing University of Chemical Technology,Beijing 100029,China

    A R T I C L E I N F o

    Article history:

    Received 24 December 2013

    Received in revised form 27 January 2014 Accepted 17 March 2014

    Available online 18 June 2014

    Methane

    Catalytic combustion

    Reverse flow reactor

    Modeling

    This work is focused on the performance prediction of pilot scale catalytic reverse flow reactors used for combustion of lean methane-air mixtures.An unsteady one-dimensional heterogeneous model for the reactor was established to account for the influence of the reactor wall on the heat transfer.Results of the simulation indicate that feed concentration,switch time and compensatory temperature impose important influence on the performance of the reactor.The amount of the heat extracted from the mid-section of the reactor can be optimized via adjusting the parameters mentioned above.At the optimal operating conditions,i.e.switching time of 400 s, feed concentration of1%(by volume),and insulation layer temperature of 343 K,the axial temperature of the reactor revealed a comparatively symmetrical“saddle”distribution,indicating a favorable operating status of the catalytic reverse flow reactor.

    ?2014 The Chemical Industry and Engineering Society of China,and Chemical Industry Press.All rights reserved.

    1.Introduction

    Global warming caused by greenhouse gases(GHG)in the atmosphere is currently a major concern[1].Primary GHG includes carbon dioxide,methane and other hydrocarbons,and oxides of nitrogen[2]. A significant source of methane is fugitive emissions released during the production and transportation of fossil fuels,for example,from underground coal mines[3,4].Methane has a global warming potential 23 times than carbon dioxide on a mass basis,therefore the complete combustion of methane will reduce equivalent carbon dioxide emissions by a factor of about 20,even though carbon dioxide is produced at the meantime[5].

    The combustion of fugitive emissions,in particular lean methane emissions,presents some problems[6].Fugitive emissions of relatively high concentration can be burnt in a fl are,however,the methane concentration in many fugitive emission streams is low(0.3%-1.0%,by volume)[4].These emissions are unbeneficial for both safety and climate change reasons,and also represent waste of energy which may be used as fuel in certain applications if captured[7].However, homogeneous combustion reactions are unfavorable at such low concentration.Catalytic combustion is an effective method of burning methane,which offers many advantages over conventional combustion, especially for lean mixtures[8,9].

    The main problem for catalytic combustion of lean methane mixtures is the necessity of a high reactor temperature.For the feed at ambient temperature,some pre-heating equipments would be required to achieve ignition temperature.Catalytic reverse flow reactor(RFR)is a very interesting alternative to overcome this problem[10,11].Catalytic RFR consists of a catalyst bed in which the feed flow direction is periodically reversed.Therefore,if the adequate switching time was ensured,the combustion heat would be restrained inside the reactor[11].

    In such a case,a quasi-periodic state operation may be achieved in which the reactor temperature pro file has a maximum value near the center of the reactor,which slowly oscillates within the reactor as the feed is switched periodically.Because the temperature in the center of the reactor can greatly exceed the insulation temperature,energy can be withdrawn from the center both as a means of control and as a source of useful energy[11].For example,RFR was used for catalytic decontamination of waste gases,as reported in[12],where a heat exchanger was added in the center of the reactor to remove heat.This ensured that the reactor did not over heat and deactivate the catalyst or damage the reactor.The heat removed from the reactor may be used for utilities such as heating a building or driving a small turbine.

    Computer modeling was a valuable tool as an aid in understanding the processes and assisting with design.Models specific to RFR have been reported for a variety of applications,and was used to investigate a variety of operating parameters,either from a theoretical viewpoint or coupled with experiments.Many models are related to one dimensional(1D),and either pseudo-homogeneous or heterogeneous models[13,14].Two-dimensional models have also been developed [15,16].

    Fig.1.Frame scheme for the reverse flow reactor with a middle heat exchanger[12].

    In this work,a scaled-up catalytic reverse flow reactor with heat exchanger added in the center of the reactor to remove heat for the combustion of fugitive emissions of methane were simulated by means of a 1D heterogeneous model,in which the intrinsic kinetics were measured independently in this work.Compared with the simulations of many years of research,the influences of heat extraction and insulation layer temperature which is used as compensation for the loss in the reactor wall were investigated.The aim of this work was to investigate the influence of operational conditions on the steady-state performances with high gaseous hourly space velocity(GHSV)and explore the feasibility of heat extraction from the middle of the reactor. Most of the relevant operating conditions,such as feed concentration, pressure drop,and switch time and insulation layer temperature were considered.Maximum temperature reached in the reactor,critical or maximum switch time allowable for stable operation,and favorable insulation layer temperature were determined for the reactor.

    2.Model Descriptions

    2.1.Physical model

    A RFR concept appropriate for relatively large flow rates was illustrated in Fig.1.The distribution of axial-radial reactor bed for RFR was shown in Fig.2.

    The reactor system has two separate stainless steel reactors with insulation layer installed around the external wall and connected by a heat exchanger at the bottom.The physical dimension of the reactor was shown in Table 1.The insulation layers are electrically heated and temperature controlled,so as to adjust the wall heat transfer condition for better RFR performance.Two three-way valves and associated transfer piping allow for either forward or reverse flow operating modes. Each reactor consisted of inert and catalyst sections and empty spaces in between.The inert sections were packed with blank Al2O3particles, whose physical properties were shown in Table 2.Reaction heat was extracted from the heat exchanger.

    Table 1 Geometry and physical properties for the RFR

    Regarding the solid bed,physical properties changed along the reactor because it was formed by two consecutive sections:upper inert and lower active.The active part was full of particles which was the mixture of catalyst(Pd/Al2O3)and inert spheres.The high activity of the catalyst was so adequate for the RFR operation that the central section was formed by a mixture of catalyst(50%,by volume)and inert particles (i.e.refractory ceramic spheres of the same diameter).The inert part was formed only by spherically shaped refractory ceramic material.

    In this work,the non-steady model was established for the combustion of methane in RFR in a heterogeneous and one-dimensional manner.The heterogeneous model can clarify the temperature difference between solid and fluid,which will be large in highly exothermic reactions.The only dimension taken into account in the model is the axial one,because for the reactor operating at the conditions of this work, temperature and velocity radial profiles are almost fl at(except in the vicinity of the reactor wall),and hence the difference in the results between one-dimensional(much easier to solve)and two-dimensional models is very small.

    The temperature of the wall insulation layer remained constant,in order to complement heat losses to the wall.So that the temperature radial profiles were almost fl at;and the reactor was assumed to be one-dimensional.Butit was necessary to take into account for the influence of the wall it self in RFR.That is to say,part of heat can be transferred axially along the wall to reach the inert part of bed.This heat transfer has a negative effect on the stability of the reactor,which had been demonstrated by research on a reactor of0.05 m diameter[17]: the reactor stability increases as the reactor diameter increases and the wall thermal conductivity decreases.Furthermore,the behavior of pilot-plant reactors,with relative small diameter and a high thermal conductivity wall(for example,constructed in steel),are strongly affected by this phenomenon,which must be taken into account in the simulation model.

    In the case of methane combustion,homogeneous reaction(combustion in gas phase)only occurs at very high temperature(>800°C). In this work,the reactor runs at much lower temperature so that the homogeneous combustion can be neglected.Ideal gas behavior for the gas phase was also assumed.

    Fig.2.Axialallocation of-reactor sections and heat exchanger for RFR.

    Table 2 Solid physical properties for particle beds[18]

    2.2.Basic equations

    The mathematical model for RFR was constituted by several partial differential and algebraic equations.The differential equations were obtained from mass and energy conservations,and applied separately to the gas and solid phases.The algebraic equations were the ones used to estimate the physical and transport properties.

    Mass balance for the gas phase:

    Energy balance for the gas phase:

    Mass balance for the solid phase with a 1.5 order reaction:

    Energy balance for the solid phase:

    Energy balance for the wall:

    The right hand side of Eq.(6)contained the axial heat conduction along the reactor wall,heat transfer between solid phase and reactor wall and that to reactor wall from insulation layer.

    Energy balance for the heat exchanger:

    2.3.Boundary and initial conditions

    The initial conditions were given as follows:

    Danckwerts boundary conditions were adopted.

    At the inlet:

    At the outlet:

    Fig.3.Diagram of the equipment used for kinetic experiments of lean methane catalytic combustion.1—air compressor;2—pressure maintaining valve;3—three-way valve;4—mass flow meter;5—one-way valve;6—methane bottle;7—pressure reducing valve;8—pressure gage;9—commingler;10—temperature display;11—temperature control instrument;12—reactor;13—condensator;14—gas-liquid separator;15—cut-off valve;16—GC.

    The outlet conditions for the reactor conform to Eq.(13)at the first half of the cycle.At the second half of the cycle,Eq.(14)is followed.

    2.4.Kinetic expressions

    In the development of the model,pseudo-1.5th order kinetics for the catalytic combustion of methane was assumed[18].This assumption had been experimentally checked on the sphere-shaped Pd/Al2O3catalyst which was a sample of catalyst for organic waste gas purification [HPA(KMK)]bought from Zhejiang University.[HPA(KMK)]catalyst is prepared by using palladium and platinum supported onγ-Al2O3. The related properties of the catalysts were shown in Table 2.The equipment used for kinetic experiments was shown as Fig.3.

    Before the experiments,the catalyst was ground to 100-250μm in order to neglect internal diffusion resistance.The influence of external diffusion was investigated by means of changing the loading amount of the catalyst particles and gaseous flow while the ratio of gaseous flow to catalyst loading being maintained constant.The constancy of estimated rate constant k versus change of GHSV indicated that the external diffusion influence could be neglected under the condition of GHSV(gas hourly space velocity)≥4600 ml·(g cat)-1·h-1.

    The required methane concentration was 1.0%(by volume),and temperature remained below 600°C.Forty-five experimental data for methane conversion were obtained under different space velocities and temperatures as listed in Table 3.

    The experimental data of k were correlated with a simple 1.5 order model.In this case the intrinsic model of lean methane catalytic combustion on the Pd/Al2O3catalyst was given as rA=kcA1.5.The rate parameter was expressed as k=k0exp(-E/RT),where k0was prexponential factor,and E is activation energy.

    The partialleast-square(PLS)and the Nelder-Mead simplex search methods were employed to estimate the parameters of kinetic model. The Arrhenius plot for experimental data of this model was shown in Fig.4,from which the estimated values of k0and E were obtained.

    The resulting apparent activation energy for this 1.5-order model was 82.64 kJ·mol-1and pre-exponential factor was 4.17×106L0.5· mol-0.5·s-1.The reaction heat for methane combustion(ΔHR)was -802.5 kJ·mol-1.Deactivation of the Pd/Al2O3catalyst was observed when the temperature was higher than 600°C.Therefore,the behavior of the catalyst cannot be adequately described by these parameters over 600°C.

    2.5.Transport properties

    The transport coefficients used to simulate the particle bed reactor were experimentally validated and reported in a previous work[19]. The effective dispersion coefficient(Deff),effective thermal conductivity (kGeff),gas-solid mass transfer coefficient(KG),gas-solid heat transfer coefficients(h),and Ergun equation to calculate the pressure drop (Δp),were listed in Table 4.

    It was assumed that the reaction heat was transported from the particle to the insulation layer through the reactor wall with the thickness of1.5 mm and thermal conductivity of19.51 W·m-1·K-1.The thermal conductivity 0.042 W·m-1·K-1of the insulation layer was much lower than the wall of reactor,which could prevent heat diffusion from the reactor wall to the ambient environment.

    Table 3 Conversion of methane at different temperatures and space velocities

    Fig.4.Arrhenius plots of the estimated values of kinetic parameters k for the power law model.

    2.6.Solution of mathematical model

    The present mathematical model consisted of5 partial differential equations and additional algebraic equations,in which physical,chemical and transport properties were calculated as mentioned above.The discretization of the mathematical model was performed by using the‘upwind scheme’.We modify the criterion of the second derivative control with a criterion based on a combination of the magnitude of second-order derivatives and that of first-order derivatives.With this criterion,the spacing distribution of grid can be well controlled.It may be better to put a little more points in the region of high first-order derivatives besides in the region of high second-order derivatives.The criterion can be described as[19]

    whereδ=0.02 and k1=1 and k2=0.5 were two selected nonnegative constant for catalysts and inert bed.The previous studies[20]proved that the simulation grid was independent with the parameterδ≤0.5. In this paper,to insure the simulation accuracy,δwas setas 0.02.

    Table 4 Correlations used for estimating internal transfer parameters in the simulation

    Table 5 Five primary variables D,F,E,and R used in Eq.(19)

    Fully implicit schemes were generally stable since a large time step can be chosen to save computational time.But completely implicit schemes usually required a series of iteration to ensure convergence. The weakly implicit algorithm described in[21]could avoid some unnecessary iteration,and could choose a large time step(0.1 s)if necessary.

    In a fully implicit scheme,all spatial derivatives are expressed using discrete variables at the n+1 time step,and only the temporal derivatives are involved variables at step n.The first-order derivatives were discretized by‘upwind scheme’[Eq.(16)]and the second-order derivatives were discretized by central difference method[Eq.(17)]:

    Fig.5.Schematic diagram of reverse flow reactor for methane catalytic combustion in experimental scale.1—methane;2—commingler;3—gas-liquid separator;4—mass flow meter; 5—pressure maintaining valve;6—electromagnetic valve;7—reactor;8-heat exchanger.A inert,B catalyst.

    Fig.6.Comparison of experimental and simulation results.

    Table 6 Operating parameters used in the simulation of base case

    Fig.7.Development of the center line axial temperature profile without heat extraction. 1—after the 25th cycle;2—after the 15th cycle;3—after the 1st cycle.

    Generally,a partial differential equation for RFR can be described as [20]

    After substitution,it becomes a series of algebraic equations of uin+1for the time instant tn+1,i.e.,

    For our 5 primary variables,functions D,F,E,and R are listed in the following Table 5.

    Also,the Crank-Nicolson technique and the corresponding predictor corrector method can be used to get higher order accuracy.In this work, MATLAB had been used to solve the ordinary differential equations for a particle bed reactor with feed flow periodically reversed by changing the spatial coordinates.The convergence residual of the simulation results was less than 10-5.

    2.7.Model validation

    The bench scale experiment was took place in the reactor shown in Fig.5.The internal diameter of the reactor was0.05 m,the total length of each reactor was 2.0 m,with the catalysts filling of 1 l to make a 0.5 m packed section(and the rest3.5 m packed with inert).The solid physical properties for particle beds were shown in Table 2.The operational condition was controlled at 1.0%(by volume)of CH4,20 m3·h-1handling capacity and a switch time of400 s.

    Fig.8.Methane concentration pro files in the reactor without heat extraction.1—after the 25th cycle;2—after the 15th cycle;3—after the 1st cycle.

    Fig.6(a)showed the experimental axial temperature pro file after operation of25 cycles with switch time of400 s when the reactor operated at a pseudo-steady state.Fig.6(b)showed the axial temperatureprofile by simulation at the same moment of operation.The comparison of them showed that the trend of temperature profile from the simulator fairly well matched the experimental results over a long period of time.Overall,the comparison results indicated that the model is able to predict the thermal energy trends in the RFR.

    Fig.9.Variation of the center line axial temperature profile after 13 cycles with heat extraction.1—before heat extraction;2—after heat extraction.

    3.Results and Discussion

    Numerical simulation was used to investigate the influence of operating parameters,as a prior step in the development of optimization and control strategies.From our previous experimental and limited modeling work reported in the literature[15-17],a number of variables were found to be important to the operation,such as thermal properties of the system,insulation layer temperature,fluid velocity,switch time,and feed concentrations,which were addressed in this work.By regulating these parameters,the reactor could reach the best performance,such as no catalyst bed‘run-away’,reactor extinction being avoided,high conversion rate of CH4,and more heat extraction.

    3.1.The base case

    Firstly,a base case without heat extraction was given,with the specific parameters shown in Table 6.

    The operation starts from an initial temperature profile without heat extraction from the heat exchanger at the bottom.The progression of the axial temperature profile was shown in Fig.7.The maximum temperature continued to increase every cycle before the 15th cycle. While the periodic flow reversal continues,the maximum temperature slightly decreases,and a stabilized state was obtained at the end of the 25th cycle.During the reaction,the reaction heat was accumulated in the reactor,the maximum temperature of reactor was high to 740 K (after 25th cycle)compared to 705 K(after 1st cycle).Fig.8 presented corresponding methane concentration profiles.It is visible from the concentration profiles that the combustion of methane was nearly completed at the position of reactor bottom while much more cycles were taken.As a result,more reaction heat was released to the system,producing a sharper gas temperature profile and higher maximum temperature at the end of the pseudo-steady state.

    With the operating parameters showed in Table 6,when the pseudo-steady state without heat extraction reached,the cooling system situated at the middle of the reactor started to operate.The cooling system that consumes heat requires that the water temperature at the outlet should be above 90°C and at the inlet below 25°C at the flow rate of30 kg·h-1.Figs.9,10 and 11 showed the axial gas temperature profile,methane concentration profile and the solid temperature profile,respectively,after 13 cycles with heat extraction.It appeared that the maximum temperature decreased from 750 K to 690 K,while the maximum solid temperature decreased from 720 K to 690 K after heat extraction,and a comparatively symmetrical“saddle”distribution was found in the axial gas temperature profile.Moreover,in Fig.9,because some released heat was extracted in the heat exchange range(2.4 m-3.9 m),the temperature of reactor was reduced to 610 K;the cooled gas recommenced reaction in the second stage reactor(3.9 m-4.3 m). Due to the lower gas temperature,the reaction rate was decreased, and only a small portion of reactant was combusted.In Fig.10,it can be observed that the outlet concentration of methane was increased slightly with heat extraction.This is explained considering the lower reaction rate by the heat removal.However,the concentration of methane was less than 0.3%(by volume)with some of heat extraction.Therefore, a favorable running status was achieved.In Fig.11,the solidtemperature profile over catalyst and inert sections was shown.At the entrance,because the heat transferred from gas to solid,the solid temperature showed a trend of decreasing firstly,and then increasing.All the simulated results indicate that it is possible to extract some of heat produced from the combustion in the reactor through the middle heat exchanger,with the reactor running properly at the same time.

    Fig.10.Variation of methane concentration profiles after 13 cycles with heat extraction.1—before heat extraction;2—after heat extraction.

    Fig.11.Variation of the solid temperature profile after 13 cycles with heat extraction.1—before heat extraction;2—after heat extraction.

    3.2.Effects of switch time

    Switch time was an important operating parameter,namely the time for finishing a forward flow and a reverse fl ow.The temperature profiles at a pseudo-steady state with different switch times without heat extraction were shown in Fig.12(a).In the simulation,it was observed that with longer switch time,it would take much more cycles to achieve the pseudo-steady state.For example,for the shortests witch time of200 s,the reactor reached a pseudo-steady state over 5 cycles. While the switch time was as long as 1000 s,a pseudo-steady state was achieved over 47 cycles.

    It can be seen from Fig.12(a)that the maximum temperature at the stabilized operation increased as the switch time increased without heat extraction.However,this trend was reversed when the heat was extracted away from the reactor,as shown in Fig.12(b).This effect can be explained by the increased reaction heat loss from the bed carried away by exhaust gas if the switch time is too long.But when the switch time was only 200 s,there was almost no heat to be extracted. Because the residence time of gas was less than the switch time,a larger fraction of the reaction energy was carried away by the effluent gas at the outlet,leading to a lower average gas temperature over the total reactor at the end of pseudo-steady state,as shown in Fig.13.

    Therefore,there should exist a critical switch time.In our case,the approximate switch time was within the range of 400 s to 800 s.In this range of switch time,the reactor could run properly without catalysts‘run-away’and‘reactor extinction’,while more combustion heat was able to transfer.

    Fig.12.Axial temperature profile at pseudo-steady state with different switch time(at the end of a cycle with flow along the z direction).

    Fig.13.Average temperature of the reactor with different switch time.

    Fig.14.Center line axial temperature profiles at different feed concentrations after heat extraction.1—feed concentration is 1%;2—feed concentration is 0.5%.

    3.3.Effectoffeed concentration

    The effect of feed concentration was illustrated in Fig.14,after25 full cycles when the inlet concentration of methane reduces from 1%to 0.5% by volume.It can be seen that when the feed concentration reduced to 0.5%while other operating condition including heat removal was maintained,the reaction energy was almostre moved,and thus temperature is not high enough to maintain the reaction.It was expected that the combustion would stop finally,if no other measure was taken to maintain the energy balance.

    Reducing the heat removal is an effective way to make the reaction proceeding with lower methane concentration.As shown in Fig.15, when the cooling system that the water temperature at the outlet should be below 70°C and at the inlet above 40°C at the flow rate of 25 kg·h-1;the methane concentration at outlet can decrease to 0.35%,closed to the status of inlet concentration of methane 1%.

    Above all,when the feed concentration was lower,the most effective way is to minimize the heat extraction as much as possible,in order to prevent the reactor extinction.

    Fig.16.Axial temperature profiles at different insulation layer temperatures.1—before heat extraction,Ta=323 K;2—before heat extraction,Ta=343 K;3—before heat extraction,Ta=373 K;4—before heat extraction,Ta=500;5—after heat extraction, Ta=343 K.

    3.4.Effect of insulation layer temperature(Ta)

    Fig.15.Concentration of methane profiles at different feed concentrations with cooling system changed.1—feed concentration is 1%;2—feed concentration is 0.5%.

    It had been reported[15,21-23]that the insulation layer played a significant role in the reactor operation from a scale-up standpoint. The insulation layer supplied energy for the system.As another point of view,Tacan be interpreted as insulation layer temperature.The effectof the insulation layer temperature can be shown independently by adjusting the temperature in the base case.For operation convenience, the insulation layer temperature was adjusted to an appointed value. Fig.16 illustrated the result of man-made increasing the insulation layer temperature(Ta)from 323 K to 500 K with switch time of 600 s, and feed concentration of1.0%.Simultaneously,the temperatures of reactor wall with different insulation layer temperatures were shown in Fig.17.From the graphs,it was found that at a temperature larger than 373 K,there was no significant difference in the gas temperature profile,while the reactor wall temperature increased obviously under the base case operating conditions.This condition leads to high energy loss in the reactor wall.

    Fig.17.Reactor wall temperature profiles at different insulation layer temperatures.1—before heat extraction,Ta=323 K;2—before heat extraction,Ta=343 K;3—before heat extraction,Ta=373 K;4—before heat extraction,Ta=500 K.5—after heat extraction,Ta=343 K.

    Besides,when the insulation layer temperature was lower than 323 K,the reaction will be extinguished if the insulation layer cannot supply enough energy for the system,as shown in Fig.18.During the operation,much of the reaction energy was transferred to the reactor wall,with the result that the self-heating of the reaction cannot be preserved.

    4.Conclusions

    In this work,the effect of operating parameters on the performance of a reverse-fl ow reactor had been investigated for the catalytic combustion of lean methane/air mixtures.The simulated calculation results demonstrate that under the suitable operating conditions,the system can achieve a favorable running status,without catalyst bed‘runaway’,reactor extinction,and more heat energy extraction.On the other hand,it was found that the extraction of thermal energy from the mid-section of the reactor was feasible.The temperature profile remained smooth and continuous in the case of desired energy extracting from heat exchanger.

    Fig.18.Methane concentrations at different insulation layer temperatures.1—Ta=323 K;2—Ta=343 K.

    The proposed simulation was appropriate for the optimization of RFR at a pilot scale,especially for the system with heat extraction. Actually,RFRwas a very intricately complex system.To fully understand allof the processes in the system,more experimental research must be done in the future work.

    Nomenclature

    Cpheat capacity,J·kg-1K-1

    c molar concentration,mol·m-3

    Deffaxial dispersion coefficient,m2·s-1

    DRinternal diameter of the reactor,m

    DR′internal diameter of the heat exchanger

    dpparticle bed particle diameter,mm

    dWreactor wall thickness,m

    E activation energy for heterogeneous reaction,J·mol-1

    ΔHRcombustion enthalpy,J·mol-1

    h gas-solid heat transfer coefficient,W·m-2.K-1

    KGgas-solid mass transfer coefficient,W·m-2.K-1

    k thermal conductivity,W·m-1K-1

    kShepre-exponential factor of heterogeneous reaction

    kGeffaxial effective thermal dispersion coefficient,W·m-1K-1

    R ideal gas constant,J·mol-1K-1

    LRreactor length,m

    rAheterogeneous reaction rate per catalyst surface

    Δp pressure drop,Pa

    T temperature,K

    u superficial velocity,m·s-1

    z axial length along the reactor,m

    αvsurface-volume ratio,m2·m-3

    δpositive number was given to judge the step size

    εbbed porosity

    ρdensity,kg·m-3

    Superscripts

    dinertpacking

    he heterogeneous reaction

    ′heat exchanger

    Subscripts

    A reaction compound(methane)

    a insulation layer

    d inert stuffing

    G gas phase

    S solid phase

    W reactor wall

    0 reactor inlet conditions

    [1]J.P.Guo,C.D.Zhou,Greenhouse gas emissions and mitigation measures in Chinese agroecosystems,Agric.For.Meteorol.142(2-4)(2007)270-277.

    [2]K.Ito,Y.Uchiyama,T.Takeshita,H.Hayashibe,Study on GHG control scenarios by life cycle analysis—world energy outlook until 2100,Energy Convers.Manag.38 (7)(1997)607-614.

    [3]C.?.Karacan,F.A.Ruiz,M.Cote,S.Phipps,Coalmine methane:a review of capture and utilization practices with benefits to mining safety and to greenhouse gas reduction,Int.J.CoalGeol.86(2-3)(2011)121-156.

    [4]Y.P.Cheng,L.Wang,X.L.Zhang,Environmental impact of coal mine methane emissions and responding strategies in China,Int.J.Greenhouse Gas Control5(1) (2011)157-166.

    [5]R.E.Hayes,Catalytic solutions for fugitive methane emissions in the oil and gas sector,Chem.Eng.Sci.59(2)(2004)4073-4080.

    [6]S.Andrea,S.B.Paola,L.Gianluca,P.Raffaele,R.Gennaro,Combustion of methane hydrogen mixtures on catalytic tablets,Chem.Eng.J.154(1-3)(2009)315-324.

    [7]R.E.Hayes,S.T.Kolaczkowski,Introduction to Catalytic Combustion,Gordon and Breach,UK,1997.323-325.

    [8]N.Russo,P.Palmisano,D.Fino,Pd substitution effects on perovskite catalyst activity for methane emission control,Chem.Eng.J.154(1-3)(2009)137-141.

    [9]K.Xu,Z.L.Liu,H.He,S.Y.Cheng,C.F.Ma,Experimental study on emission control of premixed catalytic combustion of natural gas using preheated air,Chin.J.Chem.Eng. 15(1)(2007)68-74.

    [10]Y.S.Matros,G.A.Bunimovich,Reverse-fl ow operation in fixed bed catalytic reactors, Catal.Rev.38(1)(1996)1-67.

    [11]P.Marin,S.Ordonez,F.V.Diez,Procedures for heatrecovery in the catalytic combustion oflean methane-air mixtures in a reverse flow reactor,Chem.Eng.J.147(2-3) (2009)356-365.

    [12]K.Gosiewskia,K.Warmuzinski,Effect of the mode of heat withdrawal on the asymmetry of temperature profiles in reverse-fl ow reactors.Catalytic combustion of methane as a test case,Chem.Eng.Sci.62(2)(2007)2679-2689.

    [13]B.Liu,R.E.Hayes,M.D.Checkel,M.Zheng,E.Mirosh,Reversing flow catalytic converter for a natural gas/diesel dual fuel engine,Chem.Eng.Sci.56(2)(2001) 2641-2685.

    [14]A.G.Miguel,O.Hevia Salvador,V.D.Fernando,Effect of the catalyst properties on the performance of a reverse flow reactor for methane combustion in lean mixtures, Chem.Eng.J.129(1)(2009)1-10.

    [15]S.Salomons,R.E.Hayes,M.Poirier,H.Sapoundjiev,Modelling a reverse flow reactor for the catalytic combustion of fugitive methane emissions,Comput.Chem.Eng.28 (1)(2004)1599-1610.

    [16]R.Litto,R.E.Hayes,H.Sapoundjiev,A.Fuxman,F.Forbes,B.Liu,F.Bertrand, Optimization of a flow reversal reactor for the catalytic combustion of lean methane mixtures,Catal.Today 117(1)(2006)536-542.

    [17]P.Marin,M.A.G.Hevia,S.Ordonez,F.V.Drez,Combustion ofmethane lean mixtures in reverse flow reactors:comparison between packed and structured catalyst beds, Catal.Today 105(3-4)(2005)701-708.

    [18]Y.H.Chin,D.E.Resasco,Catalytic oxidation of methane on supported palladium under lean conditions:kinetics,structures and properties,Catalysis 14(1)(1999) 1-39.

    [19]X.K.Niu,Modeling of Reverse Flow Catalytic Combustion Process for Air Purification, Ph.D.Thesis Beijing University of Chemical Technology,China,2003.(in Chinese).

    [20]J.Y.Huang,C.Y.Li,A robust adaptive algorithm for flow reversal fi xed-bed flow reversal fixed-bed reactor models,Process 2nd Joint China/USA chemical engineering Conference,Chem.Ind.Press,Beijing,1997,pp.266-270.

    [21]S.Salomons,R.E.Hayes,M.Poirier,H.Sapoundjiev,Flow reversal reactor for the catalytic combustion of lean methane mixtures,Catal.Today 83(1)(2003)59-69.

    [22]A.Kushwaha,R.E.Hayes,M.Poirier,H.Sapoundjiev,Effect of reactor internal properties on the performance of a flow reversal catalytic reactor for methane combustion,Chem.Eng.Sci.59(1)(2004)4081-4093.

    [23]A.Kushwaha,M.Poirier,R.E.Hayes,H.Sapoundjiev,Heat extraction from a flow reversal reactor in lean methane combustion,Chem.Eng.Res.Des.83(1)(2005) 205-213.

    ☆Supported by the National High Technology Research and Development Program of China(2006AA030201).

    *Corresponding author.

    E-mailaddress:lijw@mail.buct.edu.cn(J.Li).

    夜夜爽夜夜爽视频| 在线观看三级黄色| av卡一久久| 国产精品国产三级国产专区5o| 精品亚洲成国产av| 国产精品一区二区在线观看99| 丝袜人妻中文字幕| 亚洲av福利一区| 最近中文字幕高清免费大全6| 亚洲人与动物交配视频| 国产一区有黄有色的免费视频| 精品酒店卫生间| 大香蕉久久成人网| 黄片播放在线免费| 国产精品久久久久久精品古装| 亚洲伊人色综图| 免费大片黄手机在线观看| 亚洲精品久久久久久婷婷小说| 大片电影免费在线观看免费| 免费日韩欧美在线观看| 精品一区二区三卡| 日本免费在线观看一区| 国产一区二区在线观看av| 91在线精品国自产拍蜜月| 黄网站色视频无遮挡免费观看| 成年动漫av网址| 天堂中文最新版在线下载| 亚洲伊人久久精品综合| 国语对白做爰xxxⅹ性视频网站| 亚洲精品国产av蜜桃| 亚洲久久久国产精品| 亚洲精品久久久久久婷婷小说| 欧美国产精品一级二级三级| 国产片特级美女逼逼视频| 国产亚洲av片在线观看秒播厂| 五月开心婷婷网| 成人黄色视频免费在线看| 在线观看国产h片| 水蜜桃什么品种好| 赤兔流量卡办理| 日韩人妻精品一区2区三区| 亚洲综合色网址| 一本大道久久a久久精品| 精品亚洲成a人片在线观看| 成人免费观看视频高清| 日本vs欧美在线观看视频| 久久韩国三级中文字幕| 中文字幕免费在线视频6| 亚洲在久久综合| 日本与韩国留学比较| 国产视频首页在线观看| 久久久精品区二区三区| 亚洲,欧美精品.| 亚洲精品日韩在线中文字幕| 久久人人爽av亚洲精品天堂| 国产无遮挡羞羞视频在线观看| 蜜桃在线观看..| 亚洲一区二区三区欧美精品| 制服人妻中文乱码| 婷婷成人精品国产| videosex国产| 一本—道久久a久久精品蜜桃钙片| 18禁在线无遮挡免费观看视频| 九九爱精品视频在线观看| 黄色一级大片看看| av免费在线看不卡| 热99国产精品久久久久久7| 丁香六月天网| 亚洲欧美精品自产自拍| 国产精品欧美亚洲77777| 国产 精品1| 久久精品aⅴ一区二区三区四区 | 人人妻人人澡人人看| 中文字幕av电影在线播放| 久热久热在线精品观看| 国产男女超爽视频在线观看| 高清av免费在线| 国产精品99久久99久久久不卡 | 国产亚洲精品第一综合不卡 | a级毛色黄片| 日韩在线高清观看一区二区三区| 男人操女人黄网站| 久久午夜福利片| 国产深夜福利视频在线观看| 久久99一区二区三区| 一级片'在线观看视频| 欧美日韩av久久| 另类亚洲欧美激情| 一区在线观看完整版| 看免费av毛片| 老女人水多毛片| av.在线天堂| 国产精品人妻久久久影院| 欧美国产精品va在线观看不卡| 亚洲第一av免费看| 只有这里有精品99| 国产成人免费无遮挡视频| 在线亚洲精品国产二区图片欧美| 国产av精品麻豆| 成人国产av品久久久| 伊人久久国产一区二区| 男女免费视频国产| 欧美另类一区| 妹子高潮喷水视频| 久久久久久久久久人人人人人人| 一级片'在线观看视频| 九色亚洲精品在线播放| 欧美日本中文国产一区发布| 成人国产av品久久久| 久久久久久久精品精品| 熟妇人妻不卡中文字幕| 久久久久久久久久人人人人人人| 欧美精品国产亚洲| 丝袜脚勾引网站| 国产又爽黄色视频| 大香蕉久久成人网| 午夜福利影视在线免费观看| 在线观看免费视频网站a站| 亚洲一区二区三区欧美精品| 精品少妇内射三级| 秋霞伦理黄片| av福利片在线| 国产精品久久久久久精品电影小说| 免费看不卡的av| 午夜日本视频在线| 观看美女的网站| 国产欧美日韩一区二区三区在线| 狠狠婷婷综合久久久久久88av| 中文精品一卡2卡3卡4更新| 午夜激情av网站| 免费在线观看完整版高清| 成人黄色视频免费在线看| 王馨瑶露胸无遮挡在线观看| 久久精品aⅴ一区二区三区四区 | 爱豆传媒免费全集在线观看| 97人妻天天添夜夜摸| 两个人免费观看高清视频| 国产一区二区激情短视频 | 国产成人91sexporn| 最近手机中文字幕大全| 日韩成人伦理影院| 亚洲国产精品专区欧美| 亚洲五月色婷婷综合| 亚洲在久久综合| 亚洲综合色惰| 欧美成人午夜精品| 久久国产精品大桥未久av| 国产精品人妻久久久影院| 久久久久国产精品人妻一区二区| 成人免费观看视频高清| 久久ye,这里只有精品| 欧美亚洲日本最大视频资源| 久久人妻熟女aⅴ| 一级毛片电影观看| 高清在线视频一区二区三区| 一级片免费观看大全| 国产日韩一区二区三区精品不卡| 乱码一卡2卡4卡精品| 欧美精品一区二区大全| 一二三四中文在线观看免费高清| 午夜福利乱码中文字幕| 国产国语露脸激情在线看| 一本大道久久a久久精品| 蜜桃在线观看..| 免费播放大片免费观看视频在线观看| 中文字幕精品免费在线观看视频 | 一二三四中文在线观看免费高清| 欧美变态另类bdsm刘玥| av免费观看日本| 亚洲精品国产av成人精品| 考比视频在线观看| 国产成人91sexporn| 最近2019中文字幕mv第一页| 精品国产一区二区三区久久久樱花| 亚洲内射少妇av| 欧美丝袜亚洲另类| 亚洲精品aⅴ在线观看| 黑人欧美特级aaaaaa片| 一级爰片在线观看| 91国产中文字幕| 精品久久久精品久久久| 国产av一区二区精品久久| 免费观看a级毛片全部| 午夜福利影视在线免费观看| 国产一区二区三区综合在线观看 | 免费看不卡的av| 国产极品粉嫩免费观看在线| 99久国产av精品国产电影| 国产精品不卡视频一区二区| 一本色道久久久久久精品综合| 免费黄色在线免费观看| 色94色欧美一区二区| 又粗又硬又长又爽又黄的视频| 热re99久久精品国产66热6| 国产男女超爽视频在线观看| 视频中文字幕在线观看| 亚洲国产看品久久| 国产精品秋霞免费鲁丝片| 曰老女人黄片| 免费高清在线观看日韩| xxxhd国产人妻xxx| 久久久精品免费免费高清| av免费观看日本| 亚洲精品日本国产第一区| 99国产精品免费福利视频| 色吧在线观看| 久久精品久久久久久久性| 久久国产精品大桥未久av| 狂野欧美激情性xxxx在线观看| 亚洲美女搞黄在线观看| 欧美日韩成人在线一区二区| 国产精品熟女久久久久浪| 久久精品国产亚洲av涩爱| 亚洲国产欧美在线一区| 欧美成人午夜免费资源| 欧美日韩视频精品一区| 熟女人妻精品中文字幕| 99视频精品全部免费 在线| 国产亚洲最大av| 久久精品国产亚洲av天美| 精品人妻偷拍中文字幕| 我要看黄色一级片免费的| 久久亚洲国产成人精品v| 中国三级夫妇交换| 亚洲国产成人一精品久久久| 日本与韩国留学比较| 大香蕉97超碰在线| 不卡视频在线观看欧美| 成人影院久久| 久久精品熟女亚洲av麻豆精品| 日韩中字成人| 国产老妇伦熟女老妇高清| 久久久久精品人妻al黑| 中文字幕最新亚洲高清| 国产色婷婷99| 日日撸夜夜添| 久久久久久久大尺度免费视频| 永久免费av网站大全| 久久97久久精品| 80岁老熟妇乱子伦牲交| 亚洲av男天堂| 三级国产精品片| 一级毛片黄色毛片免费观看视频| 久久久久久久精品精品| a级片在线免费高清观看视频| 亚洲国产欧美日韩在线播放| 最近中文字幕2019免费版| 国产精品女同一区二区软件| 亚洲精品久久成人aⅴ小说| 少妇人妻 视频| 少妇人妻精品综合一区二区| 亚洲成人av在线免费| 中国国产av一级| 国产成人精品一,二区| 免费观看a级毛片全部| 熟女人妻精品中文字幕| 美女脱内裤让男人舔精品视频| 免费大片18禁| a级片在线免费高清观看视频| 黑人高潮一二区| 黄色 视频免费看| 人成视频在线观看免费观看| 国产成人a∨麻豆精品| 18+在线观看网站| 桃花免费在线播放| 亚洲国产精品国产精品| 国产一区二区在线观看日韩| 国产精品一二三区在线看| 男女免费视频国产| 国产精品久久久久成人av| 欧美日本中文国产一区发布| 观看美女的网站| 成人午夜精彩视频在线观看| 亚洲国产欧美在线一区| 少妇精品久久久久久久| 看非洲黑人一级黄片| 51国产日韩欧美| 插逼视频在线观看| 午夜91福利影院| 亚洲,欧美,日韩| 2018国产大陆天天弄谢| 久久久a久久爽久久v久久| 久久精品久久久久久噜噜老黄| 三级国产精品片| 天天操日日干夜夜撸| 精品久久久精品久久久| 亚洲成人av在线免费| av在线播放精品| 尾随美女入室| 亚洲国产精品一区三区| 国产精品不卡视频一区二区| 狂野欧美激情性xxxx在线观看| 国产xxxxx性猛交| 国产精品秋霞免费鲁丝片| 美女视频免费永久观看网站| 国产精品偷伦视频观看了| 黄色视频在线播放观看不卡| 天天操日日干夜夜撸| 一本大道久久a久久精品| 久久久国产精品麻豆| 亚洲情色 制服丝袜| 在线观看免费视频网站a站| 成人国产av品久久久| 国产精品嫩草影院av在线观看| 伦理电影大哥的女人| 精品第一国产精品| 中国美白少妇内射xxxbb| 午夜激情久久久久久久| 天天影视国产精品| 在线观看三级黄色| 国产免费又黄又爽又色| 少妇的丰满在线观看| 内地一区二区视频在线| 精品国产乱码久久久久久小说| 激情视频va一区二区三区| av网站免费在线观看视频| 性色avwww在线观看| 女性生殖器流出的白浆| 美女xxoo啪啪120秒动态图| 免费人成在线观看视频色| 18禁裸乳无遮挡动漫免费视频| 精品99又大又爽又粗少妇毛片| 只有这里有精品99| 丰满饥渴人妻一区二区三| 午夜日本视频在线| 欧美成人午夜精品| 国产精品偷伦视频观看了| 少妇 在线观看| av网站免费在线观看视频| 日本欧美国产在线视频| 国产免费一级a男人的天堂| 亚洲欧美清纯卡通| 亚洲第一区二区三区不卡| 18禁在线无遮挡免费观看视频| 在线天堂中文资源库| 在线观看美女被高潮喷水网站| 欧美老熟妇乱子伦牲交| 久久影院123| 日韩欧美一区视频在线观看| 男女国产视频网站| 热99久久久久精品小说推荐| av一本久久久久| 久久精品aⅴ一区二区三区四区 | 黄色毛片三级朝国网站| videosex国产| 成人无遮挡网站| 欧美最新免费一区二区三区| 一区二区三区乱码不卡18| 看十八女毛片水多多多| av国产久精品久网站免费入址| 久久久久国产网址| 中文字幕av电影在线播放| www.av在线官网国产| 人人妻人人添人人爽欧美一区卜| 欧美激情国产日韩精品一区| 天美传媒精品一区二区| 久久精品国产亚洲av涩爱| 天堂8中文在线网| 亚洲欧美成人精品一区二区| 视频中文字幕在线观看| 九草在线视频观看| 亚洲婷婷狠狠爱综合网| 全区人妻精品视频| 成年人午夜在线观看视频| 久久人妻熟女aⅴ| 大香蕉97超碰在线| 欧美精品一区二区大全| 亚洲精品成人av观看孕妇| 亚洲激情五月婷婷啪啪| 韩国高清视频一区二区三区| 街头女战士在线观看网站| 欧美成人午夜精品| 亚洲熟女精品中文字幕| 如何舔出高潮| 亚洲av成人精品一二三区| 国产欧美日韩一区二区三区在线| 日韩伦理黄色片| 18禁观看日本| 国产精品秋霞免费鲁丝片| 亚洲伊人色综图| 美女国产视频在线观看| 日本欧美视频一区| 黑人猛操日本美女一级片| 精品人妻在线不人妻| 有码 亚洲区| 午夜福利乱码中文字幕| 青春草亚洲视频在线观看| 欧美激情 高清一区二区三区| 丰满迷人的少妇在线观看| 亚洲精品一二三| 久久ye,这里只有精品| 久久亚洲国产成人精品v| 人人澡人人妻人| 捣出白浆h1v1| 久久午夜福利片| 熟女电影av网| 曰老女人黄片| 免费人妻精品一区二区三区视频| 超碰97精品在线观看| 日本色播在线视频| 熟妇人妻不卡中文字幕| 中国三级夫妇交换| 亚洲精品久久成人aⅴ小说| 男女边吃奶边做爰视频| 久久av网站| 久久人人爽av亚洲精品天堂| 亚洲中文av在线| 99九九在线精品视频| 免费观看av网站的网址| 亚洲欧美一区二区三区国产| 91aial.com中文字幕在线观看| 晚上一个人看的免费电影| 搡女人真爽免费视频火全软件| 精品亚洲乱码少妇综合久久| 国产成人精品婷婷| 午夜免费男女啪啪视频观看| 久久人人爽av亚洲精品天堂| 久久久国产精品麻豆| 一级毛片我不卡| 日本免费在线观看一区| 亚洲人与动物交配视频| 日本色播在线视频| 新久久久久国产一级毛片| 自拍欧美九色日韩亚洲蝌蚪91| 亚洲国产看品久久| 欧美性感艳星| 久久久国产精品麻豆| 夜夜骑夜夜射夜夜干| 中文天堂在线官网| 免费看av在线观看网站| 精品一区二区三区四区五区乱码 | 亚洲精品中文字幕在线视频| 亚洲欧美成人综合另类久久久| 亚洲精品一区蜜桃| xxx大片免费视频| 国产av一区二区精品久久| 国产精品女同一区二区软件| 亚洲成色77777| 欧美 亚洲 国产 日韩一| 伊人亚洲综合成人网| 国产乱来视频区| 国产 一区精品| 纯流量卡能插随身wifi吗| 七月丁香在线播放| 成年动漫av网址| 国产日韩一区二区三区精品不卡| 性高湖久久久久久久久免费观看| 欧美97在线视频| 香蕉国产在线看| 在线亚洲精品国产二区图片欧美| 亚洲精品乱久久久久久| 国产精品人妻久久久久久| av国产精品久久久久影院| 性高湖久久久久久久久免费观看| av不卡在线播放| 母亲3免费完整高清在线观看 | 精品国产国语对白av| 日韩,欧美,国产一区二区三区| √禁漫天堂资源中文www| 十八禁高潮呻吟视频| 欧美亚洲 丝袜 人妻 在线| 欧美人与性动交α欧美软件 | 看免费av毛片| 高清视频免费观看一区二区| 午夜久久久在线观看| 国产高清三级在线| 日韩不卡一区二区三区视频在线| 亚洲欧美一区二区三区国产| 免费女性裸体啪啪无遮挡网站| 在线观看免费高清a一片| 人妻一区二区av| 国产亚洲午夜精品一区二区久久| 涩涩av久久男人的天堂| 日本与韩国留学比较| 久久久久国产精品人妻一区二区| 寂寞人妻少妇视频99o| 美女福利国产在线| 国产毛片在线视频| www.熟女人妻精品国产 | 久久久国产一区二区| 国产又爽黄色视频| 久久综合国产亚洲精品| 国产成人精品在线电影| 精品一区二区三区视频在线| av播播在线观看一区| 国产在线一区二区三区精| 国产男女内射视频| 亚洲久久久国产精品| 一级毛片我不卡| 1024视频免费在线观看| 内地一区二区视频在线| 精品久久国产蜜桃| 不卡视频在线观看欧美| 2021少妇久久久久久久久久久| 国产永久视频网站| 久久久久久久国产电影| 九色成人免费人妻av| 亚洲精品美女久久久久99蜜臀 | 日产精品乱码卡一卡2卡三| 女的被弄到高潮叫床怎么办| 大码成人一级视频| 亚洲性久久影院| 亚洲成人手机| 日韩av在线免费看完整版不卡| 国产成人欧美| 成人国产av品久久久| 99久久人妻综合| 亚洲欧洲日产国产| 女性生殖器流出的白浆| 熟女人妻精品中文字幕| 视频区图区小说| 在线观看免费视频网站a站| av国产精品久久久久影院| 亚洲欧美一区二区三区黑人 | 精品少妇黑人巨大在线播放| 欧美 亚洲 国产 日韩一| av播播在线观看一区| 精品亚洲成国产av| 两个人看的免费小视频| 中文字幕人妻丝袜制服| 三级国产精品片| 日日撸夜夜添| 久久av网站| 春色校园在线视频观看| 免费高清在线观看视频在线观看| av在线app专区| av天堂久久9| 波野结衣二区三区在线| 秋霞在线观看毛片| 国产成人免费观看mmmm| 免费在线观看完整版高清| av电影中文网址| 国产精品偷伦视频观看了| 91成人精品电影| 亚洲情色 制服丝袜| 国产精品久久久久久久电影| 另类精品久久| 日韩三级伦理在线观看| 精品人妻熟女毛片av久久网站| 免费日韩欧美在线观看| 最黄视频免费看| 久久国产精品男人的天堂亚洲 | 你懂的网址亚洲精品在线观看| 建设人人有责人人尽责人人享有的| 极品人妻少妇av视频| 国产午夜精品一二区理论片| av在线播放精品| 免费av不卡在线播放| 日韩免费高清中文字幕av| 亚洲av日韩在线播放| 秋霞在线观看毛片| 99久久中文字幕三级久久日本| 人妻一区二区av| 国产精品熟女久久久久浪| 婷婷色综合大香蕉| 日韩一区二区视频免费看| 国产成人免费无遮挡视频| 一级,二级,三级黄色视频| 成年av动漫网址| 国产亚洲最大av| 国产有黄有色有爽视频| 婷婷成人精品国产| 国产激情久久老熟女| 黄色一级大片看看| 少妇 在线观看| 一级毛片 在线播放| 欧美激情国产日韩精品一区| 精品一区二区三区视频在线| 亚洲欧美成人综合另类久久久| 国产成人午夜福利电影在线观看| av电影中文网址| 国产 一区精品| 久久婷婷青草| 午夜福利视频精品| 欧美亚洲 丝袜 人妻 在线| 又黄又爽又刺激的免费视频.| 久久99热这里只频精品6学生| 男女免费视频国产| 亚洲美女搞黄在线观看| 青春草亚洲视频在线观看| 永久网站在线| 18禁观看日本| av天堂久久9| 国产免费福利视频在线观看| 飞空精品影院首页| 99久久中文字幕三级久久日本| 成年动漫av网址| 国精品久久久久久国模美| 91午夜精品亚洲一区二区三区| 国产深夜福利视频在线观看| 大片电影免费在线观看免费| 亚洲性久久影院| 亚洲综合精品二区| 99香蕉大伊视频| 精品人妻在线不人妻| 午夜福利影视在线免费观看| 观看av在线不卡| 久热这里只有精品99| 王馨瑶露胸无遮挡在线观看| 国产一区二区三区av在线| h视频一区二区三区| 9191精品国产免费久久| 亚洲熟女精品中文字幕| 国产av一区二区精品久久| 又黄又粗又硬又大视频| 婷婷色综合大香蕉| 蜜桃在线观看..| 老司机影院成人| 人人妻人人澡人人看| 国产成人一区二区在线| 人妻少妇偷人精品九色| 99久久中文字幕三级久久日本| 国产精品人妻久久久久久| 一区二区三区四区激情视频| 免费女性裸体啪啪无遮挡网站|