首页    期刊浏览 2024年10月01日 星期二
登录注册

文章基本信息

  • 标题:Thermal state and hydrodynamics of evaporating hydrocarbon droplets. 1. A possibility of natural circulation of the liquid in the droplet/Garuojanciu skysto angliavandenilio laseliu termine busena ir hidrodinamika. 1. Savaimines skyscio cirkuliacijos kilimo galimybe.
  • 作者:Miliauskas, G. ; Talubinskas, J. ; Adomavicius, A.
  • 期刊名称:Mechanika
  • 印刷版ISSN:1392-1207
  • 出版年度:2012
  • 期号:March
  • 语种:English
  • 出版社:Kauno Technologijos Universitetas
  • 摘要:In boilers, internal combustion engines and rocket engines liquid fuel (generally liquid hydrocarbons or their mixtures) is combusted in the sprayed form. Droplet evaporation speed is one of the factors essentially influencing the efficiency of the combustion process. The evaluation of evaporation speed using direct measurements in the flaming torch is very complicated. Thus, a lot of theoretical droplet evaporation research methods have been developed and applied [1-3]. A wide range of heat and mass transfer process research tasks with the experimentally determined values of the critical Rayleight criterion for the spherical water volume [4] in sprayed liquid systems is reflected by a variety of possible research methods. All these methods have a common starting point--the modelling of the evaporated liquid droplet, usually known as the "droplet task". The "droplet task" is integrated and consists of "internal" and "external" tasks. Both of them depend on multiple heat and mass transfer processes in the droplet and its surroundings. Their peculiarities should be taken into account when the mathematical models for the "droplet task" are being prepared. Generally, in the "external droplet task", heating of the droplet by thermal convection and radiation simultaneously needs to be evaluated. Thus its model includes heating and evaporation models. The correlation between the processes of heat convection and radiation in the gases surrounding the droplet is negligible, so both processes can be described independently. The influence Stefan's hydrodynamic flow has on the speed of droplet evaporation and intensity of heating by convection, however, can not be neglected. When evaporation models are based on the methods of similarity theory, the influence of Stefan's hydrodynamics flow on droplet's evaporation speed is evaluated by the Spalding mass transfer parameter [B.sub.M] function [Sh.sub.m] [equivalent to] Shf([B.sub.M]) [2, 5, 6]. In analytical models of the Stefan hydrodynamic flow, the influence on the speed of droplet evaporation is evaluated by logarithmic function of liquid vapour pressure difference between the surface and the surroundings of the droplet [7, 8].
  • 关键词:Evaporation;Fluid dynamics;Hydrocarbons;Hydrodynamics;Hydrofoil boats;Nuclear radiation;Thermodynamics

Thermal state and hydrodynamics of evaporating hydrocarbon droplets. 1. A possibility of natural circulation of the liquid in the droplet/Garuojanciu skysto angliavandenilio laseliu termine busena ir hidrodinamika. 1. Savaimines skyscio cirkuliacijos kilimo galimybe.


Miliauskas, G. ; Talubinskas, J. ; Adomavicius, A. 等


1. Introduction

In boilers, internal combustion engines and rocket engines liquid fuel (generally liquid hydrocarbons or their mixtures) is combusted in the sprayed form. Droplet evaporation speed is one of the factors essentially influencing the efficiency of the combustion process. The evaluation of evaporation speed using direct measurements in the flaming torch is very complicated. Thus, a lot of theoretical droplet evaporation research methods have been developed and applied [1-3]. A wide range of heat and mass transfer process research tasks with the experimentally determined values of the critical Rayleight criterion for the spherical water volume [4] in sprayed liquid systems is reflected by a variety of possible research methods. All these methods have a common starting point--the modelling of the evaporated liquid droplet, usually known as the "droplet task". The "droplet task" is integrated and consists of "internal" and "external" tasks. Both of them depend on multiple heat and mass transfer processes in the droplet and its surroundings. Their peculiarities should be taken into account when the mathematical models for the "droplet task" are being prepared. Generally, in the "external droplet task", heating of the droplet by thermal convection and radiation simultaneously needs to be evaluated. Thus its model includes heating and evaporation models. The correlation between the processes of heat convection and radiation in the gases surrounding the droplet is negligible, so both processes can be described independently. The influence Stefan's hydrodynamic flow has on the speed of droplet evaporation and intensity of heating by convection, however, can not be neglected. When evaporation models are based on the methods of similarity theory, the influence of Stefan's hydrodynamics flow on droplet's evaporation speed is evaluated by the Spalding mass transfer parameter [B.sub.M] function [Sh.sub.m] [equivalent to] Shf([B.sub.M]) [2, 5, 6]. In analytical models of the Stefan hydrodynamic flow, the influence on the speed of droplet evaporation is evaluated by logarithmic function of liquid vapour pressure difference between the surface and the surroundings of the droplet [7, 8].

