Next Article in Journal
Resonance Suppression Method Based on Hybrid Damping Linear Active Disturbance Rejection Control for Multi-Parallel Converters
Previous Article in Journal
Enhanced Hydrolysis of Carbonyl Sulfide in Coking Oven Gas Utilizing an Efficient Ca-Ba-γ-Al2O3 Catalyst
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Enhanced Transverse Dispersion in 3D-Printed Logpile Structures: A Comparative Analysis of Stacking Configurations

by
Leon R. S. Rosseau
,
Martijn A. A. van Aarle
,
Egbert van Laer
,
Ivo Roghair
and
Martin van Sint Annaland
*
Department of Chemical Engineering and Chemistry, Eindhoven University of Technology, P.O. Box 513, 5600 MB Eindhoven, The Netherlands
*
Author to whom correspondence should be addressed.
Processes 2024, 12(10), 2151; https://doi.org/10.3390/pr12102151
Submission received: 9 September 2024 / Revised: 28 September 2024 / Accepted: 29 September 2024 / Published: 2 October 2024
(This article belongs to the Section Chemical Processes and Systems)

Abstract

:
Three-dimensionally printed logpile structures have demonstrated the tunability of the transverse dispersion behavior, which is relevant in the context of chemical reactor design. The current modeling study aims to further investigate the trade-offs in such structures, extending the range of geometries investigated and addressing the limitations associated with the pseudo-2D nature of previous experiments. The relative transverse dispersion coefficient and pressure drop were determined using computational fluid dynamics simulations in OpenFOAM for structures with different stacking configurations, porosities and scaling of the structures’ unit cell along the secondary transverse axis. The latter could not be varied in previous experiments, but the current results demonstrate that this limitation suppresses vortex shedding in structures with high porosity. These vortices significantly enhance the transverse dispersion. By using a staggered stacking configuration on both transverse axes, an earlier onset of this phenomenon could be realized. Importantly, operation in this regime could be achieved without an equivalent increase in pressure drop, offering a favorable operating trade-off. The findings demonstrate that at low Reynolds numbers, the studied structures consistently outperform randomly packed beds of spheres, highlighting their potential for chemical process intensification.

1. Introduction

Adequate wall-to-bed heat transfer is of utmost importance for the optimal performance of catalytic reactors, as deviations from the intended operating temperature can lead to lower conversion and potentially lower product selectivity [1,2]. This decreases the overall process efficiency, as downstream separation requirements are increased.
In commonly used fixed bed reactors, heat transfer occurs through a combined process of solid-phase conduction and fluid-phase dispersion. However, the low thermal conductivity of catalyst powder and the presence of point contacts between particles in randomly packed bed reactors create significant resistance to heat transfer [3,4,5]. Thus, solid-phase conduction alone is insufficient for proper heat management. Fluid-phase heat transfer is a combined process of conduction and convection. The conductive contribution is governed by the thermal conductivity of the fluid, which is generally low for gases. The convective contribution, on the other hand, can be related to the transverse dispersion coefficient. This coefficient is a measure of the rate at which fluid is transported perpendicular to the flow direction, lumping the effects of convective flow phenomena on mass transport [6]. The transverse dispersion in randomly packed beds of spherical particles has been investigated extensively, revealing a linear relationship between the relative transverse dispersion coefficient and the Péclet number for gas-phase systems, as shown in Equation (1) [7]. The definition of the Péclet number is shown in Equation (2). Combining the high packing density of randomly packed beds with the relatively poor scaling with the Péclet number results in relatively low rates of transverse dispersion in these systems at typical operating conditions.
D t D m = 1 τ + 1 12 P e m
P e m = U d p D m ϵ
In randomly packed beds of spherical particles, the combined process of solid-phase conduction and fluid-phase dispersion does not sufficiently facilitate heat transfer for strongly endothermic or exothermic reactions. To achieve the necessary heat management in such systems, a multi-tubular configuration with low-diameter tubes (in the order of several inches) is typically employed, which is relatively expensive [8,9]. However, tailored catalyst geometries can enhance the transverse dispersion and offer a promising solution, potentially leading to reduced capital and operating expenditures thanks to improved heat management [10].
With the advent of additive manufacturing technologies, it is now possible to shape catalyst structures via 3D printing. Currently, the most mature in this context are the so-called logpile structures, produced via direct ink writing (DIW). This process involves extruding a viscous catalyst slurry through a circular nozzle in a programmed pattern [11,12]. As the name suggests, logpile structures consist of stacked cylinders arranged in alternating layers, forming a lattice-like structure. Tailored structures can be produced by adjusting the nozzle size, cylinder spacing to control porosity and stacking configuration [13,14,15,16]. Compared to randomly packed beds of catalyst particles, logpile structures offer more design freedom to navigate the reactor design trade-offs. At the same time, the catalyst holdup is significantly higher compared to washcoated reactor internals such as monoliths or foams. When using an aligned stacking configuration (e.g., the same pattern on all stacked layers), these structures resemble honeycomb monoliths. However, due to the stacking method employed, there are gaps in the transverse direction which may introduce transverse dispersion. Alternatively, a staggered configuration can be utilized, where every other stacked layer is offset from the previous one. This configuration offers a more tortuous path for the fluid and promotes vortex shedding, contrasting the aligned configuration where vortices may dissipate when encountering downstream cylinders (depending on the exact geometry) [17,18]. Because of these characteristics, the staggered configuration is expected to further enhance the transverse dispersion coefficient.
In a recent experimental work from our group, various designs of pseudo-2D logpile geometries were studied using a non-invasive infrared transmission method [19]. As hypothesized, it was found that there was some tunability in the transverse dispersion coefficient, specifically by varying the spacing in staggered structures. This behavior was correlated to the Péclet number following a similar relationship as described in Equation (1). Interestingly, the observed transverse dispersion coefficient exhibited a superlinear dependence, being a function of P e m 1.67 . However, a lower pre-factor was fitted, resulting in transverse dispersion coefficients lower than those in packed beds of monodisperse spheres. The underlying mechanisms behind the observed behavior still lack fundamental explanations. This particularly relates to the influence of the pseudo-2D nature of the experiments. The pseudo-2D character was an artifact of the experimental setup, as a flat column was required to capture the relevant phenomena in the structure with the infrared camera. As a result of this, the column depth was fixed and could not be scaled as the unit cell dimensions were varied.
The objective of this study is to investigate the transverse dispersion behavior and pressure drop in staggered logpile structures in more detail through computational fluid dynamics (CFD) simulations using OpenFOAM. Specifically, the impact of different stacking configurations in the secondary transverse direction is investigated. Although this design parameter could not be varied in the experimental setup, it plays an important role in translating results from pseudo-2D experiments into practical data for actual (3D) structured reactors and is of special interest since 3D printing allows for independent variations in both transverse directions. Furthermore, due to the relatively low expected pressure drop in these structures, it was not possible to obtain statistically significant experimental data for the pressure drop. However, it is essential to quantitatively assess the influence of catalyst geometry on the pressure drop for effective reactor design considerations. The current scope does not include the aligned (or straight) stacking configuration, as previous experimental efforts have indicated a relatively low transverse dispersion coefficient in such structures, which limits their practical applicability [19].
This paper is structured as follows. In Section 2, the design variables of the structures under study are introduced and the simulation setup is described. This is followed by considerations on the extraction of a representative transverse dispersion coefficient from tracer concentration profiles. In Section 3, a comparison is made between the transverse dispersion coefficients obtained from simulations and experiments from the literature. Subsequently, the influence of scaling of the structures’ depth, boundary conditions and stacking configuration on the transverse dispersion coefficient and pressure drop is investigated and discussed, alongside a comparison to a packed bed of monodisperse spheres. Finally, a conclusion is presented in Section 4.

