LBPM Tutorial, Step 8. Simulating Water Flooding (Part I)

In the previous tutorial initial conditions were established based on morphological drainage. We now consider how to use this initial configuration of fluids to simulate a water-flooding process. In experimental setting, water is pumped into the sample at a particular flow rate to displace oil from the pore space. We will use a flux boundary condition to mimic this same basic approach (see McClure et al. for details).

The flux boundary condition is determined by using a consistency condition derived for lattice Boltzmann methods to apply a constant pressure over the inlet such that a desired volumetric flow rate is obtained. Measured at the inlet, the volumetric flow rate is exact, since the boundary condition is derived analytically. The unit test TestFluxBC verifies that this an exact match is achieved at the boundary. However, if spurious currents arise in the vicinity of the boundary the volumetric flow rate measured through the system will tend not match the boundary value. This will also effect other velocity boundary conditions for the LBM. The best way to overcome this is by using a larger reservoir at the inlet. An inlet layer that is 16-voxels deep would not be an unreasonable choice. For this tutorial we will use a 4-voxel layer, which is the minimum allowed by LBPM when running simulations with non-periodic boundary conditions.

To simulate the water-flooding, we will use a single GPU on a smaller sub-domain, simply so that the system will be flooded more quickly. The input file input.db for this case will be

Domain {
   Filename = "mask_water_flooded_water_and_oil.raw.morphdrain.raw"  
   ReadType = "8bit"  // data type
   nproc = 1, 1, 1     // process grid
   n = 300, 297, 500   // sub-domain size
   N = 300, 297, 500  // size of original image
   voxel_length = 7.0  // voxel length (in microns)
   ReadValues = -1, 0, 1, 2  // labels within the original image
   WriteValues = -1, 0, 1, 2 // associated labels to be used by LBPM
   BC = 4                // boundary condition type (0 for flux)
   Sw = 0.35
}
MRT {
   tau = 1.0
   F = 0.0, 0.0, 1.0e-5
   timestepMax = 20000
   tolerance = 0.01
}
Color {
    tauA = 0.7;             // relaxation time for fluid A (labeled as "1")       
    tauB = 0.7;             // relaxation time for fluid B (labeled as "2") 
    rhoA   = 1.0;           // density for fluid A (in lattice units)
    rhoB   = 1.0;           // density for fluid B (in lattice units)
    alpha = 1e-3;           // controls the surface tension
    beta  = 0.95;           // controls the interface width 
    F = 0, 0, 0             // controls the external force
    Restart = false         // initialize simulation from restart file?
    timestepMax = 1000    // maximum number of timesteps to perform before exit
    ComponentLabels = 0, -1     // number of immobile component labels in the input image
    ComponentAffinity = -1.0, -0.9  // wetting condition for each immobile component
    flux = -10.0               // volumetric flux at the z-inlet in voxels per timestep
}
Analysis {
    analysis_interval = 1000         // Frequency to perform analysis
    visualization_interval = 100000  // Frequency to write visualization data
    restart_interval = 1000000       // Frequency to write restart data
    restart_file = "Restart"         // Filename to use for restart file (will append rank)
    N_threads    = 4                 // Number of threads to use for analysis
    load_balance = "independent"     // Load balance method to use: "none", "default", "independent"
}
Visualization {
}

In this case, changes have been made to the Domain section of the input database so that results from the morphological drainage are used as input. Several lines must be updated to use the data written by lbpm_morphdrain_pp instead of the original image:

• Filename = "mask_water_flooded_water_and_oil.raw.morphdrain.raw"
• ReadType = "8bit"
• ReadValues = 0, 1, 2
• N = 300, 297, 500

Note that the convention for flux is relative to the flux of oil. This means that a negative flux corresponds to water injection and positive flux corresponds to oil injection.

For the water-flooding simulation performed by lbpm_color_simulator, note that we have chosen timestepMax = 1000. This is not enough time for any interesting dynamics to occur. It is sufficient to

1. Evaluate the performance
2. Verify that you haven't made a mistake

In this case we see the performance reported in million fluid lattice updates per second (MFLUPs)

********************************************************
Running Color LBM	
********************************************************
MPI rank=0 will use GPU ID 0 / 4 
voxel length = 7.000000 micron 
voxel length = 7.000000 micron 
Read input media... 
Initialize from segmented data: solid=0, NWP=1, WP=2 
Media porosity = 0.130172 
Initialized solid phase -- Converting to Signed Distance function 
Domain set.
Create ScaLBL_Communicator 
Set up memory efficient layout, 6289409 | 6289440 | 45329596 
Allocating distributions 
Setting up device map and neighbor list 
Component labels: 1 
   label=0, affinity=-1.000000, volume fraction==0.870766
Initializing distributions 
Initializing phase field 
********************************************************
No. of timesteps: 1000 
Affinities - rank 0:
Main: 0
Thread 1: 16
....inlet pressure=1.014592 
-------------------------------------------------------------------
********************************************************
CPU time = 0.035938 
Lattice update rate (per core)= 175.009294 MLUPS 
Lattice update rate (total)= 175.009294 MLUPS 
********************************************************

On a V100 GPU, 175 MFLUPS is not especially fast for lbpm_color_simulator -- one can achieve as high as 550 MLUPS per V100 GPU provided that sufficiently large sub-domains are used. In this case the sub-domain is probalby large enough. However, getting full performance often requires an exercise in tuning the flags provided to mpirun. Frequently data locality or memory affinity issues are the cause of any loss in performance. Binding processes to cores can sometimes resolve this issue (-bind-to core in openMPI) but sometimes more advanced settings are needed (e.g. setting MCA parameters in OpenMPI). Since the optimal settings are particular to your environment, we leave this exercise to the reader.

To verify that the simulation has not failed, read the output file timelog.csv into python or manually inspect in from the command line

cat timelog.csv

Increasing the number of timesteps will allow for simulation of a complete water flooding. To do this, you will need to

• get the porosity of the sample (printed out by LBPM)
• multiply the porosity by the number of voxels in your simulation to get the size of a pore volume (in number voxels)
• based on flux, determine the estimated timestepMax needed to inject one pore volumes into the system
• update input.db and re-run the simulation.

Proceed to next step in the tutorial