The influence Stefan's hydrodynamic flow has on the intensity of the evaporating droplet heating by convection is evaluated by multiplying the convective solid particle heating criterion by the Spalding heat transfer corrective coefficient [B.sub.T] function. The corrective function in the models of equilibrated evaporation of a droplet, heated by convection is based on the classical expression of Spalding's heat transfer parameter. Equilibrated evaporation is the case when all heat transferred to the droplet from its surroundings is used up to evaporate the liquid. In the case of combined heating of the droplet, the surrounding heat transfer by radiation must be taken into account. When nonequilibrated evaporation processes are being modelled, the intensity by which the liquid inside the droplet is being heated needs to be evaluated [6, 8]. Because of that, some additional modelling problems and tasks typical to the "internal droplet task" are necessary to be solved. One of these problems is that the heat transfer inside of the droplet is combined. The correlation between heat transfer by conduction, convection and radiation in the liquid can not be neglected. Modelling of these interactions is made more complicated by the need to take spectral radiation characteristics of semitransparent liquid absorption and light's optical effects on the surface of the droplet into account. These effects are assessed in the spectral radiation models [8-11]. They are effective in the case of integrated heat transfer inside of the droplet by thermal conductivity and radiation. However, the use of such models when heat convection in the liquid occurs is a very delicate matter. That leaves the problem of stability assessment in the droplets of evaporating liquid. Spontaneous circulation in the nonisothermal liquid droplets may be the effect of Archimedes forces, while forced circulation inside of the droplet can cause friction on the surface of it when the speed of the moving droplet is different than that of the surrounding gasses. Thus, when the heating and evaporation of the fuel droplet is being modelled, the hydrodynamic conditions inside of it should be evaluated.

A liquid gravitational film forms on the channel's surface when a two-phase flow flows in it [12]. The interaction between the liquid film and steam, in turn, changes the course of condensation [13]. In order to model heat and mass transfer processes "droplet" and "film" models should be aggregated. The fluid hydrodynamic regime problem is relevant to both of them. Therefore, the hydrodynamic stability evaluation received in the "droplet task" is essential to the development of a complete two-phase flow and heat exchange model.

In this paper, the results of modeling of liquid n-decane unsteady temperature evaporated droplet heated by radiation are presented. The simulation results of spontaneous circulation are summarized by Rayleight's criterion. The possibility of spontaneous circulation in a nonisothermal hydrocarbon droplet was evaluated by comparing the results of modelling with experimental data of symmetrical heating of spherical water volume [4].

2. Research methodology

The intensity of heat and mass transfer in the system of gas and dispersed liquid droplets is determined by the processes of transfer inside and outside of the droplets. Heating of the droplets influences a number of phenomenons, including the liquid self-circulation caused by changes of the density. The intensity of natural circulation inside of the liquid droplet can be expressed by the criterion of Rayleight which is calculated by multiplying Grashof's and Prandtl's criteria

[Ra.sub.l] = [Gr.sub.l][Pr.sub.l] (1)

The criterion of Grashof takes into account Archimedes forces acting in the nonisothermal droplet and the Prandtl criterion assesses the peculiarities of liquid physical properties

[Gr.sub.l] = [beta][gl.sup.3]/[v.sup.2][DELTA][T.sub.l]; [Pr.sub.l] = v/a (2)

Before Grashof's and Prandtl's criteria can be calculated, the geometric shape of the analysed volume and temperature distribution in it should be taken into account. The assumption that the spherical droplet is being heated symmetrically has to be made (Fig. 1). In this case, heating a droplet by heat flux from its surroundings is uniform on the whole surface of it and can be expressed by this equation

