Hydrodynamics in a randomly packed bed of spheres

Multitubular reactors are commonly used in industry for processes involving highly exothermic chemical reactions. This reactor type consists of individual tubes, with a small diameter compared to the particle size. These slender beds facilitate heat management, but also give rise to flow maldistribution, which decreases the reactor efficiency. The aim of this article is to validate particle-resolved simulations using Magnetic Resonance Imaging experiments while focusing on the flow maldistribution. The packing structure used in the simulations is reconstructed from the experimental images to facilitate a one-to-one comparison. A good match between experiments and simulations is found for the averaged flow profile, probability density function of the velocity in axial direction and even the local velocity distributions. However, a correction of the experimental results for magnetic susceptibility artifacts is necessary to obtain a similar match in the probability density functions and the local profiles.


| INTRODUCTION
Catalytic slender packed bed reactors are commonly used in the chemical, biochemical, and petrochemical industries because of their high wall-to-bed heat transfer coefficients facilitating efficient removal of heat produced by highly exothermic reactions.These reactors consist of a randomly packed bed of solid catalyst particles.Because of the low ratio of the bed diameter to the particle diameter in the slender packed bed reactors, flow maldistribution can easily arise due to the influence of the tube walls on the packing structure.Therefore, understanding fluid flow maldistribution is essential for optimal packed bed reactor design.To study the fluid flow in packed beds, both experiments and simulations have been performed in literature for a wide range of particle Reynolds numbers (Re p ), which are summarized in a table in Appendix A. In the next two subsections, the current state-of-art on this topic for experiments and simulations will be discussed in more detail.

| Experiments
The noninvasive measurement techniques used to gain insight in flow maldistribution in packed bed reactors can be divided Amirhossein Eghbalmanesh and Noah Romijn share first authorship.
into optical techniques and Magnetic Resonance Imaging (MRI).
The optical techniques include laser doppler velocimetry, [1][2][3][4] particle image velocimetry [5][6][7][8] and photoluminescent volumetric imaging. 9,10To obtain optical access required for these techniques in the generally opaque systems, the refractive indices of the solid particles and the fluid should be matched resulting in a transparent system.As optical access is not required in MRI, it is a suitable technique to study hydrodynamics in opaque systems, for example, packed beds. 11,12Using MRI, Sederman et al. 13,14 reported on structure-flow correlations within the inter-particle space of a random packed bed of glass ballotini.It was found that approximately 8% of the pores are responsible for 40% of the total volumetric flow rate indicating the significant heterogeneity in the flow patterns.It clearly shows that the local packing structure of the bed affects the flow in the bed.Johns et al. 15 found that the transition from creeping to inertial flow in isolated pores depends on the local Reynolds number (0:84 < Re p < 14:52).
Furthermore, the features of the flow pattern scale with the overall volumetric flow rate, which implies that there is an enhancement of the inherent heterogeneity in the flow field with an increase of the flow rate.Suekane et al. 16 studied a simplified pore geometry (Re ¼ 12:17 À 204:74) enabling the determination of the main factors influencing flow inhomogeneities, for example, pore size and intrinsic flow instabilities.Ren et al. 17  Based on the flow paths of fluid elements, the fastest motion is observed near the tube wall.Lovreglio et al. 18 showed that the structure and the flow in packed beds of spheres with different sizes can be described.In this article, we will apply MRI measurements to determine the packing and flow inside a packed bed of spheres and will provide a one-to-one comparison with Particle Resolved Computational Fluid Dynamics (PR-CFD).

