# Comparison of Flow and Transport Experiments on 3D Printed Micromodels with Direct Numerical Simulations

Francesca Watson, Julien Maes, Sebastian Geiger, Eric Mackay, Mike Singleton, Thomas Mcgravie, Terry Anouilh, T. Dawn Jobe, Shuo Zhang, Susan Agar, et al.

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Francesca Watson, Julien Maes, Sebastian Geiger, Eric Mackay, Mike Singleton, et al. Comparison of Flow and Transport Experiments on 3D Printed Micromodels with Direct Numerical Simulations. Transport in Porous Media, 2019, 129 (2), pp.449-466. 10.1007/s11242-018-1136-9. hal-02301846

HAL Id: hal-02301846

https://ifp.hal.science/hal-02301846v1

Submitted on 30 Sep 2019

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /></figure>

# Comparison of Flow and Transport Experiments on 3D Printed Micromodels with Direct Numerical Simulations

Francesca Watson<fnref n="1" /> · Julien Maes<fnref n="1" /> · Sebastian Geiger<fnref n="1" /> · Eric Mackay<fnref n="1" /> · Mike Singleton<fnref n="1" /> · Thomas McGravie<fnref n="1" /> · Terry Anouilh<fnref n="2" /> · T. Dawn Jobe<fnref n="3" /> · Shuo Zhang<fnref n="3" /> · Susan Agar<fnref n="3" /> · Sergey Ishutov<fnref n="4" /> · Franciszek Hasiuk<fnref n="4" />

Received: 27 April 2018 / Accepted: 11 August 2018 / Published online: 25 August 2018
© The Author(s) 2018

## Abstract

Understanding pore-scale flow and transport processes is important for understanding flow and transport within rocks on a larger scale. Flow experiments on small-scale micromodels can be used to experimentally investigate pore-scale flow. Current manufacturing methods of micromodels are costly and time consuming. 3D printing is an alternative method for the production of micromodels. We have been able to visualise small-scale, single-phase flow and transport processes within a 3D printed micromodel using a custom-built visualisation cell. Results have been compared with the same experiments run on a micromodel with the same geometry made from polymethyl methacrylate (PMMA, also known as Perspex). Numerical simulations of the experiments indicate that differences in experimental results between the 3D printed micromodel and the Perspex micromodel may be due to variability in print geometry and surface properties between the samples. 3D printing technology looks promising as a micromodel manufacturing method; however, further work is needed to improve the accuracy and quality of 3D printed models in terms of geometry and surface roughness.

**Keywords** 3D printing · Pore-scale flow · Micromodels · Imaging

## 1 Introduction

3D printing, or additive manufacturing, is a growing method for creating custom-built physical models quickly and cost effectively. The potential uses of 3D printing within the geosciences are vast, with applications ranging from on-site manufacturing of replacement

<figure><img src="image_4.png" /></figure>

<footnote>Francesca Watson
f.watson@hw.ac.uk</footnote>

<footnote n="1">¹ Institute of Petroleum Engineering, Heriot-Watt University, Edinburgh EH14 4AS, United Kingdom</footnote>

<footnote n="2">² Institut Français du Pétrole, 232 Avenue Napoléon Bonaparte, 92852 Rueil-Malmaison, France</footnote>

<footnote n="3">³ Aramco Research Center, 16300 Park Row Drive, Houston, TX 77084, USA</footnote>

<footnote n="4">⁴ Department of Geological and Atmospheric Sciences, Iowa State University, 253 Science Hall, 2237 Osborn Drive, Ames, IA 50011, USA</footnote>
parts for drilling rigs, to printing granular samples from real sand grains for geomechanical testing (Ishutov et al. 2018b). This paper is concerned with the use of 3D printing to investigate flow and transport properties of rocks at the pore-scale.

Pore-scale flow and transport phenomena influence macroscopic flow properties within reservoirs (Blunt 2017). One way to investigate pore-scale flow processes is by visualising flow inside micromodels, small-scale models designed to analyse specific fluid flow processes (sometimes referred to as microfluidics). Current micromodel fabrication processes can be costly and time consuming and often require the use of specialist laboratories. Several previous studies have created micromodels by etching pore structures onto silicon wafers either to use directly (Buchgraber et al. 2012) or to use as a mould for making the micromodel out of polydimethylsiloxane (PDMS) (e.g. Gerami et al. 2016; Karadimitriou et al. 2013; Zhou et al. 2010). Porter et al. (2015) created micromodels by etching pore-scale fracture patterns taken from 3D computer microtomography scans of shales into samples of shale. Doing this ensures the physical and chemical properties of the models are realistic, but it cannot overcome the inherent heterogeneity between shale samples.

3D printing lends itself well to the manufacture of micromodels due to the fact that it is easy, quick and relatively inexpensive to repeatably produce 3D printed models with a wide range of geometries. One of the major advantages of 3D printing is that it allows much greater freedom when designing micromodel geometry. In particular, 3D printed micromodels can have variable cross-sectional geometries which is very difficult to achieve using existing methods. 3D printing can be used to create channels with varying thickness as well as shape. This could potentially allow us to control surface roughness by printing surfaces with specified roughness.

The manufacturing process for 3D printed micromodels is often much less complicated than that of other methods. Zhao et al. (2016) had to create several moulds and go through several curing processes before getting a final micromodel to be used in their experiments, whereas the micromodels used in this study only require one printing step and a small amount of post-processing (flushing with ethanol, drying and then curing under ultraviolet light). See Ishutov et al. (2018a) for more details. This simplified manufacturing process is particularly useful for studies where several models with different geometries are required.