[q.sup.+.sub.[SIGMA]]([tau]) = [Q.sup.+.sub.[SIGMA]]/A([tau]) [equivalent to] [Q.sup.+.sub.c]([tau]) + [Q.sup.+.sub.r]([tau])/4[pi][R.sup.2]([tau]) (3)

The time dependant temperature field of warming droplet is uniquely determined by the droplet radial coordinates and the time function T(r,[tau]), which must satisfy the condition of the droplet temperature field symmetry center-wise

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (4)

The symmetrically heated droplet temperature range can be defined by the difference of temperature between the surface and center of the droplet

[DELTA][T.sub.l]([tau]) [equivalent to] [T.sub.R]([tau]) - [T.sub.c]([tau]) (5)

In this way, the droplet characteristic dimension in Grashof's criterion is the droplet radius: l [equivalent to] R([tau]).

Physical properties of the liquid in Prandtl's criterion are selected according to the nonisothermal droplet mass average temperature

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (6)

Rayleigh's criterion for the symmetrically heated spherical liquid droplet is described as

[Ra.sub.l,m]([tau]) [equivalent to] [[beta].sub.m]([tau])g[R.sup.3]([tau])/ [v.sub.m]([tau])[a.sub.m]([tau])[DELTA][T.sub.l]([tau]) (7)

The function T([eta],[tau]) of the temperature field of the heated droplet can be expressed as an infinite integral equations line [6]

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (8)

Until functions f ([eta],[tau]) and [f.sub.n]([tau]) are defined, expression (8) is formal and the temperature field of the droplet can not be calculated. The specific shape of the functions f([eta],t) and [f.sub.n]([tau]) is obtained when the initial conditions are formulated and appropriate assumptions about physical properties of the liquid are selected [8, 14, 15].

When defining the function [f.sub.n]([tau]), the following must be taken into account: the droplet surface temperature variation rate, the absorbed radiation heat flux and the initial distribution of droplet's temperature

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (9)

When defining the function f([eta],[tau]), the droplet surface temperature [T.sub.R]([tau]), of the selected control point in time ([tau] [equivalent to] t), and the potential impact of the initial temperature distribution [T.sub.0]([eta]) is evaluated [14]

f ([eta],[tau]) [equivalent to] [T.sub.0]([eta])[[T.sub.R](t)/[T.sub.R,0]] (10)

For the initial state of isothermal droplet f([eta],[tau]) [equivalent to] [T.sub.R]([tau]).

The change of temperature on the surface of the evaporating droplet is determined by the heat flux affecting it. Assuming that the transfer processes are quasistationary, it is stated that the instantaneous speed of unsteady heat and mass transfer processes corresponds to the stationary processes, which take place in the selected instantaneous boundary conditions. Because of that, the function [T.sub.R]([tau]) of the aforementioned droplet surface temperature can be defined on the heat flux balance condition

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (11)

In Eq. (11) the vapour mass flux density on the surface of the droplet is defined by the expression [7]

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (12)

while the gradient of temperature in the droplet is found by differentiating expression (8)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (13)

The local radiation flow in the semitransparent droplet is calculated by [8] methodology. The function [T.sub.R]([tau]) for temperature on the surface of the droplet is determined by iterative calculations. The fastest descent method was used for heat flux balance on the surface functional (11).

For the calculation of droplet's evaporation dynamics, the effect of phase transformations on the droplet surface and the impact caused by thermal expansion of the heated liquid on the change of droplet's size are taken into account by the following equation

[[rho].sub.l]([tau])[d[R.sup.3]([tau])/d[tau]] = [R.sup.3]([tau])[d[[rho].sub.l,m]([tau])/d[tau]] - 3[R.sup.2]([tau])[m.sub.v]([tau]). (14)

2. Research results

Higher saturation temperature is typical to hydrocarbons with longer chain molecules [16]. Until the droplets of such liquid will reach the state of equilibrium conditioned by evaporation, more favorable conditions for higher temperature gradients inside the droplet will be reached [17].