| Simulations
To compute local flow and pressure fields in a packed bed, PR-CFD can be used.In this type of modeling, the flow field is resolved using a resolution significantly higher than the particle size.The two main PR-CFD methods are: (i) Body-conformal grid methods and (ii) structured grid methods.0][21][22][23][24][25][26] Although such a mesh can represent the interface accurately, the mesh needs to be extremely fine at the contact point of the solid particles, which leads to high computational costs. 27To solve these issues at the contact points, different approaches are proposed, for example, shrinking of the particle [28][29][30][31][32] and using a cylindrical bridge between the particles. 30,33However, these methods change the particle size and thus the porosity of the bed. 34,31ructured grid methods resolve the entire domain, that is, in the solid and fluid phases.A large advantage of these methods is that the grid generation in these methods is simplified compared to the unstructured grid, also in the vicinity of particle contacts. 355][36][37][38][39] The different methods can be classified as continuous forcing methods, which locally explicitly enforce the no-slip boundary by adding an extra force density in the calculation of the Navier-Stokes equations, and direct forcing methods, which directly enforce the no-slip boundaries at the level of the discretised Navier-Stokes equations.By enforcing the no-slip boundary directly in the discretized Navier-Stokes equations, the direct forcing methods create a very sharp definition of the fluid-solid interface, which avoids the use of calibration of the particle size in these methods. 35Therefore, we have chosen to use the direct forcing method of Deen et al., 40 which uses a second-order interpolation method to enforce the no-slip boundary condition on the particle surface.

| Objectives
The aim of this article is to experimentally validate simulations by demonstrating a good correspondence between measurement and simulation.To this end, the MRI results on the fluid flow will be used to validate the PR-CFD results.Using the reconstructed experimental packing of spherical particles in the flow simulations, the flow field inside the packed bed will be studied.Besides the comparison of the averaged axial velocity, also the velocity in the longitudinal and transverse slices will be compared.Moreover, we build upon the previous work of Lovreglio et al., 18 where the authors showed how MRI can be used to obtain data for benchmarking of flow simulations.In this article, we extended the qualitative comparison of Lovreglio et al. 18 to a quantitative comparison of results.To enable an one-to-one comparison between the MRI results and the PR-CFD, an image postprocessing method will be introduced that reduces the effects of the magnetic susceptibility issues leading to the questionable probability density function (PDF) of the velocity as obtained by Lovreglio et al. 18 Finally, the implementation of the image postprocessing method also allows for the comparison of the local velocity profiles over any line in the domain, which has, to the best of the authors' knowledge, not been reported before.
The next section describes the PR-CFD method used for simulating the experimentally obtained packed bed.This article continues with an explanation of the experimental setup, including the packing configuration and the MRI settings.Then, the method for analyzing the MRI images and the postprocessing of the images will be discussed in detail.In the results and discussion section (Section 4), the results of MRI and PR-CFD are compared for the local and average velocity fields.

| NUMERICAL METHOD
In the PR-CFD model, the Navier-Stokes (Equation 1) and continuity (Equation 2) equations are solved for an incompressible, Newtonian fluid with constant physical properties.
To determine the velocity and pressure field, a projectioncorrection method is used on a staggered grid.In this method, the velocity field is calculated by the Navier-Stokes equations using firstorder time discretization.The convective and source terms are treated explicitly.The viscous term is evaluated semi-implicitly where the implicit part is chosen such that each velocity component can be treated separately while the explicit part contains mainly the mixed derivatives of other velocity components.The diffusion term is discretized using a central differencing scheme and the convective term is discretized with a second-order flux-delimited Barton scheme.After the projection step, the velocity field is corrected to ensure mass conservation via the continuity equation.The resulting linear equations for the velocity field calculation and the pressure correction are solved using an in-house modified BiCGSTAB2 solver using the Incomplete-LU decomposition and Algebraic Multigrid (AMG) preconditioners of Trilinos, 41 respectively.
To ensure the no-slip boundary condition at the particle surface and the column wall, the second-order implicit IBM of Deen et al. 40 is applied.In this method, all the fluid cells that have one or more neighbours inside a solid object will be treated differently.After discretization, velocities are required for the neighboring cells inside the solid, which are determined using a second-order extrapolation of the velocity on the liquid-solid interface and the velocities in the liquid phase.It should be noted that if the distance between the central cell of the discretization stencil and the solid is small (less than 0.05 times the grid spacing), the extrapolation will be performed with a firstorder fit to prevent any singularities.Finally, if a liquid cell has two neighboring solid cells in the same direction a different treatment is required which includes both solid boundary conditions.The used method was already validated with the experimental results of Suekane et al. 16 in previous work. 42