2. Methods

2.1. Geometries

The structures examined in this study are based on a grid consisting of 10 by 13 cylinders with a diameter of 1.5   m m in the xz plane. Fully 3D structures include additional cylinders in the yz plane. In the simulated domain, the z-direction represents the axial direction (i.e., the flow direction), the x-direction represents the transverse direction and the y-direction represents a single logpile unit cell of varying configuration in the secondary transverse direction. Staggered structures are investigated, where every odd axial layer is offset by 50% in the transverse direction (i.e., the xz plane). This leads to two half-cylinders on the either side of the domain in these layers. Consistent with previous work, an axial stacking offset of 80% is assumed, which is common for structures produced by DIW [19]. As a result, all structures have an axial length of 3.12   c m . The basic dimensions of an exemplary structure are visualized in Figure 1.
Three main geometrical variations are considered. The first variation is the spacing between cylinders in the xz plane, with values of 0.75   m m , 1.5   m m or 3 m m (resulting in domain widths ranging from 2.25   c m to 4.5   c m ). This range of spacing values is considered reasonable within the generally expected constraints of the DIW technology [12]. These variations are denoted as S0.5, S1 and S2, respectively, with the number representing the ratio between the spacing and cylinder diameter. The second variation is the stacking configuration of the single unit cell in the yz plane. Four variations are studied, as visualized in Figure 1: 2D, side, center and stag. The 2D geometry does not have additional cylinders placed in the yz plane. The side variation features two half-cylinders on each side of every axial layer, resembling an aligned stacking configuration that closely corresponds to the pseudo-2D structures studied experimentally. The center variation, also similar to the aligned stacking configuration, has centered whole cylinders on each axial layer instead of half-cylinders on either side. Both side and center variations represent the same structure in the absence of bounding walls, but such walls were present in previous experimental work. Finally, the stag geometry represents staggered stacking in the yz plane, alternating between half-cylinders on each side and a centered cylinder. The final geometrical variation is the dimension of the yz unit cell, which can either be a fixed width (3 m m , as in the experiments) or scaled according to the spacing ( 2.25   m m for S0.5 and 4.5   m m for S2), also shown in Figure 1. This is indicated by either the fixed or scaled suffix in the geometry identifier. The variations and resulting porosities are presented in Table 1.

2.2. Simulation Setup

To study the hydrodynamics, the rhoReactingFoam solver from the OpenFOAM-v2106 software (OpenCFD Ltd., Bracknell, UK) was employed. This transient solver is based on the pressure implicit method with the splitting of the operator for pressure-linked equations (PIMPLE) and intended for the simulation of compressible multicomponent flows. The conservation of (species) mass and momentum is described by Equations (3)–(5), using the general Newtonian form of the stress tensor and the ideal gas law as an equation of state [20]. Since the focus of this study is solely on the transverse dispersion and pressure drop induced by the various geometries, the solid phase is not explicitly modeled other than imposing a no-slip boundary condition at its wall. Validation of the rhoReactingFoam solver is provided in Appendix A.
ρ t + · ρ u = 0
ρ u t + · ρ u u = p · τ
ρ Y i t + · ρ u Y i = · D i Y i
To study the transverse dispersion in the different geometries, transient simulations were conducted with a continuous tracer injection of fixed width at the inlet. The tracer gas properties were set to match those of nitrogen, which was the background gas. The diffusion coefficient of the tracer in the background was then equal to the molecular diffusion coefficient of nitrogen. The simulations were isothermal, maintained at a temperature of 300 K , and an outlet pressure of 1 atm was applied. By varying the superficial velocity, thirteen variations of the Péclet number were studied, ranging from approximately 10 to 175. It should be noted that, due to the different porosities, a comparison between different structures at the same Péclet number is by definition at different Reynolds numbers, and this may influence vortex shedding phenomena. Across all structures, Reynolds numbers ranging from approximately 4 to 142 were studied. A no-slip velocity boundary condition was imposed on the cylinders, while either a slip or a no-slip velocity boundary condition was applied to the walls. This boundary condition was varied to study the influence of bounding walls on the transverse dispersion, with the no-slip condition accounting for the bounding walls that were present in the experimental setup. The time step was determined using the Courant number, with a maximum value of 0.8. The simulated time was adjusted to ensure that steady state was reached, as discussed in Section 2.4. At higher velocities, temporal fluctuations were observed due to vortex shedding, necessitating the use of a transient solver and subsequent time-averaging of the results to obtain accurate performance estimates. Although these fluctuations could potentially be enhanced by turbulent effects, it was observed that the use of a Reynolds-averaged Navier Stokes (RANS) k ω model did not yield significantly different outcomes compared to the laminar case setup, despite the associated increase in computational cost. These considerations can be found in Appendix B.
OpenFOAM meshes were generated by the blockMesh and snappyHexMesh utilities. Initially, blockMesh was used to generate a background mesh of cubic cells. Subsequently, using a stereolithography (.stl) file containing the geometry as generated by the Blender 3.0 software (Blender Foundation, Amsterdam, The Netherlands), snappyHexMesh removes the geometry from the background mesh. In this work, snappyHexMesh was used with level two cell splitting to locally refine the mesh around the cylinders. This local refinement accurately captures the boundary layers near the structure but increases both the number of cells and the non-orthogonality significantly, necessitating the use of corrected surface-normal gradient schemes.
The refinement of the background mesh is the key variable in balancing accuracy and computational cost. A suitable value was determined through a grid size study conducted on the S12D structure, encompassing both the lower and upper limits of the Péclet number range investigated in this study. The results of the grid size study, presented in Figure 2, show the relative error in the calculated transverse dispersion coefficient and the relative (computational) cost, normalized with respect to the finest refinement. Notably, the error is relatively large for P e m = 169 , due to the relatively large coefficient of variation (COV, the standard deviation divided by the mean value) associated with temporal fluctuations. This can be explained by the vortex shedding that occurs at these conditions ( R e p = 128 ). For reference, these data exhibit a COV of 0.08, contrasting the value in the order of 10 5 for P e m = 8 . Further analysis of these uncertainties is presented in Section 2.4, but considering the high COV, the 150 μ m refinement is considered adequate for the current purpose. A visualization of the mesh refinement is shown in Appendix C.
The simulations in this work were conducted on the Dutch supercomputer ‘Snellius’. Running the S1stag case at the highest Péclet number for a single residence time required 1253 s of wall clock time using 256 cores, equal to approximately 90 core hours on AMD EPYC 7H12 CPUs. In contrast, running the S12D case at the lowest Péclet number for a single residence time only took 32 s of wall clock time on 32 cores, equal to approximately 0.28 core hours. The significant difference in computational time can mainly be attributed to the substantial difference in the number of cells between the two cases: 6,821,440 for the former and merely 167,060 for the latter.

2.3. Data Evaluation