For the analysis discussed in this paper, n-decane, a widely applied fuel in liquid fuelled rockets was used. Warming and evaporation of a droplet of this fuel's are simulated in different cases of heating. Heating conductivity is modelled supposing that the droplets are carried out by a dry airflow without slipping. In case of combined heating, the black body with outside air temperature was stated as radiation source.

Disregarding the method of heating, in the initial stage of intensive droplet warming, a sufficient non-stationary temperature field inside of the droplets is observed (Fig. 1). In Fig. 2 the aforementioned time dependent temperature field is represented by the surface of the droplet, the centre of mass and the average mass temperature.

[FIGURE 1 OMITTED]

[FIGURE 2 OMITTED]

The evaporation time and thermal state of the droplet depends on the size of it (Fig. 3, a). Therefore, real-time scale analysis of heat and mass transfer parameters of the dispersed liquid droplets of various sizes is quite tricky. A dimensionless form is convenient for analysing the unsteady temperature field

T(r,[tau]) [equivalent to] [bar.T]([eta],Fo)

where [bar.T] [equivalent to] T([eta],Fo)/[T.sub.0,l] [eta] [equivalent to] r/R(Fo), Fo [equivalent to] ([a.sub.0]/[R.sup.2.sub.0])[tau].

[FIGURE 3 OMITTED]

When the temperature for the sprayed liquid and its surroundings is known and partial steam pressure is defined, then the function [bar.T]([eta],Fo) is the same for all conductively heated droplets [18]. When time is expressed as a Fourier criterion, temperature of the droplet surface, its centre and the average mass temperature [[bar.T].sub.x](Fo) are independent from the diameter of the droplet (Fig. 3, b).

In the initial stage of droplet evaporation, the different rate of heating on the surface and the central part of the droplet may cause significant differences of temperature between them (Fig. 4). However, at the initial stage of droplet evaporation, the method of heating does not have significant influence on the temperature difference between the surface and the centre of the droplet [DELTA][bar.T] [equivalent to] [[bar.T].sub.R] - [[bar.T].sub.C] (Fig. 5).

Initially temperature in the surface layer changes much more intensively than in the centre (Fig. 4). While heating rate of the droplet surface declines consistently, warming rate of the central layer accelerates until reaching its peak. After that, the temperature difference between the surface and the centre of the droplet begins decreasing. From this moment, thermal conditions inside the droplet are heavily influenced by the method by which it is being heated (Fig. 5). The thermal conditions are distinctive to the size of the droplet when it is being heated conductively and presented in a real time scale (Fig. 6, a), but when time is expressed by the Fourier criterion, its distribution is universal (Fig. 6, b). Therefore, [DELTA][[bar.T].sub.k+r](Fo) deviation from the conductively heated droplet universal [DELTA][[bar.T].sub.k](Fo) curve in the case of combined heating shows the influence of radiation and droplet's size (Fig. 5).

[FIGURE 4 OMITTED]

For smaller than 100 micrometres n-decane droplets the radiation influence is negligible, but for large droplets with a diameter larger than 500 micrometers, the influence of radiation is quite strong. The absorbed radiation flux inside of them causes the overheating of the inner droplet layers and creates preconditions for the occurrence of a negative temperature field gradient. This causes the second [DELTA][bar.T](Fo) peak observed in larger droplets of n-decane (Fig. 5, 4, 5 curves).

[FIGURE 5 OMITTED]

[FIGURE 6 OMITTED]

The change of temperature field causes changes in the liquid droplet density, which brings Archimedes forces into effect. At the initial stage of the evaporation process, the increase of Archimedes forces is highlighted by Rayleight's criterion, which, in this case, reaches its peak (Fig. 7).

[FIGURE 7 OMITTED]

The internal circulation inside the droplets caused by Archimedes forces begins if Rayleigh's criterion reaches a critical value [Ra.sub.kr]. According to experimental results of the natural circulation of water in spherical volume [4] values of [Ra.sub.kr] are highly dependent on heating conditions. In the case of symmetric heating, the critical value of the Rayleigh's criterion for water is about two thousand, but when it is heated asymmetrically it is equal to a few hundred [4]. Diameter of the volume was taken as a characteristic dimension in the criteria when experimental data was summarised [4]. In expression (7) for symmetric heating of a fuel droplet the droplet's radius was taken as characteristic dimension. Therefore, for the thermal state assessment of the n-heptane droplet from the point of view of liquid self-circulation possibility [Ra.sub.kr] of few hundred was stated.