| Simulation setup
In the PR-CFD simulations, the packed bed reconstructed from the MRI images (see Section 3) is used that includes both the spherical particles and the cylindrical column wall.To dampen the eddies that form behind the bed an outlet region is added to the bed (0.01 m).
However, this extension is not able to dampen the eddies for Re p ¼ 50.Therefore, the outlet region in this case is modeled with a gradual increase of the viscosity to 100 times the set viscosity.The other parameters used in the simulation are given in Table 1.Based on a previously conducted grid convergence tests, 42   To ensure imaging at the particle Reynolds number (Re p ¼ ρvind μ ) of 10, 30 and 50, the volumetric flow rate was set to 69.3, 207.8 and 346.4 mL min À1 resulting in a superficial velocity, v in , of 0.33, 1.00 and 1.67 cm s À1 , respectively.The viscosity, μ, and density, ρ, for water are obtained from the literature at 20 C. 43

| Image acquisition
Three-dimensional (3D) images of the packed bed in stagnant water were acquired with a Fast-Low-Angle-Shot sequence (FLASH), using the image parameters in Table 2.Because of the limitations of the receiver field, the column was imaged in three different sections: the top, middle and bottom sections.The slices were taken with the image plane perpendicular to the central axis of the flow column.
The velocity image data were acquired by a Phase Contrast Gradient Echo (PCGRE) sequence. 44The velocity encoding was performed using a bipolar gradient pulse along the slice select axis, which leads to a phase shift proportional to the velocity in the main flow direction.
It should be noted that a reference scan with flow compensation is subtracted to eliminate the background phase shift.The image parameters of these experiments are also provided in Table 2.

| Structural images
To obtain a single 3D volume image of the packed bed, a volume reconstruction was performed in Mathematica 12 which combines the acquired images.As it is mentioned in Section 3.2, the column is imaged in three different sections of 256 slices perpendicular to the main direction of the flow (all with a 128 Â 128 Â 256 acquisition matrix), these different sections need to be 'stitched' together.To perform an accurate reconstruction of the entire column, the overlapping regions were determined based on reference points in the volume reconstruction of these sections.These reference points were determined with the "ImageCorrespondingPoints" function.The three sections were then aligned with each other by applying the required image transformation which was obtained with the "FindGeometric-Transform" function.Only translation and rotation were allowed.
Based on the volume reconstruction, the center of mass of each particle was determined using the spherical Hough algorithm in Matlab. 45To correct for any overlap obtained between particles and the particles and the column wall, an in-house swelling algorithm was applied, which grows the particles from the obtained centers of mass at a constant rate until they touch a neighboring particle or the wall.
Upon further swelling, touching particles push each other such that no overlap is created.This approach allows the particle position to be altered slightly during swelling.Based on a synthetic packing, the mean absolute error in the particle positions obtained by the spherical Hough transform was estimated to be 0.13 mm.

| Velocity images
Near liquid-solid interfaces with different magnetic susceptibility artifacts might occur, which causes local field inhomogeneity. 46To remove these artifacts, this article applies boundary postprocessing (BPP) which treats the fluid voxel next to the solid particles slightly differently.When such a voxel features a negative velocity, the velocity inside the voxel was set to zero.Both the results with and without BPP are corrected to ensure mass conservation over the different slices.The correction consists of imposing the same superficial velocity, that is known from the imposed flow rate, in all slices.This is achieved by adding a constant value to all velocities in a slice.Finally, to correct for the geometric distortion present in the image data, an affine transformation was applied.In this section, the velocity profiles of the packed bed of spherical particles, which are acquired from the MRI experiment and PR-CFD simulations will be compared with each other.Firstly, the average velocity and flow profiles and then the local velocity values in the bed will be compared.To determine the fluid flow maldistribution, the velocity data is also represented in Figure 2 by the cumulative flow profiles.Cumulative flow profiles explicitly show the fraction of the total flow that passes through a radial section of the tube.Very good correspondence between PR-CFD and MRI measurements is found.In addition, these graphs clearly show that the fluid distribution is almost independent of the Re p .There are two main flow paths in these cumulative flow profiles one close to the wall for S < 0:5 and an inner region with S between 0.5 and 1, which is different compared to our previous results. 42These differences are caused by the differences in the packings.The packing in this work is consolidated while the packings in Fathiganjehlou et al. 42 are not densified.This results in a more ordered packing in this research leading to two preferential pathways each having 30% for S < 1 of the total flow compared to the pathway which constitutes to 50% of the total flow at the same distance from the wall in Fathiganjehlou et al. 42

