2. Methods
An integrated parallel plate microreactor configuration is simulated. The oxidation and reforming reactions are carried out in alternate channels separated by walls. The stoichiometric methanol-air flow in the oxidation channel is co-current with respect to the methanol-steam flow in the reforming channel. Both the fluid streams enter the device at room temperature and exit at atmospheric pressure. Using the inherent symmetry of the geometry, only one-half of each channel and the connecting wall are simulated. The fluid flow as well as the heat and mass transfer equations are solved in the fluid phase and the energy balance is explicitly accounted for in the solid wall. A steady-state solution of the problem is sought using computational fluid dynamics. More specifically, computational fluid dynamics is used to solve the steady-state continuity, momentum, energy, and species conservation equations in the fluid phase and the heat equation in the solid phase using a finite volume approach. An adaptive meshing scheme is used for the discretization of the differential equations. The computational mesh is initialized with 200 axial nodes, 200 radial nodes for the oxidation channel and the wall sections, and 200 radial nodes for the reforming channel. This discretization translates into a total of about 80000 nodes. This initial mesh is adapted and refined during a calculation to increase the accuracy of the solution in regions of high gradients. Specifically, additional nodes are introduced to refine the mesh using the tools built in the computational software. Specifically, mesh refinement before achieving complete convergence reduces the computational effort, but a too early refinement, namely in a few iterations, may lead to refinement in wrong regions. Such an adaptive meshing strategy, starting with a relatively coarse initial mesh followed by refinement in regions of large gradients, achieves an adequate balance between accuracy and computational effort.
The fluid density is calculated using the ideal gas law. The individual properties of various gaseous species, such as thermal conductivity, are calculated using the kinetic theory of gases, whereas the specific heats are determined as a function of temperature using polynomial fits from the thermodynamic database. Mixture properties, such as specific heat and thermal conductivity, are calculated from pure component values based on the mass-fraction weighted mixing law. Binary species diffusivities are determined and then are used to calculate the multicomponent mixture diffusivities. For the solid wall, a constant specific heat and an isotropic thermal conductivity are specified. Given that material conductivity varies with temperature and more importantly with the material chosen, simulations are carried out over a wide range of conductivities. To permit modeling of fluid flow and related transport phenomena in industrial equipment and processes, various useful features are provided. These include porous media, lumped parameter streamwise-periodic flow and heat transfer, swirl, and moving reference frame models. The moving reference frame family of models includes the ability to model single or multiple reference frames. A time-accurate sliding mesh method, useful for modeling multiple stages in turbomachinery applications, for example, is also provided, along with the mixing plane model for computing time-averaged flow fields. Another very useful group of models is the set of free surface and multiphase flow models. These can be used for analysis of gas-liquid, gas-solid, liquid-solid, and gas-liquid-solid flows. For all flows, conservation equations for mass and momentum are solved. For flows involving heat transfer or compressibility, an additional equation for energy conservation is solved. For flows involving species mixing or reactions, a species conservation equation is solved or, if the non-premixed combustion model is used, conservation equations for the mixture fraction and its variance are solved. Additional transport equations are also solved when the flow is turbulent.
The flow of thermal energy from matter occupying one region in space to matter occupying a different region in space is known as heat transfer [29, 30]. Heat transfer can occur by three main methods: conduction, convection, and radiation. The net transport of energy at inlets consists of both the convection and diffusion components. The convection component is fixed by the inlet temperature specified. The diffusion component, however, depends on the gradient of the computed temperature field [31, 32]. Consequently, the diffusion component and therefore the net inlet transport is not specified a priori. When heat is added to a fluid and the fluid density varies with temperature, a flow can be induced due to the force of gravity acting on the density variations. Such buoyancy-driven flows are termed natural-convection or mixed-convection flows. Multiple simultaneous chemical reactions can be modeled, with reactions occurring in the bulk phase and on wall surfaces, and in the porous region [33, 34]. For gas-phase reactions, the reaction rate is defined on a volumetric basis and the rate of creation and destruction of chemical species becomes a source term in the species conservation equations. For surface reactions, the rate of adsorption and desorption is governed by both chemical kinetics and diffusion to and from the surface [35, 36]. Wall surface reactions consequently create sources and sinks of chemical species in the gas phase, as well as on the reacting surface. Reactions at surfaces change gas-phase, surface-adsorbed and bulk species. On reacting surfaces, the mass flux of each gas specie due to diffusion and convection to and from the surface is balanced with its rate of consumption and production on the surface. The model boundary conditions are described next. Boundary conditions consist of flow inlets and exit boundaries and wall boundaries. Symmetry boundary conditions are used since the physical geometry of interest, and the expected pattern of the flow and thermal solution, have mirror symmetry. The axis boundary type is used as the centerline of an axisymmetric geometry. Both gases enter the channels at room temperature with a uniform, flat flow velocity. The reactor exits are held at a fixed pressure and the normal gradients of species and temperature, with respect to the direction of the flow, are set to zero. Symmetry boundary condition is applied at the centerline of both channels, implying a zero normal velocity and zero normal gradients of all variables. No-slip boundary condition is applied at each wall-fluid interface. Overall, the device is adiabatic, namely no heat losses occur through the side walls. Continuity in temperature and heat flux is applied at the fluid-solid interfaces. Neither heat-transfer nor mass-transfer correlations are employed since detailed transport within the solid and fluid phases is explicitly accounted for.
The pressure-based solver traditionally has been used for incompressible and mildly compressible flows [37, 38]. The density-based approach, on the other hand, is originally designed for high-speed compressible flows [39, 40]. Both approaches are now applicable to a broad range of flows, from incompressible to highly compressible, but the origins of the density-based formulation may give it an accuracy, namely shock resolution, advantage over the pressure-based solver for high-speed compressible flows [41, 42]. Two formulations exist under the density-based solver: implicit and explicit. The density-based explicit and implicit formulations solve the equations for additional scalars, for example, turbulence or radiation quantities, sequentially [43, 44]. The implicit and explicit density-based formulations differ in the way that they linearize the coupled equations. In both methods the velocity field is obtained from the momentum equations [45, 46]. In the density-based approach, the continuity equation is used to obtain the density field while the pressure field is determined from the equation of state. On the other hand, in the pressure-based approach, the pressure field is extracted by solving a pressure or pressure correction equation which is obtained by manipulating continuity and momentum equations. In both cases a control-volume-based technique is used that consists of division of the domain into discrete control volumes using a computational grid, integration of the governing equations on the individual control volumes to construct algebraic equations for the discrete dependent variables, such as velocities, pressure, temperature, and conserved scalars, and linearization of the discretized equations and solution of the resultant linear equation system to yield updated values of the dependent variables. The two numerical methods employ a similar discretization process, but the approach used to linearize and solve the discretized equations is different. The full computational fluid dynamics problem is solved via a segregated solver using an under-relaxation method. Convergence of the solution is monitored through the residuals of the governing equations and the norm of successive iterations of the solution. The coupling of the heat equation in the wall and the reacting flow equations makes the problem stiff due to the disparity in thermal conductivity between the gases and the wall. Typical computational fluid dynamics simulation times varies from about several hours to a few days on a single processor depending on the stiffness of the problem. Natural parameter continuation is employed to study the effect of various operating parameters.