The size of a droplet and the temperature gradient are the key elements that affect the possibility of Archimedes forces to influence heat transfer inside the droplet. While heating rate of the droplet surface layer slows down, and the rate of heating of central layers increases, the temperature field gradient inside of the droplet decreases. The change in the aforementioned factors influences the effect of Archimedes forces to sharply decrease (Fig. 7). Along with rising temperature, the evaporation process intensifies and the size of the droplet decreases rapidly (Fig. 8).

[FIGURE 8 OMITTED]

3. Conclusions

1. The thermal state of the liquid hydrocarbon droplet changes during evaporation. The peculiarities of this change are influenced by the method of droplet heating, droplet diameter and liquid fuel molecular weight.

2. A numerical study of the n-decane droplet heating and evaporation showed that in the case of heating the droplet conductively, the thermal state change expressed using Fourier's criterion is universal and can be stated with typical droplet surface, centre and average mass temperature describing curves.

3. When time is expressed by Fourier's criterion, other conductively heated droplet heat and mass transfer parameters can be described and presented in a non-dimensional form. In the case of combined droplet heating, a deviation in the heat and mass transfer parameters regarding universal curves obtained from conductive heating helps evaluate how the radiation affects them.

4. According to the results of the provided numerical investigation of evaporating pure hydrocarbon n-decane droplet thermal conditions, it is clear that Archimedes forces in high molecular weight hydrocarbon droplets caused by internal temperature gradient are insufficient to cause spontaneous fluid circulation in the common for thermal technologies range of droplet size.

5. The temperature difference between the surroundings of the droplet and the droplet itself affects the thermal and hydrodynamic state of the evaporating hydrocarbon droplet. Therefore it is important to investigate the process of droplet evaporation in different initial temperature conditions of the droplet and its surroundings.

Nomenclature

A--droplet surface area, [m.sup.2]; a--thermal diffusivity, [m.sup.2]/s; [B.sub.M]--Spalding mass transfer number; [B.sub.T]--Spalding heat transfer number; [c.sub.p]--specific heat, J/(kgK); Fo--Fourier number; Gr--Grashof number; L--latent heat of evaporation, J/kg; l--characteristic dimension, m--vapour mass flux density, kg/(s x [m.sup.2]); n--number of the term in infinite sum; q--heat flux density, W/[m.sup.2]; p--pressure, Pa; R--radius of the droplet, m; r--coordinate of the droplet, m; Ra--Rayleigh number; Nu--Nusselt number; T--temperature, K; [beta]--thermal volumetric expansion coefficient, [K.sup.-1]; [eta]--droplets dimensionless coordinate; [lambda]--thermal conductivity, W/(m K); [rho]--density, kg/[m.sup.3]; [tau]--time, s.

Subscript: c--heating by convection; g--gas; f--phase transformation; k--conductive; l--liquid; m--mass average; r--heating by radiation; R--droplet surface; 0--initial state; [infinity]--far from a droplet.

Superscript: m--modified; + --external side of droplet surface; - --internal side of droplet surface.

Received May 11, 2011

Accepted April 12, 2012

References

[1.] Semenov, S.; Starov, V.M.; Velargde, M.G.; Rubio, R.G. 2011. Droplets evaporation: Problems and solution, Eur. Phys. J. Special Topic 197, 265-278. http://dx.doi.org/10.1140/epjst/e2011-01468-1.

[2.] Sazhin S.S. 2006. Advanced models of fuel droplet heating and evaporation, Progress in Energy and Combustion Science 32: 162-214. http://dx.doi.org/10.1016/j.pecs.2005.11.001.

[3.] Merouane, H.; Bounif, A.; Abidat, M. 2011. Simulation of evaporating isopropyl alcohol droplets injected into a turbulent flow, Mechanika 17(6): 594-600. http://dx.doi.org/10.5755/j01.mech.17.6.1001.

[4.] Gershuni G.Z.; Zhukhovitsky E.M. 1972. Convective stability of an incompressible liquid, Moscow: Nauka. 392 p. (in Russian).