| PDF of the velocity
To determine the distribution of velocities over the entire column the PDFs of the normalized velocity are presented in Figure 3.The figure shows the PDFs before and after application of the BPP analysis.Earlier works reported large deviations between PDFs from MRI and PR-CFD. 18,47The figure also clearly shows significant differences in the shape of the PDFs.After application of the BPP, a better match is found between the experimental and numerical PDFs.An excellent match can be seen at Re p ¼ 50.With decreasing Re p , the peak position in the BPP-MRI PDF shifts to lower values and the negative tail becomes larger, which is not observed for the PR-CFD.These shifts in the MRI data are caused by the larger fraction of velocities close to zero at lower Re p , which are not captured well by the MRI measurement.where V is the velocity and N the number of compared pixels in the MRI image.

| Local velocity comparison
The RMSE of the average radial porosity shows that there is a small error in the reconstruction of the packing, which is explained in Section 3.3.1.The table clearly shows an improvement of the RMSE values when the BPP method is applied, which implies the importance of addressing the magnetic susceptibility artifacts in the results.The table shows a smaller effect of the BPP on the cumulative flow profile compared to the average radial velocity profile, which indicates that the effects occur at the center of the bed as the surface area is small and the velocity is low due to the low porosity.This is also confirmed by comparing the differences obtained in the radial flow profiles in Figure 2A, C, E and the cumulative flow profiles in Figure 2B, D, F.
The RMSE values of the PDF shows a significant improvement of the RSME when the BPP method is used.This is also clear when the results presented in Figure 3A,C,E are compared to the results of Lovergilo et al. , 18 as it results in a significant decrease in the negative tail of the pdf is observed.In addition to the improvements in the results for the average velocity profiles, the RMSE of the local velocity is also improved by applying the BPP method.
Finally, the results given in Table 3 show that the RMSE values decrease with increasing Re, which was also noticeable in the comparisons shown in Figures 2-5.The higher difference at low Re p is explained by the decrease in the velocity signal-to-noise ratio for the measurement of low Re p .

| DATA AVAILABILITY AND REPRODUCIBILITY STATEMENT
This section is divided into two parts: MRI experiments and PR-CFD simulations.The purpose of these sections is to define the settings and postprocessing procedures for each method separately.The data used and required scripts to plot

| Velocity profiles
The velocity data obtained with the experimental and postprocessing procedure described in Sections 3.1,3.

| Procedure error analysis
The measurement error for the flow experiments was determined via the superficial velocity.The deviation in the superficial velocity w.r.t. the set superficial velocity could be determined based on the condition that the superficial velocity in the z-direction (along the central axis of the column) needs to remain constant to provide for mass The RMSE of the average radial porosity, the average radial flow profile, the cumulative flow rate, the normalized velocity PDF and the normalized velocity for both the comparison of the PR-CFD to the MRI images without BPP (no BPP) and with BPP (BPP). is the average velocity and ϵ the porosity.The error analysis was carried out with a 95% confidence interval.It should be noted that the error analysis could only be performed on the raw velocity data.

| PR-CFD simulations
The PR-CFD simulations were conducted on a slender packed bed of spherical particles reconstructed from MRI structural images, as explained in Section 3.3.The exact geometrical information of the packed bed is provided in the Data S1 inside the folder /Reconstructed packing.The packed beds were imported into an in-house C++ code, and the simulations were carried out using the settings described in Table 4.It is important to note that the center of the cylindrical column should be adjusted according to the center of the packing, which is determined from the MRI simulation.The validation and grid dependency of the code are extensively discussed in Fathiganjehlou et al. 42 The simulations were run in parallel mode using 64 mpi-threats.
To post-process the results, data extracted from all subdomains were combined using the MATLAB script provided in the Data S1 folder named "Matlab script for PR-CFD."These MATLAB scripts were used to generate the PR-CFD data presented in Figures 2 and 3.The raw data used to generate the graphs is also available in the Data S1 folder with the same name.The results presented in Figures 2 and 3were calculated for the height range between 0.015875 and 0.040875 of the bed (determined according to the slice range used in MRI postprocessing).It is important to note that all velocity data extracted from the simulation were normalized with the inlet velocities specified in Table 4.
Table 5 shows the parameters that are used for Figure 3.
The PR-CFD portions of Figures 4,5 and Figure 2 in Appendix B were calculated using a Python script that can be found in the supplementary materials folder named "Python script for PR-CFD" (the data that are used for Figure 1 of Appendix B is similar to Figure 4).The raw data of the simulation used in the calculation can also be provided upon request.As the resolution of the PR-CFD simulations was higher This rotation may result in porosity profile values exceeding 1 or falling below 0, which were manually set to 1 or 0, respectively.