Using the built-in OpenFOAM sampling tool, transverse concentration, flow and pressure profiles were extracted at various axial positions. To describe the concentration profile on the reactor scale rather than the individual cells in the mesh, the simplified convection–diffusion equation in Equation (6) was used. This equation was derived by applying the same assumptions as in previous experimental work [19]:
  • The system is at steady state and temporal fluctuations due to vortex shedding are averaged. The resolution of averaging for representative estimates is discussed in Section 2.4.
  • A uniform axial velocity profile is assumed from the reactor scale perspective. The influence of transverse velocity profiles is not considered explicitly, as their effect is incorporated into the dispersion coefficient. It is assumed that the dispersion coefficient remains constant throughout the structure. This assumption is valid in the absence of hydrodynamic entrance effects, which is discussed in Section 2.5.
  • The concentration profiles in the secondary transverse direction (the y-direction) were condensed into a single concentration profile in the primary transverse direction (the x-direction). This is justified since the tracer is injected along the entire length of secondary transverse axis. Furthermore, it was observed that concentration gradients along this axis were negligible compared to the profiles obtained along the x-axis.
u C z = D t 2 C x 2
As in previous experimental work, the appropriate analytical solution in Equation (7), with the corresponding standard deviation given by Equation (8), was derived by applying the boundary conditions specified in Table 2 [21]. These conditions specify finite values for the injection width and tracer concentration at the inlet and impose a zero-gradient condition on the side walls so that mass is contained within the domain boundaries.
C = C 0 2 n = erf Δ x 2 n L + x σ 2 + erf Δ x + 2 n L x σ 2
σ = 2 D t z u
The dispersion coefficient can be obtained by rewriting Equation (8) for D t . However, when considering the hydrodynamic entrance length, it is important to account for the adjusted reference point. For this purpose, the dispersion of a tracer in the entrance zone is first calculated by fitting a reference standard deviation using the concentration profile sampled after the entrance zone. This is then subtracted from the fitted standard deviation at the intended outlet position, isolating the transverse dispersion that occurred in the region bounded by the inlet and outlet. This is shown in Equation (9).
D t = σ structure 2 u 2 z = σ outlet 2 σ inlet 2 U 2 z ϵ
Finally, to determine the pressure drop, the pressure was sampled on the transverse plane at specified inlet and outlet positions. The obtained profiles were then averaged, after which the outlet pressure was subtracted from the inlet pressure to obtain the pressure drop.

2.4. Steady State and Time-Averaging

As mentioned previously, a transient solver was used in this study to accurately capture temporal fluctuations due to vortex shedding at higher Reynolds numbers. However, before such data can be obtained, a (pseudo-)steady state must be reached. It can be expected that this would take at least one residence time, as initially there is no tracer present in the domain, but settling effects may require additional time. In Figure 3, the calculated relative transverse dispersion coefficient is plotted as a function of the dimensionless time (the time divided by the residence time) for the S22D structure with P e m = 169 (the highest value currently under investigation). The plot demonstrates a sharp decrease in the relative transverse dispersion coefficient until t * = 2 , after which the profile stabilizes, exhibiting a repetitive pattern characteristic of a pseudo-steady state after approximately t * = 4 . Hence, for all subsequent simulations, it was assumed that the (pseudo-)steady state was achieved after four residence times had elapsed.
Once the steady state was reached, the transverse dispersion coefficient and the COV of time-averaging were determined based on a discrete number of time steps sampled during a final residence time. These parameters were evaluated for the S12D structure with P e m = 169 , considering the entire structure. By using either 5, 10 or 20 time steps during the final residence time, a relative transverse dispersion coefficient of 11.48, 11.53 or 11.55 was obtained, respectively. Regardless of sampling frequency, a COV of approximately 0.032 was obtained. Considering the increasing memory requirements associated with higher sampling frequencies, a sampling frequency of 10 was deemed acceptable.

2.5. Hydrodynamic Entrance Effects

It is important to distinguish between bulk structure transverse dispersion and transverse dispersion induced by the settling of the velocity profile at the inlet and outlet of the structure. To quantify the hydrodynamic entrance length, the development of the velocity profile along the axial coordinate is investigated in structures with additional axial layers. The settling of the velocity profile is evaluated by calculating the mean absolute error (MAE) between subsequent uneven layers, using Equation (10). Figure 4 illustrates the development of this error along the axial coordinate for six exemplary 2D geometries (in the first 4 c m of the structure). To provide further context, both the COV of time-averaging and the pressure drop are presented in Table 3.
MAE = 1 N u ( x , z ) u ( x , z + 1 ) u ( x , z + 1 )
For all structures with P e m = 8 , it is evident that a steady error is achieved within less than 1 c m , or four axial layers. The relative positions of the different curves can be attributed to the corresponding pressure drops, as shown in Table 3. The pressure drop influences the velocity magnitude, causing a deviation in each layer compared to the preceding layer. The curves with P e m = 169 exhibit higher error values and a larger degree of variation as a function of the axial position. The higher pressure drop at these conditions accounts for the relatively high error, and the larger degree of variation can be explained by the high COV of time-averaging. Taking into account these influences, it can be seen that a relatively stable error is also present in all of these simulations after approximately 2 c m , or eight axial layers. This value is in agreement with entrance length estimates from the literature for similar structures [22,23].
To verify the observed settling behavior, the relative transverse dispersion coefficient was calculated for different points of reference. Figure 5 shows this behavior as a function of the relative injection width for an exemplary S12D structure. It is evident that after eight layers, the calculated relative transverse dispersion coefficient becomes independent of axial position. Thus, it can be concluded that using eight layers as an inlet zone, and excluding the outlet but instead using the penultimate layer, allows for accurate determination of the bulk transverse dispersion coefficient. Additionally, employing an injection width at least twice as large as the cylinder size yields a representative solution.

3. Results

