README.md

partFilt README

Dan Puddephatt 2016-12-10

The R library partFilt is used to calculate particulate filtration based on the methods proposed by J.E. Tobiason and C.R. O'Melia in Physicochemical Aspects of Particle Removeal in Depth Filtration (1988). The library can be installed with the following code:

library(devtools)
install_github("dpphat/partFilt")

I will not go into the details of the physics and chemistry involved in particle filtration that are discussed in Tobiason and O'Melia (1988) suffice it to say that suspended solids interact with soil particles (collector particles) through size exclusion, gravitational attraction (sedimentation), and chemical attraction and repulsion. The chemical interaction term (ALPHA) is a fitted parameter that needs to be investigated through calibration. The units of measurement in partFilt are based on those used in the paper and include "cm", "g", "s". There is an option to specify the units of hydraulic conductivity to be either "cm/s" or "ft/day". Length units can also be specified as either "cm", "m", or "ft".

The single-collector efficiency (eta_parm) is the efficiency of a single soil particle to attract a suspended solid from flowing groundwater. The eta_parm() function is used to calculate the single-collector efficiency. Figure 4 from Tobiason and O'Melia (1988) is prepared using the following code:

library(scales)
ETA <- eta_parm(Kx = 0.14 * 0.4, 
                TEMP = 25, 
                GRAD = 1, 
                POR = 0.4, 
                DIA_SOIL = 0.04, 
                PART_DENS = 1.05, 
                DIA_SUSP = 10^(seq(-2, 2, 0.1) - 4))

ETA <- tibble(PART_SIZE = 10^(seq(-2, 2, 0.1)), 
              ETA = ETA)

ETA %>% ggplot(aes(x = PART_SIZE, y = ETA))+
        geom_line(size = 1) +
        scale_x_log10(name = "Particle Size (µm)", 
                      limit = c(0.01, 100), 
                      breaks = 10^(-10:10), 
                      minor_breaks = c(1:9 %o% 10^(-10:10)), 
                      labels = scales::trans_format("log10", math_format(10^.x))) +
        scale_y_log10(name = expression(eta), 
                      limit = c(10^-4, 1), 
                      breaks = 10^(-10:10), 
                      minor_breaks = c(1:9 %o% 10^(-10:10)), 
                      labels = scales::trans_format("log10", math_format(10^.x))) +
        theme(
                      plot.title          = element_text(face = "bold", hjust = 0, size = 10),
                      legend.position     = "none", 
                      strip.background    = element_rect(colour = "white", fill = "white", size = 0),
                      strip.text          = element_text(face = "bold"), 
                      panel.border        = element_rect(colour = "black", size = 1, fill = NA), 
                      panel.background    = element_rect(fill = NA), 
                      panel.grid.major    = element_line(colour = rgb(217, 217, 217, 255, maxColorValue = 255), size = 0.25), 
                      panel.grid.minor    = element_line(colour = rgb(217, 217, 217, 255, maxColorValue = 255), size = 0.25), 
                      axis.text.x         = element_text(size = 9),
                      axis.title.x        = element_text(size = 9),
                      strip.text.x        = element_text(size = 9, face = "bold", vjust = 0, hjust = 0),
                      strip.text.y        = element_text(size = 9),
                      axis.text.y         = element_text(size = 9), 
                      axis.title.y        = element_text(size = 9)
                      )

When there is no chemical interaction (i.e., ALPHA == 1) the minimum collector efficiency occurs when the suspended sediment is:

ETA %>% filter(ETA == min(ETA))
## # A tibble: 1 × 2
##   PART_SIZE          ETA
##       <dbl>        <dbl>
## 1         1 0.0007073822

The calculated reduction values presented on Figure 13 when the ALPHA parameter is 0.01 (i.e., there is a chemical repulsion acting on suspended particles) can be replicated with the following:

CC0 <- cc0_filt(LENGTH = 60, 
         LENGTH_UNIT = "cm", 
         alpha = 0.01, 
         TEMP = 25, 
         POR = 0.4, 
         Kx = 5 * 100 / 3600 * 0.4, 
         GRAD = 1, 
         PART_DENS = 1.055, 
         DIA_SOIL = 0.05, 
         DIA_SUSP = 10^(seq(-2, 3, 0.01) - 4))

CC0 <- tibble(CC0 = CC0, 
              DIA_SUSP = 10^(seq(-2, 3, 0.01)))      