[5.] Sirignano, W.A. 1999. Fluid dynamics and transport of droplets and sprays, Cambridge: Cambridge University Press. 309 p. http://dx.doi.org/10.1017/CBO9780511529566.

[6.] Renksizbulut, M.; Bussmann, M.; Li, X. 1992. A droplet vaporization model for spray calculations, Part. Part. Syst. Charact. 9: 59-65. http://dx.doi.org/10.1002/ppsc.19920090110.

[7.] Kuzikovskij, A.V. 1970. Dynamic of spherical particle in powerful optical field, Izv. VUZ Fizika 5: 89-94 (in Russian).

[8.] Miliauskas, G. 2001. Regularities of unsteady radiative-conductive heat transfer in evaporating semitransparent liquid droplets, Int. J. Heat Mass Transfer 44: 785-798. http://dx.doi.org/10.1016/S0017-9310(00)00127-7.

[9.] Tseng, C.C.; Viskanta, R. 2006. Enhancement of water droplet evaporation by radiation absorption, Fire Safety J. 41: 236-247. http://dx.doi.org/10.1016/j.firesaf.2006.01.001.

[10.] Lage, P.L.C.; Rangel, R.H. 1993. Single droplet vaporization including thermal radiation absorption, J. of Thermophysics and Heat Transfer 7: 502-509. http://dx.doi.org/10.2514/3.446.

[11.] Tuntomo, A.; Tien, C.L.; Park, S.H. 1991. Internal distribution of radiant absorption in a spherical particle, J. Heat Transfer 113: 407-412. http://dx.doi.org/10.m5/L2910576.

[12.] Sinkunas, S.; Kiela, A. 2011. Heat transfer for film condensation of vapour, Mechanika 17(2): 144-148. http://dx.doi.org/10.5755/j01.mech.17.2.329.

[13.] Saffari, H.; Dalir, N. 2011. Effect of friction stress of droplets with film on prediction of pressure changes in condensing tubes, Mechanika 17(1): 50-56. http://dx.doi.org/10.5755/j01.mech.17.L203.

[14.] Miliauskas, G. 2003. Interaction of the transfer processes in semitransparent liquid droplets, Int. J. Heat Mass Transfer 46: 4119-4138. http://dx.doi.org/10.1016/S0017-9310(03)00231-X.

[15.] Miliauskas, G.; Sinkunas, S. 2008. Influence of droplets heating manner to the sprayed liquid thermal state, CHT'08: proceedings of the International Symposium of the Advances in Computational Heat Transfer, Moroco, pp. 1-17.

[16.] Reid, R.C.; Prausnitz, J.M.; Sherwood, T.K. 1977. The properties of gases and liquid, New York: McGraw-Hill Book Company, 549 p.

[17.] Miliauskas, G.; Sinkunas, S.; Talubinskas, J.; Sinkunas, K. 2010. Peculiarities of hydrodynamics in the evaporation of hydrocarbon droplets, Advances in Fluid Mechanics VIII: WIT Transactions on Engineering Sciences 69, Boston, 283-292.

[18.] Miliauskas, G.; Sabanas, V.; Bankauskas, R.; Miliauskas, G.; Sankauskaite, V. 2008. The peculiarities of sprayed liquid's thermal state change, as droplets are heated by conduction, Int. J. of Heat and Mass Transfer 51: 4145-4160. http://dx.doi.org/10.1016/j.ijheatmasstransfer.2008.01.0 16.

G. Miliauskas, Kaunas university of technology, K.Donelaicio 20, 44239 Kaunas, Lithuania, E-mail: gimil@ktu.lt

J. Talubinskas, Kaunas university of technology, K.Donelaicio 20, 44239 Kaunas, Lithuania, E-mail: _julius. talubinskas@gmail.com

A. Adomavicius, Kaunas university of technology, K.Donelaicio 20, 44239 Kaunas, Lithuania, E-mail: arado@ktu.lt

E. Puida, Kaunas university of technology, K.Donelaicio 20, Kaunas, LT-44239, Lithuania, E-mail: epuida@ktu. lt

http://dx.doi.org/ 10.5755/j01.mech.18.2.1573
联系我们|关于我们|网站声明
国家哲学社会科学文献中心版权所有