Reported resolutions of 3D printers can vary between 0.001 and 100 s of µm (Ishutov et al. 2018b) depending on the type of printer used. This resolution can however also depend on the geometry of the model being printed. For microfluidic studies, the dimensions of models created are typically on the order of 100 µm.

In terms of cost, 3D printers can range in price from 1000 to 100,000s of USD. However several studies, including this one, have created micromodels using lower end printers (~1000s USD which have produced low-cost 3D printed samples (< 1 USD) that have worked well (e.g. Chan et al. 2015; Comina et al. 2014; Shallan et al. 2014). In this study, we used a Formlabs Form 2 printer which cost 4500 USD and each 3D printed sample cost around $0.20 USD (see Sect. 2.2 for more details).

Most work to date that uses 3D printed micromodels to investigate flow and transport processes focusses on biomedical applications. These applications usually involve controlled mixing and reaction of substances. Fluids are initially flowed through separate microchannels that join up within the micromodel, causing them to mix and in some cases react. Kitson et al. (2012) created 3D printed micromodels with various designs for the purpose of carrying out controlled chemical reactions. One of their designs included silos which were filled with reactants whilst the devices were being printed. Subsequently fluid was flowed through the silos and the reactants were dissolved. Comina et al. (2014) created a 3D printed micromodel to assess glucose levels by looking at colour change due to reaction of certain chemicals with
glucose. They also created 3D printed micromodels with geometry specifically designed to enhance mixing in order to efficiently mix certain chemicals in preparation for chemical analyses. Rogers et al. (2015) used a custom resin formulation to 3D print devices with integrated membrane-type valves. The valves were formed by printing a thin enough layer of resin such that certain chemicals could pass through it under enough pressure. The channel widths and depths in these studies ranged from 0.2 to 1 mm.

Previous 3D printed micromodel studies assess mixing and reaction processes by observing the properties of fluids at micromodel outlet (e.g. colour, chemical composition). Instead, in this paper we describe the visualisation of basic flow and transport experiments run on 3D printed micromodels where we can visualise the spatial distribution of mixing processes within the micromodel. This work has implications for visualising mixing processes at the pore-scale.

Several studies have used 3D printing to explore flow processes at the Darcy-scale. Suzuki et al. (2017) investigated single-phase flow and transport in 3D printed, discrete fracture networks. They compared tracer breakthrough curves obtained from numerical simulations with experimental results on 3D printed fracture networks. This work highlighted the importance of using numerical simulations with experimental data to understand the 3D printing process, as discrepancies between the results were attributed to incomplete removal of wax support material used during the printing process and lower than anticipated printer resolution.

Head and Vanorio (2016) used 3D printed samples to investigate the effects of changing rock microstructure on bulk permeability and porosity. They created digital models of pore space geometries from CT scans of rock samples, which were then scaled in different ways to represent the effects of dissolution and compaction. Flow experiments were run on these different models and the relationship between bulk porosity and permeability parameters and properties obtained from the digital models such as tortuosity and matrix surface area were assessed. Their work demonstrates the potential of 3D printing for systematically changing the geometry of printed models to investigate fundamental flow processes and subsequently applying the results to real-world situations.

Ishutov et al. (2017) investigated discrepancies in porosity, permeability and mean pore throat radius between sandstone samples and upscaled 3D printed versions based on CT scans of the natural rock samples. The printed samples had similar permeabilities and scaled mean pore throat radii to the original samples; however, porosity in the printed samples was higher than expected. This was attributed to the nature of the manufacturing process and post-processing of the 3D printed samples. Ishutov et al. (2018b) explored the effect of processing techniques in an attempt to reduce the porosity difference between the 3D printed samples and the natural rocks they were based on.

In contrast to other works which have investigated bulk Darcy-scale properties of 3D printed samples, our work explores the utility of 3D printing for investigating pore-scale flow and transport processes. We have performed single-phase tracer flow experiments on 3D printed micromodels and visualised the results. Image analysis techniques have been used to determine the tracer concentration through time during the experiments. Results have been compared between experiments run on Perspex micromodels and numerical simulations to enable us to identify differences that are due to the micromodel manufacturing method. Our work is important for identifying how suitable 3D printing is as a micromodel manufacturing method for the investigation of pore-scale flow processes.
# 2 Methodology

## 2.1 Conceptual Model

Our conceptual model consists of a 3D pore triplet with constant depth and with a single inlet channel, which splits into three channels of different width, which then recombine into a single outlet channel (Fig. 1). This relatively simple geometry was chosen in the first instance to allow us to investigate how well 3D printing could replicate a simple geometry with relatively simple flow dynamics.

Initially, the system is filled with fluid containing no tracer. Fluid containing a tracer is injected at the inlet at a constant flow rate and spatial variations in tracer concentration are monitored through time. The tracer moves through the system via advection and diffusion.

## 2.2 Micromodel Fabrication

The 3D printed sample (Fig. 2) was printed on a Formlabs Form 2 stereolithography printer. The printer works by using a laser to successively harden layers of resin in the desired shape. This method of printing is preferable to binder deposition printing, which involves powdered material being fused or glued together. This can create a secondary porosity within the walls of the 3D printed sample that is undesirable when flowing fluid through the sample. More details on different methods of 3D printing can be found in Ishutov et al. (2018b). The reported resolution of the printer is 0.025 mm in the x direction, 0.050 mm in the y direction and 0.100 mm in the z direction.

<figure><img src="image_5.png" /><figcaption>**Fig. 1** Original print geometry provided to the printer as an .stl file. Dimensions are in mm. The actual printed channel depth was 0.3 mm deep instead of 0.7 mm due to issues with the printing process. See text for more details</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Fig. 2 a Photograph of the 3D printed micromodel. b Photograph of the Perspex micromodel</figcaption></figure>

The sample geometry sent to the printer consisted of a 2.50 mm thick disc with the pore triplet inscribed into the surface down to a depth of 0.70 mm (Fig. 1). The final printed depth of the channel measured approximately 0.30 mm instead of 0.70 mm (Table 1), which is a greater discrepancy than the reported 0.10 mm printer resolution. This is probably due to the solidification of the channel bottom caused by the orientation of the model in relation to the build platform. This has also been seen in other studies where horizontal printing produces lower resolution channels than vertical printing (e.g. Waheed et al. 2016). We printed the model at 45 degrees to the platform to avoid warping; however, this inadvertently caused the reduction in channel depth. In the future, this could be avoided by tuning the print parameters and designing channels with the optimal channel depth that can be printed accurately.

The Perspex sample can be seen in Fig. 2b. It was cut from Perspex according to the dimensions specified in Fig. 1; however, the channels were cut to a depth of 0.30 mm instead of 0.70 mm to give a closer match to the 3D printed sample.

The dimensions of both samples were measured using scanning electron microscopy (SEM). Table 1 shows the channel depths at the locations shown in Fig. 3. The standard deviation in channel depth from the desired depth of 0.3 mm for the 3D printed sample was 0.11 mm which was greater than the standard deviation for the Perspex sample of 0.08 mm. The absolute error in depth was ± 0.2 mm for the 3D printed sample and ± 0.1 mm for the Perspex sample.

In terms of channel width (Table 2), the 3D printed channel widths had a standard deviation from the desired width of 0.08 mm and an absolute error of ± 0.11 mm whilst the Perspex sample channel widths had a standard deviation from the desired width of 0.05 mm and an absolute error of ± 0.07 mm.

Interestingly, the majority of channel widths for the 3D printed sample were smaller than the intended channel widths whereas the Perspex channel widths were mostly wider than intended. This difference could possibly be a result of the difference between an additive manufacturing process (3D printing) and a subtractive one (cutting away the Perspex).

Figure 4 shows SEM images of the surface of the channels in the samples at the locations shown in Fig. 3. The channel sides can be clearly seen in the 3D printed sample (Fig. 4a). This indicates that the channel narrows towards the bottom. Some large features can be seen on the bottom of the channel (red circle, Fig. 4a). This is the result of incomplete mixing
<table><caption>**Table 1** Measured channel depths</caption><thead><tr><th>Location</th><th>Intended depth (mm)</th><th>Measured depth 3DP model (mm)</th><th>Measured depth Per-spex model (mm)</th></tr></thead><tbody><tr><td>1</td><td>0.3</td><td>0.2</td><td>0.2</td></tr><tr><td>2</td><td>0.3</td><td>0.3</td><td>0.3</td></tr><tr><td>3</td><td>0.3</td><td>0.3</td><td>0.2</td></tr><tr><td>4</td><td>0.3</td><td>0.1</td><td>0.2</td></tr><tr><td>5</td><td>0.3</td><td>0.3</td><td>0.4</td></tr><tr><td>6</td><td>0.3</td><td>0.5</td><td>0.2</td></tr><tr><td>7</td><td>0.3</td><td>0.3</td><td>0.3</td></tr><tr><td>8</td><td>0.3</td><td>0.2</td><td>0.2</td></tr><tr><td>9</td><td>0.3</td><td>0.4</td><td>0.2</td></tr><tr><td>10</td><td>0.3</td><td>0.3</td><td>0.3</td></tr><tr><td>11</td><td>0.3</td><td>0.4</td><td>0.3</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>**Fig. 3** Location of depth and width measurements and SEM images</figcaption></figure>

<table><caption>**Table 2** Measured channel widths</caption><thead><tr><th>Location</th><th>Intended width (mm)</th><th>Measured width 3DP model (mm)</th><th>Measured width Per-spex model (mm)</th></tr></thead><tbody><tr><td>1</td><td>1.0</td><td>0.95</td><td>1.05</td></tr><tr><td>2</td><td>1.0</td><td>0.91</td><td>1.06</td></tr><tr><td>3</td><td>1.0</td><td>0.92</td><td>1.03</td></tr><tr><td>4</td><td>1.0</td><td>0.97</td><td>1.04</td></tr><tr><td>5</td><td>1.0</td><td>1.02</td><td>1.06</td></tr><tr><td>6</td><td>1.5</td><td>1.42</td><td>1.54</td></tr><tr><td>7</td><td>1.5</td><td>1.48</td><td>1.57</td></tr><tr><td>8</td><td>1.5</td><td>1.45</td><td>1.57</td></tr><tr><td>9</td><td>0.5</td><td>0.39</td><td>0.50</td></tr><tr><td>10</td><td>0.5</td><td>0.39</td><td>0.51</td></tr><tr><td>11</td><td>0.5</td><td>0.42</td><td>0.50</td></tr></tbody></table>
<figure><img src="image_8.png" /><figcaption>**Fig. 4** a SEM image of the narrow channel of the 3D printed sample. Location is shown in the red dotted rectangle in Fig. 3. b SEM image of the narrow channel of the Perspex sample. Location is shown in the red dotted rectangle in Fig. 3. c SEM image of the bottom of the largest channel of the 3D printed sample. Location is shown in the blue rectangle in Fig. 3. d SEM image of the bottom of the largest channel of the Perspex sample. Location is shown in the blue dotted rectangle in Fig. 3. See text for more details</figcaption></figure>

of excess resin close to the build platform once a layer has been printed. The resin which has not been cleared away has subsequently been cured. Some roughness and microporosity can be seen at the edges of the channel (yellow dotted circle, Fig. 4 a). The likely reason for this is incomplete curing of the resin because of attenuation of the laser beam at the channel edges. More details about this phenomenon can be found in Ishutov et al. (2018a).

The Perspex channel (Fig. 4b) has much sharper channel edges than the 3D printed sample indicating more uniform channel width. There are some small bits of Perspex stuck to the bottom and the edges of the channel which affect the surface roughness. It is also possible to see circular marks on the base of the channel. These are the marks made by the cutting machine.

The bottom of the largest channel can be seen in Fig. 4c, d for the 3D printed and Perspex samples, respectively. The small ring like structures on the bottom of the 3D printed channel (yellow arrow, Fig. 4c) could be due to remaining pieces of supports left behind during the printing process or irregular curing of the resin related to interference patterns in the lasers. The Perspex sample is relatively smooth with some small bits of Perspex stuck to the bottom (yellow circle, Fig. 4d).

The Formlabs Form 2 printer cost $4500 (current UK price £3360 Formlabs 2018). Printing one 3D printed micromodel used 1.3 g of resin which cost $150 per litre giving a price per sample of $0.20 (about £0.15) excluding the cost of the printer. In comparison, the cost of Perspex for one micromodel was £3 and 1.5 hours of time on the machine required for cutting the sample cost £75 (about $100). Therefore, the 3D printed models are individually
far cheaper than the Perspex models. Additionally, the cost of the printer would be recovered
after about 45 models if we were to print them instead of using the cutting machine.

**2.3 Experimental Setup**

The visualisation cell consists of a transparent, Perspex holder containing the micromodel with a Schott A20960.1 LED light source above and a Baumer TXG08c colour CCD camera positioned below (Fig. 5). The camera resolution is 1028 × 772 pixels. The inlet is connected to a pump, which controls the flow rate of fluid into the cell. Fluid from the outlet is allowed to flow into a waste bottle. To encourage uniform mixing of fluid at the inlet a mesh is inserted between the fluid line and the micromodel inlet.

During flow experiments, the micromodel is placed upside down into the cell, i.e. with the inscribed channels face down at the bottom of the holder (Fig. 6) in order to provide a more efficient seal between the fluid lines and the inlet to the micromodel. A Perspex lid is then screwed on top.

Flow experiments were conducted by placing the micromodel into the visualisation cell, filling the model with deionised water then flowing deionised water containing a tracer through the inlet at a rate of 0.5 mL/min. Videos of the flow experiments were recorded using the camera. The tracer chosen was Winsor and Newton red drawing ink as it fulfilled our criteria of being dark enough to see in the videos, having a distinct variation in light intensity with tracer concentration and maintaining a consistent colour over time. We tested other tracers, Cobalt (II) Chloride Hexahydrate and food dye; however, these did not give a high enough colour contrast at different concentrations and the colour of the food dye was not consistent over time. The concentration of the tracer used in the experiments was 10 mL of ink to 40 mL of distilled water. Experiments were repeated three times.

The presence of ambient light, which could differ between experiments, could affect the recorded videos and calculated tracer concentrations. To minimise the impact of ambient

<figure><img src="image_9.png" /><figcaption>**Fig. 5** Schematic diagram of the experimental setup</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Fig. 6** Orientation of the micromodel within the visualisation cell. The micromodel is placed in the cell with the channels facing downwards to prevent fluid leaking out</figcaption></figure>

light on experimental results, we enclosed the flow cell in a cover to exclude ambient light. This ensured that light entering the camera was exclusively light which had passed through the sample.

## 2.4 Image Analysis

Tracer concentration can be estimated from the intensity of light transmitted through the sample containing dyed fluid. Videos were recorded of each experiment and light intensities were converted into concentration values for comparison with the numerical simulation results. Greyscale intensities were used for simplicity. Normalised intensity was found by dividing each video frame by the background image, which was taken from the first frame of the video before solute was present. This removed the effect of any uneven illumination and eliminated background features continuously present in the images.

A Wiener filter was applied to reduce random noise in the images. The MATLAB function *wiener2* (MathWorks 2018) was used with a 20 by 20 pixel neighbourhood. A mask was also applied to mask out areas of the sample not containing flow channels. This was generated by finding the areas of the image where pixel value changed by less than 50.

Once images had been normalised and filtered, resulting intensities were converted to concentrations. Using this method we are assuming that concentration variations in the *z* direction are small as measured concentrations are integrated values across the *z* dimension.

To carry out the conversion between light intensity and concentration, a calibration curve relating both properties is required. A separate calibration curve was created for the Perspex experiments and the 3D printed experiments as the imaging conditions were different for the each sample, due to differences in how light passes through the different materials.

The flow cell was filled with fluids at different tracer concentrations and videos were recorded each time. For each video the first 10 frames, where no solute was present, were averaged to produce an image of the average background frame. The last 10 frames were averaged to produce an image of the final frame to use for the analysis. The final image was divided by the background image to reduce the effect of uneven illumination. Care was taken to ensure that the experiment had been running long enough such that tracer concentration
<table><caption>**Table 3** Calibration data for the Perspex machined sample (M) and the 3D printed sample (3DP)</caption><thead><tr><th>Tracer concentration [–]</th><th>Normalised intensity M [–]</th><th>Normalised intensity 3DP [–]</th></tr></thead><tbody><tr><td>0.00</td><td>1.00</td><td>1.00</td></tr><tr><td>0.10</td><td>0.86</td><td>–</td></tr><tr><td>0.25</td><td>0.70</td><td>0.67</td></tr><tr><td>0.50</td><td>0.54</td><td>0.52</td></tr><tr><td>0.75</td><td>0.44</td><td>0.40</td></tr><tr><td>1.00</td><td>0.38</td><td>0.35</td></tr></tbody></table>

Normalised intensity is intensity divided by the background intensity at 0.00 tracer concentration. Tracer concentration is given as a fraction of the concentration used in the final experiments 10 mL dye / 40 mL of water

was constant. Average pixel intensity of the dyed fluid within a 14 × 28 pixel rectangle was measured at 11 locations throughout the image and the mean value of these measurements was used in the calibration (Table 3).

Light intensity is related to concentration using the Beer–Lambert law:

$$ \frac{I}{I_0} = 10^{-\mu C} \qquad (1) $$

where *I* is the measured intensity, *I*₀ is the intensity at *C* = 0, *C* is the solute concentration and *μ* = *ϵ*I where *ϵ* is the molar absorption coefficient of the solute and is a constant and *I* is the length of the path containing the solute which the light travels through. Our results are not well described by the Beer–Lambert law, with the curve being more pronounced than the analytical solution would suggest for both the 3D printed and Perspex sample (Fig. 7). This is possibly due to experimental reasons such as a nonlinear response from the camera at very light and very dark intensities. Therefore, instead of using the Beer–Lambert law we used linear interpolation in between data points (blue line in Fig. 7) to obtain concentration values for all intensities.

## 2.5 Governing Equations

Generally flow in our system can be described by the Navier–Stokes equation:

$$ \frac{\partial(\rho\mathbf{u})}{\partial t} = -\rho(\mathbf{u} \cdot \nabla)\mathbf{u} + \nabla \cdot (\mu\nabla\mathbf{u}) - \nabla P + \rho\mathbf{g}, \qquad (2) $$

where *ρ* is fluid density, *u* is velocity, *t* is time, *μ* is dynamic viscosity, *P* is pressure and *g* is acceleration due to gravity.

Assuming steady-state and neglecting the effect of gravity Eq. 2 becomes

$$ \mu \nabla^2 \mathbf{u} = \rho(\mathbf{u} \cdot \nabla)\mathbf{u} + \nabla P. \qquad (3) $$

and for incompressible flow we have the continuity equation

$$ \nabla \cdot \mathbf{u} = 0. \qquad (4) $$

Solute transport is governed by the classical advection-diffusion equation

$$ \frac{\partial C}{\partial t} + \nabla \cdot (\mathbf{u}C) = \nabla \cdot (\mathbf{D}\nabla C), \qquad (5) $$
<figure><img src="image_11.png" /><figcaption>**Fig. 7** Normalised intensity [–] versus tracer concentration [–] for the Perspex machined sample (M) and the 3D printed sample (3DP)</figcaption></figure>

where *C* is the solute concentration, *D* is the diffusion coefficient and *u* is the velocity field calculated using Eqs. 3, 4.

## 2.6 Numerical Modelling Setup

Numerical modelling was carried out with OpenFOAM 2.3.0 (OpenFOAM 2017). Open-FOAM is an open-source C++ library designed for undertaking computational fluid dynamics calculations using the finite volume method. First, we used the simpleFOAM solver to determine the velocity field **u** in the model domain. SimpleFOAM solves Eqs. 3 and 4 on a collocated Eulerian grid with a predictor-corrector strategy based on the semi-implicit method for pressure linked equations (SIMPLE) algorithm (Patankar and Spalding 1972). Then we used the scalarTransportFOAM solver to solve Eq. 5 using **u** as input.

A constant flow rate *Q* with uniform pressure boundary condition is applied at the inlet, a no-slip boundary condition is applied at the solid walls and a free flow condition with constant pressure is applied at the outlet.

Input parameters are shown in Table 4. The flow rate $Q = 8.33 \times 10^{-9} \text{ m}^3/\text{s}$ is equal that imposed in the experimental setup. To determine a reference velocity in the domain, we consider the total cross section of the three channels $A = 8.1 \times 10^{-7} \text{ m}^2$. The reference velocity $u_{ref}$ can then be calculated as $u_{ref} = Q/A = 0.0103 \text{ m/s}$. Using the average channel width as a reference length $L_{ref} = 1 \text{ mm}$, we obtain a Reynolds number $Re = 10.3$. The diffusion coefficient of the tracer in water is unknown so a reasonable value $D = 2 \times 10^{-9} \text{ m}^2/\text{s}$ has been chosen. This corresponds to a Péclet number $Pe = 5150$. In the incompressible case, the absolute value of pressure is not required, only a reference value is needed. Therefore, the outlet pressure is set to 0.

To mesh the computational domain, a 3D uniform Cartesian grid is first generated, and then all cells containing solid are removed and replaced by tetrahedral cells that match the solid boundaries using the OpenFoam snappyHexMesh utility (OpenFOAM 2017).

In order to improve the match between the experimental system and the simulations, channel widths used for the simulations were taken from average SEM measurements of the channel widths shown in Table 2. For the 3D printed model, channel widths used were 1.46
<table><caption>**Table 4** Numerical simulation input parameters used in the governing equations (Sect. 2.5)</caption><thead><tr><td>Parameter</td><td>Symbol</td><td>Value</td><td>Units</td></tr></thead><tbody><tr><td>Outlet pressure</td><td>P</td><td>0</td><td>Pa</td></tr><tr><td>Flow rate</td><td>Q</td><td>8.33 × 10<sup>-9</sup></td><td>m<sup>3</sup> s<sup>-1</sup></td></tr><tr><td>Dynamic viscosity</td><td>μ</td><td>10<sup>-3</sup></td><td>Pa s</td></tr><tr><td>Density</td><td>ρ</td><td>10<sup>3</sup></td><td>kg m<sup>-3</sup></td></tr><tr><td>Diffusion coefficient</td><td>D</td><td>2 × 10<sup>-9</sup></td><td>m<sup>2</sup> s<sup>-1</sup></td></tr></tbody></table>

<figure><img src="image_12.png" /><figcaption>**Fig. 8** Computational domain and finite volume grid</figcaption></figure>

mm, 0.95 mm and 0.4 mm and for the Perspex model the average channel widths used were 1.56 mm, 1.05 mm and 0.5 mm. The channel depth used for both simulations was 0.3 mm.

A grid convergence study was performed and convergence was reached for an initial Cartesian grid comprising 300 × 150 × 5 grid blocks. The final grid for the 3D printed geometry contained 64,090 cells and the final grid for the Perspex geometry contained 71,230 cells (Fig. 8).

## 3 Results

Figure 9 shows the processed images and corresponding concentration contours at four times (37 s, 37.5 s, 38 s, 38.5 s) during the experiment on the 3D printed micromodel.

The contoured image seems to indicate very low concentrations passing through the narrowest channel. It is likely that this is an artefact of the image processing of the channel that is only a few pixels wide. Image processing requires the frames to be filtered in order to produce a sufficiently clear image to generate contours from (otherwise the image is too noisy and the contours cannot be clearly seen). The filtering includes pixels located beyond the edge of the channels. These reduce the average pixel value at the edge of the channels leading to lower reported concentrations. In the smallest channel, this filtering completely obscures the results leading to lower concentrations throughout the whole of the channel.

Results of the numerical simulations using the 3D printed geometry can be seen on the right-hand side of Fig. 9. These plots show the normalised concentration, integrated over the depth of the micromodel to obtain a 2D image similar to the normalised intensity in the experiments. The inlet concentration used through time was matched to those seen in the experimental results (Fig. 11).
<figure><img src="image_13.png" /><figcaption>**Fig. 9** Normalised, greyscale images of the flow experiment in the 3D printed micromodel through time with corresponding numerical simulation results on the right-hand side. Concentration contours are plotted. Blue $C = 0.25$, green $C = 0.30$, red $C = 0.35$</figcaption></figure>

In the experimental results, fluid in both the top and the middle channel travels at approximately the same rate with contours in both channels reaching approximately the same distance away from the inlet at the same time (see $t = 37.5$ s). However, this is not the case for the numerical results where the contours travel fastest in the middle channel. A clear boundary layer can be seen close to the outlet where fluid from the middle channel and fluid from the top channel come together ($t = 38$ s), for both the experimental and numerical results.

Results from the experiments on the Perspex micromodel can be seen in Fig. 10. As with the 3D printed micromodel, flow in the narrowest channel is much lower in the experimental results compared to the numerical simulations, most likely due to the image processing methodology. In contrast to flow in the 3D printed micromodels, the tracer appears to travel fastest through the middle channel. A result that is also seen in the simulations. A clear boundary layer is evident close to the outlet.
<figure><img src="image_14.png" /><figcaption>**Fig. 10** Normalised, greyscale images of the flow experiment in the Perspex micromodel through time with corresponding numerical simulation results on the right-hand side. Concentration contours are plotted. Blue $C = 0.25$, green $C = 0.30$, red $C = 0.35$</figcaption></figure>

There is a good correspondence between results of numerical simulations and experimental results from the Perspex micromodel (Fig. 10). This suggests that the physical dimensions of the channel width are close to the mean channel width used in the simulation. This is supported by relatively small standard deviations of channel widths from the mean channel widths for the Perspex model (Table 5). There is a larger difference between the numerical simulations and experimental results for the case of the 3D printed micromodel compared to the Perspex micromodel. The variation in channel width from the mean is greater for the 3D printed model (Table 5). Also the standard deviation of channel depth from the depth used in the simulations is greater for the 3D printed model compared to the Perspex model (0.11 mm vs. 0.08 mm, respectively).

Numerical simulations of the 3D printed experiments (Fig. 9) show contours moving faster at subsequent times compared to the speed of contours in the simulations with the
<table><caption>**Table 5** Standard deviation of channel widths from the mean</caption><thead><tr><th></th><th>3DP (mm)</th><th>Perspex (mm)</th></tr></thead><tbody><tr><td>Narrow channel</td><td>0.0195</td><td>0.0101</td></tr><tr><td>Medium channel</td><td>0.0326</td><td>0.0159</td></tr><tr><td>Large channel</td><td>0.0444</td><td>0.0112</td></tr></tbody></table>

Actual channel widths are shown in Table 2

<figure><img src="image_15.png" /><figcaption>**Fig. 11** Concentration values used as the inlet boundary condition in the numerical simulations. The 3DP curve shows the inlet concentration in the 3D printed micromodel experiment. This was used in the numerical simulation shown in Fig. 9. The M curve shows the inlet concentration in the Perspex micromodel experiment. This is used in the numerical simulation shown in Fig. 10</figcaption></figure>

Perspex geometry (Fig. 10). Experimental results also show the contours moving faster for the 3D printed micromodel compared to the Perspex micromodel. This is likely to be as a consequence of the smaller channel dimensions in the 3D printed micromodel.

In absolute terms, the blue contour ($C = 0.25$) in the Perspex model is always ahead of the blue contour in the 3D printed model. However, the speed that the blue contour is travelling at is slower for the Perspex model (note the difference between blue contours at subsequent timesteps for the same experiment). The green contour ($C = 0.30$) for the Perspex model is always behind the green contour for the 3D printed model and is also travelling more slowly. The difference in the absolute location of the blue contour between the Perspex and 3D printed model can be explained by the differing inlet concentration shown in Fig. 11.

## 4 Discussion

Our results show that there is some difference between flow and transport processes observed in a 3D printed micromodel compared to flow and transport in a Perspex micromodel made to the same specifications. This is likely due to inconsistencies between the geometry of the 3D printed micromodel and the Perspex micromodel. SEM imaging shows that there is a greater variation in channel width and depth for the 3D printed model compared to the Perspex model. SEM imaging also shows some unintended microporosity and a non-uniform
channel width with depth for the 3D printed sample. These factors will affect the flow and reduce the match between the experimental and numerical results for the 3D printed sample. Discrepancies between intended geometries and printed geometries have been also found in other studies (e.g. Suzuki et al. 2017; Head and Vanorio 2016; Wicker et al. 2005).

Differences between the experimental and numerical results could also be caused by errors in calculating the concentration contours from the images. As mentioned in Sect. 3, the image processing methodology leads to a slight reduction in concentration at the edges of the channels which has led to the flow in the smallest channel being different between the numerical and experimental results. It is also possible that this processing error could obscure an actual difference between the simulations and the experiments.

Errors in the experimental setup, such as slight variations in flow rate during the experiment, could also lead to differences between experimental results and numerical results. Experiments were repeated; however, it was difficult to compare results exactly between experiments as the time when the tracer reached the cell inlet was hard to constrain. Work is ongoing to improve this by modifying the cell design so there is more precise control at the inlet.

To gain a better understanding of how variations in geometry affect flow and transport, the next step would be to generate more detailed images of the sample using a technique such as surface profilometry. Results could be used to produce a mesh for the simulations which is much closer to the geometry of the physical model. This would give a clearer indication of how much of the variation in flow patterns was caused by inconsistencies in micromodel geometry.

For 3D printed micromodels to be more useful in the future for flow and transport applications, we need to improve the geometrical accuracy of the printed micromodels. One way of doing this could involve modifying the way prints are made in terms of the orientation in which models are created within the printer (e.g. Ishutov et al. 2018b). Ultimately though, improvements in print accuracy are likely to come about by advances in 3D printing technology. The current interest in 3D printing would indicate that this is likely to happen in the near future (Ishutov et al. 2018b).

Future work should also be undertaken into establishing how repeatable 3D printing is when printing micromodels. This would help determine how much variation in flow and transport is because of errors in micromodel geometry. Variations in geometry between prints need to be as small as possible if we are to carry out repeatable experiments. It would also be interesting to assess repeatability of prints between different types of printers.

This paper looks solely at single-phase flow. To investigate flow processes involving two or more phases, work needs to be carried out to understand other properties such as the wettability of the 3D printed micromodels which can be influenced by things such as the surface roughness and chemistry. In order to conduct repeatable experiments, variation in these properties between samples should be minimised.

# 5 Conclusion

We have conducted single-phase flow and transport experiments on a 3D printed micromodel. Comparison of our results with flow experiments on a micromodel cut from Perspex shows some discrepancies in terms of the speed and location of tracer flow in a relatively simple geometry. Numerical simulations suggest that some of this variation in flow patterns is caused by differences in channel width and depth between the 3D printed and Perspex micromodels.
There is a good match between numerical simulation results and experimental results for the Perspex micromodel. The match between experimental and numerical results for the 3D printed model is not as good as for the Perspex model. This is possibly due to more variable geometry, increased surface roughness and unintended microporosity in the 3D printed sample.

The use of 3D printed micromodels in pore-scale flow and transport investigations looks promising as our experimental results were not significantly different from the experimental results using the Perspex micromodel and numerical simulations. Using our setup, we were able to visualise the flow and observe mixing within the porespace which will be very useful for investigating pore-scale flow phenomena. However, to be able to use 3D printed micro-models we need to be able to produce the required geometry accurately and we need minimal variations in properties (e.g. geometry, wettability, etc.) between different prints. This could potentially be achieved by improving printing methodologies, such as by printing models in different orientations, or more likely by general improvements in 3D printer technology and accuracy.

In the meantime, future work should look at establishing how accurate and repeatable current micromodels are and how much variation between prints is allowable before inconsistencies make a noticeable difference to simulation results. This characterisation of 3D printed models is an important step in their use for flow and transport experiments simulating natural rocks.

**Acknowledgements** We thank Energi Simulation for supporting the chairs held by Sebastian Geiger and Eric Mackay. We thank EPSRC for supporting Julien Maes' position (Grant Number EP/P031307/1). We are grateful to Mohamed Arab and Jim Buckmann for help with the collection and interpretation of the SEM data.

**Open Access** This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

## References

Blunt, M.J.: Multiphase Flow in Permeable Media: A Pore-Scale Perspective. Cambridge University Press, Cambridge (2017)

Buchgraber, M., Al-Dossary, M., Ross, C., Kovscek, A.R.: Creation of a dual-porosity micromodel for pore-level visualization of multiphase flow. J. Pet. Sci. Eng. **86**, 27–38 (2012)

Chan, H.N., Shu, Y., Xiong, B., Chen, Y., Chen, Y., Tian, Q., Michael, S.A., Shen, B., Wu, H.: Simple, cost-effective 3d printed microfluidic components for disposable, point-of-care colorimetric analysis. ACS Sens. **1**(3), 227–234 (2015)

Comina, G., Suska, A., Filippini, D.: Low cost lab-on-a-chip prototyping with a consumer grade 3d printer. Lab Chip **14**(16), 2978–2982 (2014)

Formlabs. https://formlabs.com (2018)

Gerami, A., Mostaghimi, P., Armstrong, R.T., Zamani, A., Warkiani, M.E.: A microfluidic framework for studying relative permeability in coal. Int. J. Coal Geol. **159**, 183–193 (2016)

Head, D., Vanorio, T.: Effects of changes in rock microstructures on permeability: 3-d printing investigation. Geophys. Res. Lett. **43**(14), 7494–7502 (2016)

Ishutov, S., Hasiuk, F.J., Fullmer, S.M., Buono, A.S., Gray, J.N., Harding, C.: Resurrection of a reservoir sandstone from tomographic data using three-dimensional printing. Am. Assoc. Pet. Geol. Bull. **101**(9), 1425–1443 (2017)

Ishutov, S., Hasiuk, F.J., Jobe, D., Agar, S.: Using resin-based 3d printing to build geometrically accurate proxies of porous sedimentary rocks. Groundwater **56**(3), 482–490 (2018a)
Ishutov, S., Jobe, T.D., Zhang, S., Gonzalez, M., Agar, S.M., Hasiuk, F.J., Watson, F., Geiger, S., Mackay, E., Chalaturnyk, R.: Three-dimensional printing for geoscience: fundamental research, education, and applications for the petroleum industry. Am. Assoc. Pet. Geol. Bull. **102**(1), 1–26 (2018b)

Karadimitriou, N., Musterd, M., Kleingeld, P., Kreutzer, M., Hassanizadeh, S., Joekar-Niasar, V.: On the fabrication of pdms micromodels by rapid prototyping, and their use in two-phase flow studies. Water Resour. Res. **49**(4), 2056–2067 (2013)

Kitson, P.J., Rosnes, M.H., Sans, V., Dragone, V., Cronin, L.: Configurable 3d-printed millifluidic and microfluidic lab on a chipreactionware devices. Lab Chip **12**(18), 3267–3271 (2012)

MathWorks. https://uk.mathworks.com/help/images/ref/wiener2.html (2018)

OpenFOAM. http://www.openfoam.org/ (2017)

Patankar, S.V., Spalding, D.B.: A calculation procedure for heat, mass and momentum transfer in three-dimensional parabolic flows. Int. J. Heat Mass Transf. **15**(10), 1787–1806 (1972)

Porter, M.L., Jiménez-Martínez, J., Martinez, R., McCulloch, Q., Carey, J.W., Viswanathan, H.S.: Geo-material microfluidics at reservoir conditions for subsurface energy resource applications. Lab Chip **15**(20), 4044–4053 (2015)

Rogers, C.I., Qaderi, K., Woolley, A.T., Nordin, G.P.: 3d printed microfluidic devices with integrated valves. Biomicrofluidics **9**(1), 016501 (2015)

Shallan, A.I., Smejkal, P., Corban, M., Guijt, R.M., Breadmore, M.C.: Cost-effective three-dimensional printing of visibly transparent microchips within minutes. Anal. Chem. **86**(6), 3124–3130 (2014)

Suzuki, A., Watanabe, N., Li, K., Horne, R.N.: Fracture network created by 3-d printer and its validation using ct images. Water Resour. Res. **53**(7), 6330–6339 (2017)

Waheed, S., Cabot, J.M., Macdonald, N.P., Lewis, T., Guijt, R.M., Paull, B., Breadmore, M.C.: 3d printed microfluidic devices: enablers and barriers. Lab Chip **16**(11), 1993–2013 (2016)

Wicker, R.B., Ranade, A.V., Medina, F., Palmer, J.A.: Embedded micro-channel fabrication using line-scan stereolithography. Assembly Automation **25**(4), 316–329 (2005)

Zhao, B., MacMinn, C.W., Juanes, R.: Wettability control on multiphase flow in patterned microfluidics. Proc. Natl. Acad. Sci. **113**(37), 10251–10256 (2016)

Zhou, J., Ellis, A.V., Voelcker, N.H.: Recent developments in pdms surface modification for microfluidic devices. Electrophoresis **31**(1), 2–16 (2010)