Figure 6 presents a comparison between sidefixed geometry simulations with no-slip velocity boundary condition applied to the walls and correlations derived from experimental data from previous work [19]. The results from the two sets of data do not fully agree, showing increasing deviations for larger spacing values and higher experimental values in all cases. The mean absolute percentage error (MAPE) values associated with the comparisons are 7.06%, 11.4% and 23.7% for S0.5, S1 and S2, respectively. It is important to consider these values in the context of the MAPE obtained for fitting the experimental correlation, which was 3.74%. For all simulations in Figure 6, the COV of time-averaging was below 0.005, which is negligible, and therefore, error bars were not included. To eliminate the possible influence of the different diffusion coefficients, the simulations of the S1 geometry were repeated with the molecular diffusion coefficient of propane in nitrogen instead of the self-diffusion coefficient of nitrogen. The results were very similar, with a MAPE of only 0.79% across the range of P e m shown in Figure 6. This confirms that the molecular diffusion coefficient has minimal impact on the relationship between the relative transverse dispersion coefficient and the Péclet number, consistent with findings in the literature [7,24].
The trends observed in Figure 6 reveal a significant difference in the slopes of the different data sets. The superlinear relationship exhibited by the experimental correlation is typically associated with a transitional flow regime [7]. However, simulations using a turbulent model did not show statistically significant differences (as mentioned in Appendix B), and the particle Reynolds number is relatively low (up to 95). Potentially, the observed discrepancy from the trends in the experimental data could be attributed to non-idealities resulting from the printing process and the roughness of the structures. The 3D printing of polymeric structures with fused deposition modeling introduces this roughness due to the layered nature of the printing process. It was established by Ahn et al. [25] that the roughness is in the order of tens of micrometers for features in the order of hundreds of micrometers, which is relatively large. It can be hypothesized that this roughness induces an earlier onset of the laminar–turbulent transition than could be expected based solely on the value of the particle Reynolds number [26,27,28]. Since the polymeric structures used in previous experimental work were models for actual printed catalyst, and since the roughness of such structures is likely not representative of DIW-structured catalyst, no further attempts were made to precisely replicate the experimentally observed trends in the simulations [29,30]. The explanation for the observed deviations is supported by the increasing discrepancy as a function of the spacing, as the combination of increased hydraulic diameter and roughness-induced mixing may enable an earlier onset of a transitional regime. To conclude, while it cannot be substantiated that roughness effects are the cause of the observed discrepancy, it is potentially a relevant factor that could be evaluated in future work.
Two additional factors may contribute to the observed discrepancy. Firstly, the current method accounts for hydrodynamic entrance and exit effects that influence the observed dispersion coefficient. However, since local measurement of concentration profiles within the structure was not possible experimentally, the obtained experimental data represent an averaged value for the entire structure. Secondly, in the experimental setup, the structure was placed downstream of the tracer injection, resulting in a disperse tracer plume entering the structure. Although this was accounted for in the analysis, the potentially non-ideal tracer injection may lead to slightly different behavior compared to the simulated (ideal) injections.
Next, the influence of the pseudo-2D bounding walls that were present in the experimental setup was investigated, by changing the velocity boundary condition on these walls. Specifically, this analysis was performed for the side stacking variation with fixed width, which serves as an intermediate step between the experimental comparison and the complete results. The comparison between the results with slip and no-slip velocity boundary conditions is shown in Figure 7. The coefficient of variation for time-averaging was in the order of 10 3 for the S2 geometries and in the order of 10 4 for the other spacings. These values are negligible, and therefore, error bars were not included.
In Figure 7, it can be observed that the differences between the slip and the no-slip velocity boundary conditions are relatively small, with slightly larger variations for S2. To further support this observation, a parity plot relating the slip results to the no-slip results is presented in the right-hand side plot of Figure 7. The plot confirms the similarity of the results, which is likely attributed to the side configuration of the structure. Since the majority of the fluid travels through the center of the unit cell in the yz plane, the walls have minimal contact with the fluid, rendering the influence of the velocity boundary condition on them negligible. The larger deviations observed for the S2 variation can be attributed to the onset of vortex shedding, resulting from the larger hydraulic diameter and higher particle Reynolds number (up to 105 for S2, compared to maximum values of 80 and 92 for S0.5 and S1, respectively). Vortex shedding is likely suppressed by the no-slip boundary condition, as indicated by the corresponding values below the x = y curve in the parity plot of Figure 7. Based on these findings, it can be concluded that the presence of column walls at the front and back of the structure in the experimental setup did not significantly impact the results.
Continuing the analysis on the influence of geometrical attributes on performance, the impact of fixing or scaling the width is investigated. Figure 8 presents the results for the scaled geometries across all stacking configurations.
The most notable observation concerning the scaled geometries is the relatively high transverse dispersion exhibited by the S2 geometries at high Péclet numbers, which was not observed in previous results. These high values coincide with relatively high COV values, indicating the presence of unsteady vortices that enhance transverse dispersion. Vortex shedding is facilitated by unobstructed axial channels, and it is shown in the top view in Figure 1 that the S2scaled geometries contain these, while they are not present in any of the other geometries. This is confirmed by the visualization in Figure 9, where it is seen that the tracer plume oscillates over time for the S2 geometry, while it is stable for the S0.5 geometry. Such oscillations were not observed in previous fixed geometries, likely since the lower hydraulic diameter suppressed vortex shedding [31]. To show the relation between enhanced transverse dispersion and the unsteady character of the flow, the relative transverse dispersion coefficient is plotted as a function of the COV associated with time-averaging in Figure 10 for the S2 structures. It is clear that there is a strong correlation between these parameters, indicating that unsteady vortex shedding contributes strongly to the observed transverse dispersion.
The enhanced transverse dispersion is most pronounced for the 2D and the stag stacking configurations, as seen in Figure 8. In addition to this, it is observed that both the 2D and the stag stacking configurations exhibit higher transverse dispersion for S1 compared to S0.5, at higher Péclet number, unlike the side and center variations. In the 2D configuration, reduced obstructions in the yz plane facilitate the earlier onset of the formation of vortices. Moreover, the high porosity necessitates relatively high superficial velocities to maintain Péclet numbers equal to the operating conditions in the 3D structures. This point is demonstrated by plotting the relative transverse dispersion coefficient for the various S2 structures as a function of the particle Reynolds number rather than the Péclet number, as in Figure 10. In this graph, the increase in relative transverse dispersion coefficient thanks to the onset of vortex shedding occurs at more similar points on the x-axis compared to Figure 8. The onset of vortex shedding for these structures seems to be at particle Reynolds numbers between 60 and 80, which is in line with expectations from the literature [32,33].
In both graphs of Figure 10, the stag configuration significantly enhances transverse dispersion. In addition to the vortex shedding, the tortuous path in the yz plane likely contributes to the observed high dispersion values.
The results for the side and center stacking variations are very similar, which is demonstrated by the parity plot in Figure 11. Both variations represent aligned stacking in the yz plane, with the center variation shifted by half a unit cell length compared to the side variation. Given the slip velocity boundary condition on the walls, similar hydrodynamics are expected in both geometries, which is confirmed by the comparable transverse dispersion behavior. However, some outliers are present for the S2 spacing, where the center geometry exhibits a lower transverse dispersion coefficient at higher Péclet number. This discrepancy can again be attributed to the available space for vortex shedding, which is twice as large in the side geometry compared to the center geometry. Due to the lower hydraulic diameter, the development of vortices will occur at higher Péclet number, thus explaining the position of the outliers (at P e m = 118 and P e m = 135 ) below the ideal x = y curve.
To gain a proper understanding of the significance of these results, Figure 12 presents the pressure drop of the different scaled geometries. It is evident that the pressure drop is lower for the 2D geometries compared to the 3D geometries, while the stag geometries exhibit slightly higher pressure drop compared to the other stacking configurations. These relative trends align with expectations, as the 2D variation offers the least amount of resistance to the fluid flow, while the tortuous path of the stag configuration presents the most resistance. Interestingly, at higher velocities, the pressure drop does not increase as significantly as the transverse dispersion.
Figure 12 also includes correlations based on Equation (11), utilizing the Reynolds number as defined in Equation (12), with the spacing as characteristic length and the superficial velocity as characteristic velocity. The fitted coefficients and associated MAPE values are provided in Table 4. Vanapalli et al. [34] demonstrated the applicability of this correlation for describing the pressure drop in 2D chromatographic pillar arrays, which bear close resemblance to the structures studied here. As expected, the fitted values of β are in the same range as those found by these authors. Although the obtained MAPE values are relatively high, they are considered acceptable for this relatively simplified representation of the trends. Attempts to improve the MAPE values by utilizing different characteristic lengths or velocities did not yield better results.
Δ p h = ρ f 2 U 2 d characteristic α R e characteristic β
R e characteristic = ρ f U characteristic d characteristic μ f
To further evaluate the practical viability of these structures, their performance was compared to that of a packed bed of monodisperse spheres. This comparison was achieved by plotting the relative transverse dispersion coefficient against the pressure drop required to operate under these conditions. For the packed bed, the relative transverse dispersion coefficient was calculated using Equation (1), while the pressure drop was determined using the Ergun equation (Equation (13)). The results are presented in Figure 13. All structures currently under study are able to provide enhanced transverse dispersion at constant pressure drop compared to the packed bed. However, it is important to note that this analysis does not consider the changes in catalyst holdup, as it assumes a comparison based on a constant reactor length. Although the exclusion of catalyst holdup changes may present a limitation, it is crucial to recognize that the selection of a comparison approach is inherently subjective. Regardless of the exact basis of comparison, it is expected that especially the S2 geometries at high Péclet number will present a favorable operating trade-off.
Δ p h = 150 μ f d p 2 1 ϵ 2 ϵ 3 U + 1.75 ρ f d p 1 ϵ ϵ 3 U 2