| Error analysis
The error analysis data which are reported in Table 3 are calculated using the values which are located inside the folder /Error Analysis.

| CONCLUSION AND OUTLOOK
In this article, a comparison between PR-CFD simulations and MRI experiments is reported in order to validate the numerical prediction of single-phase flow in a randomly packed bed with spherical particles.Since the real packing configuration was reproduced, the data shown here did not only allow for a qualitative but also a quantitative comparison of the simulated and experimentally observed velocity distributions.The flow simulations and experiments reported were conducted at Re p = 10, 30, and 50.This data can be used as benchmark for flow simulation studies.
The results show a significant improvement in the experimentally resolved velocity around the particles after the application of the proposed BPP method, which sets any negative velocities in the vicinity of a liquid-solid interface to zero.At Re p ¼ 50 the simulated results are in good agreement with the experimental data.This shows that the PR-CFD model is able to predict the local flow velocity maps in the packed bed.At Re p = 10 and 30, there are larger differences between the PR-CFD results and BPP-MRI results.This can be caused by the larger experimental error in the flow velocity that exists at lower Re p .It should be noted that these differences are more apparent in the local assessment of the flow velocity, while the comparison of average velocity and flow profiles is still very good.In both the simulations and the experiments flow channelling occurred both close to the wall and between the first and second layers of particles from the column wall, which is in contrast with literature observations.This is a consequence of the packing, which was more consolidated in these measurements, indicating the strong effect of the porosity distribution on the flow maldistribution.
Although a significant improvement is found in the agreement when the BPP is used, a one-to-one comparison remains still challenging, because the BPP only corrects for the negative velocities near the solid phase while the magnetic susceptibility artefacts can both manifest as positive and negative velocities.As the velocity field is corrected to match the superficial velocity, the remaining magnetic susceptibility artifacts affect the entire flow field.Therefore, it is recommended to further investigate methods to also remove the positive velocities created by the susceptibility artifacts.For further research of this method it would also be interesting to compare the results with an magnetic susceptibility matched system or with results acquired by a spin echo sequence that is less prone to magnetic susceptibility differences between materials because of its pi refocusing pulse.This will allow for the investigation of the effect of the magnetic susceptibility artifacts.