CC0 %>% ggplot(aes(x = DIA_SUSP, y = CC0)) +      
            geom_line(size = 1) +
            scale_x_log10(name = "Particle Radii (µm)", 
                          breaks = 10^(-10:10), 
                          minor_breaks = c(1:9 %o% 10^(-10:10)), 
                          labels = scales::trans_format("log10", math_format(10^.x)), 
                          expand = c(0, 0)) +
            scale_y_log10(name = "Fraction Remaining (C / C0)", 
                          limit = c(1e-5, 1), 
                          breaks = 10^(-10:10), 
                          minor_breaks = c(1:9 %o% 10^(-10:10)), 
                          labels = scales::trans_format("log10", math_format(10^.x)), 
                          expand = c(0, 0)) +
            theme(
                          plot.title          = element_text(face = "bold", hjust = 0, size = 10),
                          legend.position     = "none", 
                          strip.background    = element_rect(colour = "white", fill = "white", size = 0),
                          strip.text          = element_text(face = "bold"), 
                          panel.border        = element_rect(colour = "black", size = 1, fill = NA), 
                          panel.background    = element_rect(fill = NA), 
                          panel.grid.major    = element_line(colour = rgb(217, 217, 217, 255, maxColorValue = 255), size = 0.25), 
                          panel.grid.minor    = element_line(colour = rgb(217, 217, 217, 255, maxColorValue = 255), size = 0.25), 
                          axis.text.x         = element_text(size = 9),
                          axis.title.x        = element_text(size = 9),
                          strip.text.x        = element_text(size = 9, face = "bold", vjust = 0, hjust = 0),
                          strip.text.y        = element_text(size = 9),
                          axis.text.y         = element_text(size = 9), 
                          axis.title.y        = element_text(size = 9)
                          )        
## Warning: Removed 68 rows containing missing values (geom_path).

We can use the command L_filt to determine the particle travel distance required to achieve a suspended solid reduction to meet standards. In the following example, I will assume that we want to reduce suspended solids concentrations to 1e-10 of the initial concentrations. I will assume that the sediment has a hydraulic conductivity value of 0.02 cm/s for coarse sands. The groundwater temperature is assumed to be 10 degrees Celcius with is common for temperate climates. Porosity will also be assigned as 0.21. We will investigate a range of hydraulic gradients that range from quite small to quite large, spanning 6 orders of magnitude.

LENGTHS <- L_filt(C_C0 = 1e-10, 
                  alpha = 0.01, 
                  Kx = 0.02, 
                  GRAD = 10^seq(-4, -2, 0.1), 
                  POR = 0.21, 
                  TEMP = 10)

LENGTHS <- tibble(GRAD = 10^seq(-4, -2, 0.1), 
                  LENGTHS = LENGTHS)

LENGTHS %>% ggplot(aes(x = GRAD, y = LENGTHS)) + 
            geom_line(size = 1) +
            scale_x_log10(name = "Hydraulic Gradients (cm / cm)", 
                          breaks = 10^(-10:10), 
                          minor_breaks = c(1:9 %o% 10^(-10:10)), 
                          labels = scales::trans_format("log10", math_format(10^.x))) +
            scale_y_log10(name = "Filtration Length (cm)", 
                          breaks = 10^(-10:10), 
                          minor_breaks = c(1:9 %o% 10^(-10:10)), 
                          labels = scales::trans_format("log10", math_format(10^.x))) +
            theme(
                          plot.title          = element_text(face = "bold", hjust = 0, size = 10),
                          legend.position     = "none", 
                          strip.background    = element_rect(colour = "white", fill = "white", size = 0),
                          strip.text          = element_text(face = "bold"), 
                          panel.border        = element_rect(colour = "black", size = 1, fill = NA), 
                          panel.background    = element_rect(fill = NA), 
                          panel.grid.major    = element_line(colour = rgb(217, 217, 217, 255, maxColorValue = 255), size = 0.25), 
                          panel.grid.minor    = element_line(colour = rgb(217, 217, 217, 255, maxColorValue = 255), size = 0.25), 
                          axis.text.x         = element_text(size = 9),
                          axis.title.x        = element_text(size = 9),
                          strip.text.x        = element_text(size = 9, face = "bold", vjust = 0, hjust = 0),
                          strip.text.y        = element_text(size = 9),
                          axis.text.y         = element_text(size = 9), 
                          axis.title.y        = element_text(size = 9)
                          )



dpphat/partFilt documentation built on May 15, 2019, 1:47 p.m.