4. Conclusions

In this work, the transverse dispersion characteristics of 3D-printed logpile structures were assessed, exploring novel variations in configuration along the secondary transverse axis. This was performed to address the limited flexibility in structural design in earlier experimental work. It can be concluded that the pseudo-2D nature of previous experimental efforts has underexposed the potential of enhanced transverse dispersion in 3D-printed logpile structures. Specifically, by using structures with a fixed depth rather than scaling this parameter as a function of the structures’ spacing, vortex shedding was suppressed. These vortices enhance mixing and in this way increase the apparent transverse dispersion for S2 structures at high Péclet number. Employing a staggered stacking configuration on both transverse axes further enhances this thanks to the increased tortuosity. It was observed that operation in this regime could be achieved without an equivalent pressure drop penalty. This makes these structures a more viable alternative to conventional randomly packed bed reactors than was previously assumed. However, it remains open for discussion whether these benefits justify the use of 3D printing technology or if this technology is better exploited by considering structures with internal variations, further enhancing the transverse dispersion. A techno-economic evaluation, along with further studies on such advanced structures, is needed to determine whether the improved performance compensates for the likely higher costs compared to conventional catalyst structuring methods.

Author Contributions

Conceptualization, L.R.S.R., I.R. and M.v.S.A.; methodology, L.R.S.R., M.A.A.v.A. and E.v.L.; software, L.R.S.R., M.A.A.v.A. and E.v.L.; validation, L.R.S.R., M.A.A.v.A. and E.v.L.; investigation, L.R.S.R. and M.A.A.v.A.; resources, L.R.S.R. and I.R.; writing—original draft preparation, L.R.S.R.; writing—review and editing, M.A.A.v.A., E.v.L., I.R. and M.v.S.A.; visualization, L.R.S.R.; supervision, I.R. and M.v.S.A.; funding acquisition, M.v.S.A. All authors have read and agreed to the published version of the manuscript.

Funding

This article is based on research undertaken in relation to a project (ZEOCAT-3D) which has received funding from the European Union’s Horizon 2020 research and innovation program, under Grant Agreement No. 814548. The authors would like to thank the EU Horizon 2020 program for this opportunity. This publication only reflects the authors’ views, and neither the funding agency nor the European Commission are responsible for any use that may be made of the information contained therein. This work made use of the Dutch national e-infrastructure with the support of the SURF Cooperative using grant no. EINF-3116 and grant no. EINF-4785.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors on request.

Conflicts of Interest

The authors declare no conflicts of interest. The funders had no role in the design of the study; in the collection, analyses or interpretation of data; in the writing of the manuscript or in the decision to publish the results.

Nomenclature

SymbolRepresentationUnit
CConcentrationmol m−3
C d Drag coefficient-
DDispersion coefficientm2 s−1
dDiameter m
hAxial length m
kTurbulent kinetic energym2 s−2
LHalf column width m
NNumber of data points-
nCounting variable-
P e Péclet number-
pPressure Pa
R e Reynolds number-
tTime s
t * Dimensionless time-
USuperficial velocitym s−1
u Velocity vectorm s−1
uInterstitial velocitym s−1
xPrimary transverse coordinate m
YMass fraction-
ySecondary transverse coordinate m
zAxial coordinate m
α Fitted coefficient-
β Fitted coefficient-
Δ p Pressure drop Pa
Δ x Half injection width m
ϵ Porosity-
μ Dynamic viscosity Pa s
ρ Densitykg m−3
σ Standard deviation m
τ Stress tensor Pa
τ Tortuosity-
ω Specific dissipation rates−1
SubscriptRepresentation
0Reference
fFluid
injInjection
mMolecular
pParticle
tTransverse

Appendix A. Solver Validation: Drag Coefficient

The solver was validated based on a comparison between the calculated drag coefficient and a correlation from the literature for flow past a single cylinder, following the same protocol as in a previous work from our group [35]. While the protocol and solver remained unchanged, there were differences in the velocity range and background mesh refinement, necessitating a repetition of the validation process. The drag coefficient was selected as the parameter for validation since it represents the flow resistance exerted by an object against the fluid flow, which is relevant to the context of this work. Based on empirical data, the relationship between drag coefficient and particle Reynolds number for flow past a single cylinder is described by Equation (A1), proposed by Cheng [36]. It is valid in the range of 0.2 < R e p < 200,000.
For this validation, a case was simulated involving nitrogen flow past a single cylinder with a diameter of 1.5   m m . The cylinder was positioned at the center of a 2D domain of 1.5   c m by 1.5   c m , and a background mesh refinement of 150 μ m was applied. To extract the drag coefficient, the built-in OpenFOAM forceCoeffs functionObject was used and the values were averaged to account for unsteady flow at higher Reynolds numbers. Thirteen different superficial velocities ranging from 0.05 m s−1 to 1 m s−1 were considered, and the results are shown in Figure A1. The agreement between the calculated values and the correlation is good, with a MAPE value of 1.7%. Despite being slightly higher than the value in previous work, this is considered satisfactory for validation, considering the broader velocity range, coarser mesh and the inherent variability present in available empirical data [35,36,37,38].
C d = 11 R e p 0.75 + 0.9 1 exp 1000 R e p + 1.2 1 exp R e p 4500 0.7
Figure A1. The drag coefficient as a function of the particle Reynolds number for flow past a single cylinder. The dashed line represents the empirical correlation (Equation (A1)), and markers represent calculated values.
Figure A1. The drag coefficient as a function of the particle Reynolds number for flow past a single cylinder. The dashed line represents the empirical correlation (Equation (A1)), and markers represent calculated values.
Processes 12 02151 g0a1

Appendix B. Validity of Laminar Flow Regime

Given the moderate particle Reynolds numbers explored in this work (up to 142), it is required to ensure the validity of a laminar model to describe the hydrodynamics. To this end, the relative transverse dispersion coefficient was determined using a k- ω SST RANS turbulence model for the S22D geometry at the highest Péclet number. Based on earlier findings, it is anticipated that any turbulent effects would be most pronounced in the 2D geometries, which are tested at higher particle Reynolds numbers due to their higher porosities. Therefore, this evaluation serves to obtain the maximum deviation that the omission of turbulence would introduce. The turbulent kinetic energy, k, was initialized with an intensity of 5% at the inlet, using the turbulentIntensityKineticEnergyInlet boundary condition. Similarly, the specific dissipation rate, ω , was initialized using the turbulentMixingLengthDissipationRateInlet boundary condition with a mixing length of 2.4   m m (equivalent to the height of two stacked cylinders). The accompanying wall functions were applied for both parameters. To assess the sensitivity of the solution to these input parameters, both the inlet intensity and mixing length were varied. However, it was observed that the influence of these variations was relatively small. For instance, using an initial intensity of 1% resulted in a 0.7% change in the transverse dispersion coefficient, while using a mixing length of 1.5   m m (equivalent to the cylinder diameter) led to a change of 0.05%.
By following the simulation procedure and conditions outlined in the Section 2, a relative transverse dispersion coefficient 0.5% lower than the laminar value was obtained. In the laminar case setup, a COV of 0.29 was calculated. Therefore, it can be concluded that the RANS model yields a slightly different value compared to the laminar model, but this discrepancy falls well within the standard deviation of the laminar simulations. Despite the similar result, conducting a RANS simulation required 35% more computational time compared to a laminar simulation. With this argument and the similar results for the conditions with the highest anticipated turbulent effects, it is concluded that the dispersion is dominated by (laminar) vortex shedding, rendering turbulence modeling unnecessary.