ACKNOWLEDGMENTS
We would like to thank SURF (https://www.surf.nl)for the support in using the National Supercomputer Snellius.In addition, we thank Claire M.Y.Claassen for providing a python script for some parts of the postprocessing of the simulation results.
applied a spin-tagging MRI to directly visualize the evolution of flow patterns in a packed bed of glass beads.
the simulations are performed with 40 grid cells per particle diameter.To ensure a correct implementation of the method, the pressure drop (ΔP=L) over 0.21 m bed at Re p ¼ 50 (3662 Pa m À1 in the packing region) was compared to experimental pressure drop (3846 Pa m À1 , using a differential pressure transmitter MDM490).Moreover, the cumulative flow rates calculated by the PR-CFD simulations have an error of 0.042% compared to the inlet mass flow rate, which shows the error of the PR-CFD simulations regarding integral mass conservation.In Section 4, the inlet section of 5 particle diameters and the outlet section of 4 particle diameters of the bed are not considered in the comparison with MRI measurements to remove any effects of the inflow and outflow regions.

3 | EXPERIMENT 3 . 1 |
SetupAll experiments were conducted on an MR Solutions Ltd MRI system with a 7 T, 17 cm bore horizontal magnet equipped with shielded gradient coils for x-, y-, and z-direction, supplying each a maximum gradient strength of approximately 0.5 T m À1 .A proton ( 1 H) transmitreceive birdcage radio frequency coil was used with an inner diameter of 35 mm and length of 60 mm, which was tuned to 298.1 MHz in order to image the 1 H signal from liquid water.Image data were recorded with the Preclinical Scan™ software.A hollow polycarbonate (LEXAN ® ) cylindrical column with an internal diameter of 21 mm and length of 61 mm was packed with polypropylene spherical particles with a diameter of 3 AE 0.05 mm.The column contains a glass (SCHOTT DURAN ® ) frit with a pore size of 90 μm at the inlet to provide a uniformly distributed liquid flow during operation.The packing was consolidated by vibration using an orbital shaker (Lab Dancer, Ika-Werke GMBH) with a preset vibration time of 10 min followed by tapping.The free space created at the outlet of the column was filled up with open-cell foam to ensure a static packing.The flow column was inserted into a plastic holder of polyoxymethylene to fix its position in the magnet.The column was filled with deionized water prior to imaging and any air bubbles were expelled from the set-up by tapping on the exterior of the column.A peristaltic pump (Qdos60, Watson-Marlow Ltd) was used to flow water through the bed where buffer-vessels were used to dampen pulsations.The flow rate was monitored by an ultrasonic volume flow meter (Cori-Flow M15-AAD-55-0-S, Bronkhorst ® ) throughout the duration T A B L E 1 The parameters that are used in the particle resolved computational fluid dynamics (PR-CFD) simulations., 0:375, 0:225 ms of the experiment.A schematic representation of the experimental setup is shown in Figure 1.

1
Schematic representation of the experimental setup used in the acquisition of image data.T A B L E 2 The parameters that are used in the magnetic resonance imaging (MRI) experiments., and 25 cm s À1 F I G U R E 2 Radial profiles obtained from MRI measurements (square markers) and PR-CFD simulations (solid lines) showing the porosity (red; A, C, E) and normalized flow velocity (black; A, C, E), and the dimensionless cumulative flow profile (B, D, F) averaged along the central axis of the column as a function of the normalized distance from the column wall for the different Re p of 10 (A, B), 30 (C, D), and 50 (E, F).

4. 1 |
Average radial fluid flow profilesThe average radial superficial velocity profiles of the experiments and the simulations are compared in Figure2A,C,E.Besides the fact that the radial velocity profile follows the trend of the radial profile of porosity, the figure shows a good match between the experiments and the simulations.The most significant differences are obtained in the center of the bed and close to the wall, although the trends in these regions are still similar.These differences are probably caused by the resolution of the MRI images and the errors in the reconstruction of the packing.To remove the inaccuracies due to the resolution of the MRI in the radial profile of porosity close to the column wall the value close to the wall has been set to one manually.

4 . 4 |
Besides the average fluid flow, the local velocities are also compared in this article on a random line at the center of the bed.To compare the local velocities, Figure4shows the two-dimensional normalized velocity maps parallel (longitudinal image) and perpendicular (transverse image) to the main flow direction at different Re p on a chosen line at the center of the bed (the comparison of the local velocity has been performed for more positions in the bed, which can be found in Appendix B).The interested reader is welcome to use the provided detailed full field data as outlined in Section 5 to perform their own reconstruction.In general, the flow patterns in the PR-CFD coincide with the BPP-MRI.The figures also show high velocities near the column wall in the voids created by the structuring of the packing, which is a clear example of flow channelling.The regions with back flow are observed in both the experiments and the simulations.The differences between the experimental results are caused by the errors in the flow imaging and the small deviations in the packing structure.The error analysis shows a measurement error of 57%, 23%, and 9% with Re p ¼ 10, 30 and 50, respectively (see Section 5.1.3for details).The one-dimensional (1D) velocity profiles for different Re p are presented in Figure 5.The figure shows large negative velocities near the solid surfaces in the MRI data before the BPP is applied.These negative values are unexpected considering the no-slip boundary condition.These negative values are obtained due to the magnetic susceptibility variations between the two materials.The BPP-MRI profiles show a better agreement with the PR-CFD simulated profiles.Although the BPP corrects for the negative velocities near the solid-liquid boundary, the negative velocities in the recirculation zones remain present.It should be noted that local inhomogeneities near the interface can also result in positive velocities, as these values fall within the expected velocity range they are difficult to identify.Although the PR-CFD simulated profiles and BPP-MRI results show many similarities, significant differences are obtained at lower Re p .For example, for Re p ¼ 10 the MRI profile shows sharp peaks, while the PR-CFD predicts more smooth peaks.In addition, both positive (Re p ¼ 10, distance ≈ À 0:6) and negative (Re p ¼ 30, distance ≈ 0:75) velocities are sometimes overmeasured by MRI, due to experimental noise.Quantitative differences between MRI and PR-CFD The previous results of the average and local velocity fields show good correspondence between the MRI and PR-CFD results.However, there are mismatches in the local normalized velocity F I G U R E 3 Probability density functions (PDFs) obtained from MRI measurements (square markers) and PR-CFD simulations (solid lines) showing the probability to find a certain normalized velocity within the bed for the different Re p of 10 (A,B), 30 (C,D) and 50 (E,F).Data are shown before (A, C, and E) and after (B, D, and F) the application of the postprocess image analysis to minimize magnetic susceptibility induced artifacts in the measured flow velocity.profileswhich occur mainly close to the solid phase indicating that they are probably caused by the magnetic susceptibility artifacts in the MRI images.To determine if the application of BPP decreases the error, Table3shows the Root Mean Square Error (RMSE) of the average radial flow profile, the cumulative flow rate, the normalized velocity PDF and the normalized velocity for both the comparison of the PR-CFD to the MRI images without BPP and with BPP.The RMSE is calculated based on Equation(3)   (which shows the calculation for the local normalized velocity), Figures 2-5 and the figures in Appendix B can be found in the Data S1: https://doi.org/10.4121/4b8d0613-5ee7-4971-b695-835f287336e8,inside the folder /Figures.5.1 | MRI experiments 5.1.1 | Packing structure The particle positions and radii obtained with the experimental and post-processing procedure described in Sections 3.1,3.2,and 3.3.1 are provided in the Data S1 inside the folder /Reconstructed packing.The raw image data (DICOM format) used to obtain these data is provided inside the subfolder /Raw data.
2, and 3.3.2 is provided in the Data S1 inside the folder /Matlab script for MRI.The script used to calculate the profiles in Figures 2-5 and the figures in Appendix B can be found inside the same folder.The raw spectrometer files (MRD format) used to obtain these data are provided inside the subfolder /Raw data.

F I G U R E 4
Flow fields showing the longitudinal (A,B,E,F,I,J) and transverse (C,D,G,H,K,L) images for both the MRI measured and PR-CFD simulated normalised velocity in the main flow direction for the different Re p of 10 (A-D), 30 (E-H), and 50 (I-L).The flow in the longitudinal images is from left to right.Positive velocities indicate flow in the direction of the main flow.The blue lines are the lines along which the onedimensional velocity profiles shown in Figure 5 are compared.F I G U R E 5 One-dimensional (1D) flow profiles showing the MRI measured (black dashed lines) and PR-CFD simulated (black solid lines) normalized velocity as a function of the normalized distance from the center of the column wall for the different Re p of 10 (A,B), 30 (C,D) and 50 (E,F).The corresponding porosity profiles obtained from the flow experiments (red dashed lines) and reconstruction (red solid lines) are also provided.(Particles are identified by zero velocity values.)See Figure 4 for the lines in the packing where these profiles are measured.Data is shown before (A, C, and E) and after (B, D, and F) the application of the postprocess image analysis to minimize magnetic susceptibility induced artifacts in the measured flow velocity.A positive velocity indicates flow in the main flow direction.than that of the MRI experiments, the graphs presented in these figures are averaged results over four grid cells of the PR-CFD simulations.This averaging applies to the longitudinal, transverse, and all 1D profiles.The simulation domain of the PR-CFD was rotated for 18, and 19 degrees for Re p ¼ 10,30, and 50, respectively, as indicated in the Python script.

T A B L E 4 5
The detailed parameters that are used in the particle resolved computational fluid dynamics (PR-CFD) simulations.Parameters that are used for Figure3.Parameter ValueNumber of bins 140Range of the normalized velocity À20 to 20