Get drawdown relationships

Calculate Dii, Dij (or PHIii, PHIij) for groups of wells

get_drawdown_relationships(wells, aquifer, group_column, weights_column)

Arguments

wells	Wells data.frame or tibble
aquifer	Aquifer object
group_column	Name of column to identify groups of wells (no quotes)
weights_column	Name of column to identify weights of well pumping (no quotes). Required -- create column equal to 1 if equal weights desired.

Value

Returns a data.frame containing the drawdown due to unit pumping for all combinations of groups in the group_column. The variable, var, is given as D_i_j (confined aquifers), or PHI_i_j (unconfined aquifers), which represent the (weighted) average drawdown at wells i due to (weighted) unit pumping from wells j.

Examples

# define aquifer
bounds_df <- data.frame(bound_type=c("CH","NF","NF","NF"),m=c(Inf,0,Inf,0),b=c(0,1000,1000,0))
aquifer_unconfined <- define_aquifer("unconfined",1e-3,bounds=bounds_df,h0=100)

# define wells and well images
library(dplyr)
set.seed(30)
wells_df <- data.frame(x=runif(8,0,1000),y=runif(8,0,1000),diam=1,R=1000) %>%
  mutate(country=factor(y>500,levels=c(FALSE,TRUE),labels=c("A","B"))) %>%
  group_by(country) %>%
  mutate(weights=1,Q=1/n()) %>% group_by()
wells <- define_wells(wells_df) %>% generate_image_wells(aquifer_unconfined)

get_drawdown_relationships(wells, aquifer_unconfined, group_column=country, weights_column=weights)
#> # A tibble: 4 x 4
#>   var       pot units     description                                           
#>   <chr>   <dbl> <chr>     <chr>                                                 
#> 1 PHI_A_A 1190. [m^2/cum… Weighted effect of pumping [cumec] in group A on disc…
#> 2 PHI_A_B  125. [m^2/cum… Weighted effect of pumping [cumec] in group B on disc…
#> 3 PHI_B_A  125. [m^2/cum… Weighted effect of pumping [cumec] in group A on disc…
#> 4 PHI_B_B  862. [m^2/cum… Weighted effect of pumping [cumec] in group B on disc…

library(ggplot2)
ggplot() +
  geom_segment(data=aquifer_unconfined$bounds,aes(x1,y1,xend=x2,yend=y2,color=bound_type)) +
  geom_abline(slope=0,intercept=500,linetype="dashed") +
  geom_point(data=wells %>% filter(wID==orig_wID),aes(x,y,fill=country),shape=21) +
  coord_equal()