Appendix C. Mesh Refinement Visualization

The mesh refinement near overlapping cylinders is visualized for reference in Figure A2.
Figure A2. Visualization of mesh refinement near overlapping cylinders.
Figure A2. Visualization of mesh refinement near overlapping cylinders.
Processes 12 02151 g0a2

References

  1. Lerou, J.J.; Froment, G.F. Velocity, temperature and conversion profiles in fixed bed catalytic reactors. Chem. Eng. Sci. 1977, 32, 853–861. [Google Scholar] [CrossRef]
  2. Kapteijn, F.; Moulijn, J.A. Structured catalysts and reactors—Perspectives for demanding applications. Catal. Today 2020, 383, 5–14. [Google Scholar] [CrossRef]
  3. Griesinger, A.; Spindler, K.; Hahne, E. Measurements and theoretical modelling of the effective thermal conductivity of zeolites. Int. J. Heat Mass Transf. 1999, 42, 4363–4374. [Google Scholar] [CrossRef]
  4. Miscke, R.A.; Smith, J.M. Thermal conductivity of alumina catalyst pellets. Ind. Eng. Chem. Fundam. 1962, 1, 288–292. [Google Scholar] [CrossRef]
  5. Jayachandran, S.; Reddy, K.S. Estimation of effective thermal conductivity of packed beds incorporating effects of primary and secondary parameters. Therm. Sci. Eng. Prog. 2019, 11, 392–408. [Google Scholar] [CrossRef]
  6. Dixon, A.G. Fixed bed catalytic reactor modelling—The radial heat transfer problem. Can. J. Chem. Eng. 2012, 90, 507–527. [Google Scholar] [CrossRef]
  7. Delgado, J.M.P.Q. A critical review of dispersion in packed beds. Heat Mass Transf. 2006, 42, 279–310. [Google Scholar] [CrossRef]
  8. Shinnar, R.; Doyle, F.J.; Budman, H.M.; Morari, M. Design considerations for tubular reactors with highly exothermic reactions. AIChE J. 1992, 38, 1729–1743. [Google Scholar] [CrossRef]
  9. Stankiewicz, A. Advances in modelling and design of multitubular fixed-bed reactors. Chem. Eng. Technol. 1989, 12, 113–130. [Google Scholar] [CrossRef]
  10. Wehinger, G.D. Improving the radial heat transport and heat distribution in catalytic gas-solid reactors. Chem. Eng. Process. Process Intensif. 2022, 177, 108996. [Google Scholar] [CrossRef]
  11. Rosseau, L.R.S.; Middelkoop, V.; Willemsen, J.A.M.; Roghair, I.; van Sint Annaland, M. Review on additive manufacturing of catalysts and sorbents and the potential for process intensification. Front. Chem. Eng. 2022, 4, 834547. [Google Scholar] [CrossRef]
  12. Lawson, S.; Adebayo, B.; Robinson, C.; Al-Naddaf, Q.; Rownaghi, A.; Rezaei, F. The effects of cell density and intrinsic porosity on structural properties and adsorption kinetics in 3D-printed zeolite monoliths. Chem. Eng. Sci. 2020, 218, 115564. [Google Scholar] [CrossRef]
  13. Huang, K.; Elsayed, H.; Franchin, G.; Colombo, P. Additive manufacturing of SiOC scaffolds with tunable structure-performance relationship. J. Eur. Ceram. Soc. 2021, 41, 7552–7559. [Google Scholar] [CrossRef]
  14. Lawson, S.; Li, X.; Thakkar, H.; Rownaghi, A.A.; Rezaei, F. Recent Advances in 3D Printing of Structured Materials for Adsorption and Catalysis Applications. Chem. Rev. 2021, 121, 6246–6291. [Google Scholar] [CrossRef] [PubMed]
  15. Lefevere, J.; Mullens, S.; Meynen, V. The impact of formulation and 3D-printing on the catalytic properties of ZSM-5 zeolite. Chem. Eng. J. 2018, 349, 260–268. [Google Scholar] [CrossRef]
  16. Beck, V.A.; Ivanovskaya, A.N.; Chandrasekaran, S.; Forien, J.B.; Baker, S.E.; Duoss, E.B.; Worsley, M.A. Inertially enhanced mass transport using 3D-printed porous flow-through electrodes with periodic lattice structures. Proc. Natl. Acad. Sci. USA 2021, 118, e2025562118. [Google Scholar] [CrossRef] [PubMed]
  17. Buyruk, E. Numerical study of heat transfer characteristics on tandem cylinders, inline and staggered tube banks in cross-flow of air. Int. Commun. Heat Mass Transf. 2002, 29, 355–366. [Google Scholar] [CrossRef]
  18. Hishikar, P.; Gaba, V.K.; Dhiman, S.K.; Tiwari, A.K. Heat transfer analysis of nine cylinders arranged inline and staggered at subcritical Reynolds number. Numer. Heat Transf. A 2024, in press. [Google Scholar] [CrossRef]
  19. Rosseau, L.R.S.; Schinkel, M.A.M.R.; Roghair, I.; van Sint Annaland, M. Experimental Quantification of Gas Dispersion in 3D-Printed Logpile Structures Using a Noninvasive Infrared Transmission Technique. ACS Eng. Au 2022, 2, 236–247. [Google Scholar] [CrossRef]
  20. Greenshields, C.; Weller, H. Notes on Computational Fluid Dynamics: General Principles; CFD Direct Ltd: Reading, UK, 2022. [Google Scholar]
  21. Crank, J. The Mathematics of Diffusion; Oxford University Press: Oxford, UK, 1975. [Google Scholar]
  22. Li, T.; Deen, N.G.; Kuipers, J.A. Numerical investigation of hydrodynamics and mass transfer for in-line fiber arrays in laminar cross-flow at low Reynolds numbers. Chem. Eng. Sci. 2005, 60, 1837–1847. [Google Scholar] [CrossRef]
  23. Gnielinski, V. Heat transfer in cross-flow around single rows of tubes and through tube bundles. In VDI Heat Atlas; Springer: Berlin/Heidelberg, Germany, 2010; pp. 725–730. [Google Scholar]
  24. Delgado, J.M.P.Q. Longitudinal and transverse dispersion in porous media. Chem. Eng. Res. Des. 2007, 85, 1245–1252. [Google Scholar] [CrossRef]
  25. Ahn, D.; Kweon, J.H.; Kwon, S.; Song, J.; Lee, S. Representation of surface roughness in fused deposition modeling. J. Mater. Process. Technol. 2009, 209, 5593–5600. [Google Scholar] [CrossRef]
  26. Puckert, D.K.; Rist, U. Experiments on critical Reynolds number and global instability in roughness-induced laminar–turbulent transition. J. Fluid Mech. 2018, 844, 878–904. [Google Scholar] [CrossRef]
  27. Kandlikar, S.G.; Joshi, S.; Tian, S. Effect of Surface Roughness on Heat Transfer and Fluid Flow Characteristics at Low Reynolds Numbers in Small Diameter Tubes. Heat Transf. Eng. 2010, 24, 4–16. [Google Scholar] [CrossRef]
  28. Brackbill, T.P.; Kandlikar, S.G. Effect of Sawtooth Roughness on Pressure Drop and Turbulent Transition in Microchannels. Heat Transf. Eng. 2007, 28, 662–669. [Google Scholar] [CrossRef]
  29. Buj-Corral, I.; Domínguez-Fernández, A.; Gómez-Gejo, A. Effect of Printing Parameters on Dimensional Error and Surface Roughness Obtained in Direct Ink Writing (DIW) Processes. Materials 2020, 13, 2157. [Google Scholar] [CrossRef]
  30. Hossain, S.S.; Lu, K. Recent progress of alumina ceramics by direct ink writing: Ink design, printing and post-processing. Ceram. Int. 2023, 49, 10199–10212. [Google Scholar] [CrossRef]
  31. Tanino, Y.; Nepf, H.M. Lateral dispersion in random cylinder arrays at high Reynolds number. J. Fluid Mech. 2008, 600, 339–371. [Google Scholar] [CrossRef]
  32. Franke, R.; Rodi, W.; Schönung, B. Numerical calculation of laminar vortex-shedding flow past cylinders. J. Wind Eng. Ind. Aerod. 1990, 35, 237–257. [Google Scholar] [CrossRef]
  33. Khalifa, Z.; Pocher, L.; Tilton, N. Regimes of flow through cylinder arrays subject to steady pressure gradients. Int. J. Heat Mass Transf. 2020, 159, 120072. [Google Scholar] [CrossRef]
  34. Vanapalli, S.; ter Brake, H.J.M.; Jansen, H.V.; Burger, J.F.; Holland, H.J.; Veenstra, T.T.; Elwenspoek, M.C. Pressure drop of laminar gas flows in a microchannel containing various pillar matrices. J. Micromech. Microeng. 2007, 17, 1381–1386. [Google Scholar] [CrossRef]
  35. Rosseau, L.R.S.; Jansen, J.T.A.; Roghair, I.; van Sint Annaland, M. Favorable trade-off between heat transfer and pressure drop in 3D printed baffled logpile catalyst structures. Chem. Eng. Res. Des. 2023, 196, 214–234. [Google Scholar] [CrossRef]
  36. Cheng, N.S. Calculation of Drag Coefficient for Arrays of Emergent Circular Cylinders with Pseudofluid Model. J. Hydraul. Eng. 2013, 139, 602–611. [Google Scholar] [CrossRef]
  37. Tritton, D.J. Experiments on the flow past a circular cylinder at low Reynolds numbers. J. Fluid Mech. 1959, 6, 547–567. [Google Scholar] [CrossRef]
  38. Qu, L.; Norberg, C.; Davidson, L.; Peng, S.H.; Wang, F. Quantitative numerical analysis of flow past a circular cylinder at Reynolds number between 50 and 200. J. Fluids Struct. 2013, 39, 347–370. [Google Scholar] [CrossRef]
Figure 1. A visualization of the geometries under study. The figure presents the front view (xz plane) of an S1center structure, indicating the relevant dimensions and illustrating walls with shaded blocks. The intended flow direction is from the bottom to the top, in the positive z-direction. Also shown is the side view (yz plane), showcasing the four stacking variations in S1 geometries, and the top view (xy plane), displaying the three spacing variations in stagscaled geometries.
Figure 1. A visualization of the geometries under study. The figure presents the front view (xz plane) of an S1center structure, indicating the relevant dimensions and illustrating walls with shaded blocks. The intended flow direction is from the bottom to the top, in the positive z-direction. Also shown is the side view (yz plane), showcasing the four stacking variations in S1 geometries, and the top view (xy plane), displaying the three spacing variations in stagscaled geometries.
Processes 12 02151 g001
Figure 2. A grid size study for the transverse dispersion coefficient of an S12D structure as a function of the blockMesh cell size. The error in the calculated transverse dispersion coefficient and the computational cost are normalized relative to the case with a 37.5   μ m background mesh. The anomalous cost value for the coarsest grid with P e m = 8 is attributed to the unfavorable 32-core decomposition for such a low amount of cells. Simulations were conducted using the optimal settings outlined in Section 2.4 and Section 2.5.
Figure 2. A grid size study for the transverse dispersion coefficient of an S12D structure as a function of the blockMesh cell size. The error in the calculated transverse dispersion coefficient and the computational cost are normalized relative to the case with a 37.5   μ m background mesh. The anomalous cost value for the coarsest grid with P e m = 8 is attributed to the unfavorable 32-core decomposition for such a low amount of cells. Simulations were conducted using the optimal settings outlined in Section 2.4 and Section 2.5.
Processes 12 02151 g002
Figure 3. The relative transverse dispersion coefficient as a function of the dimensionless time for the S22D structure with P e m = 169 . The transverse dispersion coefficient is calculated based on all thirteen layers of the structure and is plotted with a resolution of 0.1 t * .
Figure 3. The relative transverse dispersion coefficient as a function of the dimensionless time for the S22D structure with P e m = 169 . The transverse dispersion coefficient is calculated based on all thirteen layers of the structure and is plotted with a resolution of 0.1 t * .
Processes 12 02151 g003
Figure 4. MAE in axial velocity profile between subsequent uneven layers along the axial coordinate for six exemplary 2D geometries.
Figure 4. MAE in axial velocity profile between subsequent uneven layers along the axial coordinate for six exemplary 2D geometries.
Processes 12 02151 g004
Figure 5. The relative transverse dispersion coefficient as a function of the axial position in a S12D geometry with P e m = 12 , calculated with different reference positions and varying injection width. The dashed line represents the average relative transverse dispersion coefficient obtained for a relative injection width of 8.
Figure 5. The relative transverse dispersion coefficient as a function of the axial position in a S12D geometry with P e m = 12 , calculated with different reference positions and varying injection width. The dashed line represents the average relative transverse dispersion coefficient obtained for a relative injection width of 8.
Processes 12 02151 g005
Figure 6. The relative transverse dispersion coefficient as a function of the Péclet number and spacing for sidefixed geometries with a no-slip velocity boundary condition applied to the walls. Dashed lines represent the correlation obtained from previous experimental work. Only data points within the applicable range of the correlation are shown. The y-axis is scaled logarithmically.
Figure 6. The relative transverse dispersion coefficient as a function of the Péclet number and spacing for sidefixed geometries with a no-slip velocity boundary condition applied to the walls. Dashed lines represent the correlation obtained from previous experimental work. Only data points within the applicable range of the correlation are shown. The y-axis is scaled logarithmically.
Processes 12 02151 g006
Figure 7. The relative transverse dispersion coefficient as a function of the Péclet number and spacing for side-fixed geometries with either a slip (left) or no-slip (center) velocity boundary condition applied to the walls. A parity plot of the results is also shown (right), in which the ideal x = y curve is represented by the black dashed line. The y-axes and the x-axis of the right-hand side plot are scaled logarithmically, and the legend is applicable to all plots.
Figure 7. The relative transverse dispersion coefficient as a function of the Péclet number and spacing for side-fixed geometries with either a slip (left) or no-slip (center) velocity boundary condition applied to the walls. A parity plot of the results is also shown (right), in which the ideal x = y curve is represented by the black dashed line. The y-axes and the x-axis of the right-hand side plot are scaled logarithmically, and the legend is applicable to all plots.
Processes 12 02151 g007
Figure 8. The relative transverse dispersion coefficient as a function of the Péclet number and spacing for all stacking configurations of scaled geometries with a slip velocity boundary condition applied to the walls. The y-axes are scaled logarithmically, and the legend is applicable to all plots. Error bars represent the standard deviation to illustrate the variability of the data due to temporal fluctuations.
Figure 8. The relative transverse dispersion coefficient as a function of the Péclet number and spacing for all stacking configurations of scaled geometries with a slip velocity boundary condition applied to the walls. The y-axes are scaled logarithmically, and the legend is applicable to all plots. Error bars represent the standard deviation to illustrate the variability of the data due to temporal fluctuations.
Processes 12 02151 g008
Figure 9. A visualization of the tracer molar fraction in the penultimate axial layer of S0.5 (left) and S2 (right) stag geometries as a function of dimensionless time. Paraview line integral convolution (LIC) streamlines are included to visualize the local flow direction. The snapshots cover half the width of the domain, centered on the x- and y-axes.
Figure 9. A visualization of the tracer molar fraction in the penultimate axial layer of S0.5 (left) and S2 (right) stag geometries as a function of dimensionless time. Paraview line integral convolution (LIC) streamlines are included to visualize the local flow direction. The snapshots cover half the width of the domain, centered on the x- and y-axes.
Processes 12 02151 g009
Figure 10. The relative transverse dispersion coefficient as a function of the particle Reynolds number (left) and the COV (right) for all stacking configurations of S2-scaled geometries with a slip velocity boundary condition applied to the walls. The y-axes and right-hand side x-axis are scaled logarithmically, and the legend is applicable to both plots. Error bars in the left-hand side plot represent the standard deviation to illustrate the variability of the data due to temporal fluctuations.
Figure 10. The relative transverse dispersion coefficient as a function of the particle Reynolds number (left) and the COV (right) for all stacking configurations of S2-scaled geometries with a slip velocity boundary condition applied to the walls. The y-axes and right-hand side x-axis are scaled logarithmically, and the legend is applicable to both plots. Error bars in the left-hand side plot represent the standard deviation to illustrate the variability of the data due to temporal fluctuations.
Processes 12 02151 g010
Figure 11. A parity plot of the side and center results for the scaled geometries as a function of the spacing. Both axes are scaled logarithmically. The ideal x = y curve is represented by the black dashed line.
Figure 11. A parity plot of the side and center results for the scaled geometries as a function of the spacing. Both axes are scaled logarithmically. The ideal x = y curve is represented by the black dashed line.
Processes 12 02151 g011
Figure 12. Pressure drop as a function of the superficial velocity and spacing for all stacking configurations of scaled geometries with a slip velocity boundary condition applied to the walls. The y-axes are scaled logarithmically, and the legend is applicable to all plots. Fitted correlations based on Equation (11) with the parameters in Table 4 are shown as dashed lines.
Figure 12. Pressure drop as a function of the superficial velocity and spacing for all stacking configurations of scaled geometries with a slip velocity boundary condition applied to the walls. The y-axes are scaled logarithmically, and the legend is applicable to all plots. Fitted correlations based on Equation (11) with the parameters in Table 4 are shown as dashed lines.
Processes 12 02151 g012
Figure 13. The relative transverse dispersion coefficient as a function of the pressure drop and spacing for all stacking configurations of scaled geometries with a slip velocity boundary condition applied to the walls. The trade-off for a packed bed of monodispersed 1.5 mm spheres and 37% porosity is also shown as a black line. Both axes are scaled logarithmically, the legend is applicable to all plots and error bars are omitted for clarity.
Figure 13. The relative transverse dispersion coefficient as a function of the pressure drop and spacing for all stacking configurations of scaled geometries with a slip velocity boundary condition applied to the walls. The trade-off for a packed bed of monodispersed 1.5 mm spheres and 37% porosity is also shown as a black line. Both axes are scaled logarithmically, the legend is applicable to all plots and error bars are omitted for clarity.
Processes 12 02151 g013
Table 1. The porosities of the various geometries based on stacking configuration and spacing. The side, center and stag variations have identical porosities and are collectively denoted by the ‘3D’ entry.
Table 1. The porosities of the various geometries based on stacking configuration and spacing. The side, center and stag variations have identical porosities and are collectively denoted by the ‘3D’ entry.
S0.5S1S2
2D67.5%75.6%83.7%
3D—fixed47.1%54.7%62.2%
3D—scaled40.3%54.7%69.4%
Table 2. The boundary conditions for Equation (6) for the injection of a tracer with a fixed injection width of 2 Δ x and a fixed concentration of C 0 in a domain with a fixed width of 2 L .
Table 2. The boundary conditions for Equation (6) for the injection of a tracer with a fixed injection width of 2 Δ x and a fixed concentration of C 0 in a domain with a fixed width of 2 L .
z = 0 Δ x x Δ x C = C 0
z = 0 L x < Δ x Δ x < x L C = 0
z > 0 x = L x = 0 x = L C x = 0
z L < x < L C = Δ x L C 0
Table 3. The COV of time-averaging and pressure drop per layer for the six exemplary 2D geometries from Figure 4.
Table 3. The COV of time-averaging and pressure drop per layer for the six exemplary 2D geometries from Figure 4.
StructureCOV [-] Δ p [Pa]
S0.5  P e m = 8 1.52 · 10 5 8.34 · 10 2
S1  P e m = 8 6.69 · 10 6 3.06 · 10 2
S2  P e m = 8 2.11 · 10 7 1.51 · 10 2
S0.5  P e m = 169 2.54 · 10 4 4.87
S1  P e m = 169 7.32 · 10 2 2.62
S2  P e m = 169 7.29 · 10 2 1.28
Table 4. Fitted coefficients based on Equation (11) for the pressure drop per stacking configuration.
Table 4. Fitted coefficients based on Equation (11) for the pressure drop per stacking configuration.
Stacking Configuration α β MAPE
2D11.50.418.9%
side61.90.586.4%
center64.40.606.5%
stag49.80.4710%
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Rosseau, L.R.S.; van Aarle, M.A.A.; van Laer, E.; Roghair, I.; van Sint Annaland, M. Enhanced Transverse Dispersion in 3D-Printed Logpile Structures: A Comparative Analysis of Stacking Configurations. Processes 2024, 12, 2151. https://doi.org/10.3390/pr12102151

AMA Style

Rosseau LRS, van Aarle MAA, van Laer E, Roghair I, van Sint Annaland M. Enhanced Transverse Dispersion in 3D-Printed Logpile Structures: A Comparative Analysis of Stacking Configurations. Processes. 2024; 12(10):2151. https://doi.org/10.3390/pr12102151

Chicago/Turabian Style

Rosseau, Leon R. S., Martijn A. A. van Aarle, Egbert van Laer, Ivo Roghair, and Martin van Sint Annaland. 2024. "Enhanced Transverse Dispersion in 3D-Printed Logpile Structures: A Comparative Analysis of Stacking Configurations" Processes 12, no. 10: 2151. https://doi.org/10.3390/pr12102151

APA Style

Rosseau, L. R. S., van Aarle, M. A. A., van Laer, E., Roghair, I., & van Sint Annaland, M. (2024). Enhanced Transverse Dispersion in 3D-Printed Logpile Structures: A Comparative Analysis of Stacking Configurations. Processes, 12(10), 2151. https://doi.org/10.3390/pr12102151

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop