Nothing
## ---- include = FALSE---------------------------------------------------------
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>"
)
## ---- fig.width = 6, fig.height= 4, fig.align = "center", warning = FALSE-----
library(Rmbal)
library(Rrelperm)
library(pracma)
library(minpack.lm)
library(ggplot2)
library(dplyr)
library(magrittr)
p_pvt <- c(3700, 3650, 3400, 3100, 2800, 2500, 2200, 1900, 1600, 1300, 1000, 700,
600, 400)
Bo <- c(10.057, 2.417, 2.192, 1.916, 1.736, 1.617, 1.504, 1.416, 1.326, 1.268, 1.205,
1.149, 1.131, 1.093)
Rv <- c(84.11765, 84.11765, 70.5, 56.2, 46.5, 39.5, 33.8, 29.9, 27.3, 25.5, 25.9,
28.3, 29.8, 33.5) / 1e6
Rs <- c(11566, 2378, 2010, 1569, 1272, 1067, 873, 719, 565, 461, 349, 249, 218, 141)
Bg <- c(0.87, 0.88, 0.92, 0.99, 1.08, 1.20, 1.35, 1.56, 1.85, 2.28, 2.95, 4.09,
4.68, 6.53) / 1000
cw <- 3e-6
Bwi <- 10.05
Bw <- Bwi * exp(cw * (p_pvt[1] - p_pvt))
muo <- c(0.0612, 0.062, 0.1338, 0.1826, 0.2354, 0.3001, 0.3764, 0.4781, 0.6041,
0.7746, 1.0295, 1.358, 1.855, 2.500)
mug <- c(0.0612, 0.062, 0.0554, 0.0436, 0.0368, 0.0308, 0.0261, 0.0222, 0.0191,
0.0166, 0.0148, 0.0135, 0.0125, 0.0115)
muw <- rep(0.25, length(p_pvt))
pvt_table <- data.frame(p = p_pvt, Bo = Bo, Rs = Rs, Rv = Rv, Bg = Bg,
Bw = Bw, muo = muo, mug = mug, muw = muw)
p <- c(3700, 3650, 3400, 3100, 2800, 2500, 2200, 1900, 1600, 1300, 1000, 700,
600)
We <- rep(0, length.out = length(p))
Np <- c(0, 28.6, 93, 231, 270, 379, 481, 517.2, 549, 580, 675, 755, 803) *1e3
Gp <- c(0, 0.34, 1.2, 3.3, 4.3, 6.6, 9.1, 10.5, 12, 12.8, 16.4, 19.1, 20.5) * 1e9
Wp <- rep(0, length.out = length(p))
Wi <- rep(0, length.out = length(p))
wf <- rep(1, length.out = length(p))
mbal_optim_gas_lst <- mbal_optim_param_gas(input_unit = "Field", output_unit = "Field",
unknown_param = "G", aquifer_model = NULL,
phi = 0.1, swi = 0.2, Np = Np,
Gp = Gp, Wp = Wp, Wi = Wi, We = We,
pd = 3650, p = p, pvt = pvt_table, M = 0,
cf = 2e-6, wf = wf, sgrw = 0.15)
time_lst <- mbal_time(c(1:length(p)), "year")
# a number of plots will be automatically generated for quality check
optim_results <- mbal_optim_gas(mbal_optim_gas_lst, time_lst)
glimpse(optim_results)
## ---- fig.width = 6, fig.height= 4, fig.align = "center", warning = FALSE-----
mbal_results <- mbal_perform_gas(optim_results, time_lst)
mbal_results
## ----fig.align="center", fig.height=4, fig.width=6, warning=FALSE-------------
# Step I: generating a set of pseudo relative permeability curves using
# laboratory 'Kr' values
sg_lab <- c(0.05, 0.152, 0.248, 0.352, 0.448, 0.552, 0.65)
krg_lab <- c(0, 0.05, 0.09, 0.18, 0.3, 0.5, 1)
kro_lab <- c(1, 0.6, 0.35, 0.13, 0.04, 0.01, 0)
swcrit_lab <- 0.2
sgcrit_lab <- 0.05
sorgr_lab <- 0.15
fun_kr <- function(x, swcrit, sgcrit, sorg, sg, krg, kro) {
kr_est <- Rrelperm::kr2p_gl(SWCON = swcrit, SOIRG = sorg, SORG = sorg,
SGCON = sgcrit, SGCRIT = sgcrit, KRGCL = 1,
KROGCG = 1, NG = x[1], NOG = x[2], NP = 101)
krg_est_sub <- approx(x = kr_est[,1], y = kr_est[,3], xout = sg, rule = 2)$y
kro_est_sub <- approx(x = kr_est[,1], y = kr_est[,4], xout = sg, rule = 2)$y
error <- (krg - krg_est_sub) ^ 2 + (kro - kro_est_sub) ^ 2
return(error)
}
par <- c(2, 2)
opt_results <- minpack.lm::nls.lm(par = par, fn = fun_kr, swcrit = swcrit_lab,
sgcrit = sgcrit_lab, sorg = sorgr_lab, sg = sg_lab,
krg = krg_lab, kro = kro_lab,
lower = c(0.1,0.1), upper = c(10,10))
sol <- opt_results$par
sol
rel_perm <- as.data.frame(Rrelperm::kr2p_gl(SWCON = swcrit_lab, SOIRG = sorgr_lab,
SORG = sorgr_lab, SGCON = sgcrit_lab,
SGCRIT = sgcrit_lab, KRGCL = 1,
KROGCG = 1, NG = sol[1], NOG = sol[2],
NP = 101))
colnames(rel_perm) <- c("Sg", "Sl", "Krg", "Krog")
p_forecast <- p
wf_forecast <- wf
time_lst_forecast <- mbal_time(c(1:length(p_forecast)), "year")
forecast_lst <- mbal_forecast_param_gas(input_unit = "Field", output_unit = "Field",
G = 2.41e10, phi = 0.1, swi = 0.2,
pd = 3650, p = p_forecast, pvt = pvt_table,
M = 0, cf = 2e-6, wf = wf_forecast,
rel_perm = rel_perm)
glimpse(forecast_lst)
forecast_results <- mbal_forecast_gas(forecast_lst, time_lst_forecast)
forecast_results
p1 <- forecast_results %>% ggplot(aes(`P (psia)`, SOg, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = mbal_results, aes(`P (psia)`, SOg, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("Condensate Saturation Plot") +
theme_bw()
p1
p2 <- forecast_results %>% ggplot(aes(`P (psia)`, `GOR (SCF/STB)`, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = mbal_results, aes(`P (psia)`, `GOR (SCF/STB)`, color = "Field")) +
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("GOR Plot") +
theme_bw()
p2
p3 <- forecast_results %>% ggplot(aes(`P (psia)`, `RF_oil`, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = mbal_results, aes(`P (psia)`, `RF_oil`, color = "Field")) +
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("Condensate Recovery Plot") +
theme_bw()
p3
p4 <- forecast_results %>% ggplot(aes(`P (psia)`, `Liq_volume`, color = "Forecast")) +
geom_point(size = 3) +
scale_color_manual(name="Data", values=c("Forecast" = "red")) +
ggtitle("CCE Liquid Volume Plot") +
theme_bw()
p4
p5 <- forecast_results %>%
tidyr::pivot_longer(cols = c("Igd", "Ifwd"), names_to = "Drive Mechanism",
values_to = "Fraction", values_drop_na = TRUE) %>%
ggplot(aes(`P (psia)`, Fraction, fill = `Drive Mechanism`)) +
geom_area() +
ggtitle("Energy Plot") +
theme_bw()
p5
## ---- fig.width = 6, fig.height= 4, fig.align = "center", warning = FALSE-----
library(Rmbal)
library(Rrelperm)
library(pracma)
library(ggplot2)
library(dplyr)
library(magrittr)
p_pvt <- c(8000, 7500, 7280, 7250, 7000, 6500, 6000, 5500, 5000, 4500, 4000, 3500,
3000, 2500, 2000, 1500, 1000) # psia
Bo <- c(12.732, 13.044, 13.192, 1.054, 1.041, 1.018, 1.002, 0.983, 0.965, 0.947,
0.93, 0.913, 0.896, 0.877, 0.858, 0.839, 0.819) # RB/STB
Rs <- c(22527, 22527, 22527, 860, 819, 754, 704, 648, 593, 541, 490, 440,
389, 336, 280, 228, 169) #SCF/STB
Bg <- c(0.565, 0.579, 0.586, 0.587, 0.595, 0.613, 0.634, 0.661, 0.694, 0.737,
0.795, 0.817, 0.997, 1.178, 1.466, 1.963, 2.912) / 1000 # RB/SCF
Rv <- c(44.4, 44.4, 44.4, 44.3, 43.9, 40.3, 36.5, 32.9, 29.2, 25.4, 21.4, 17.6,
13.7, 10.5, 7.9, 5.8, 4.4) / 1e6 # STB/SCF
cw <- 2e-6
Bwi <- 1.0
Bw <- Bwi * exp(cw * (p_pvt[1] - p_pvt))
muo <- c(0.049, 0.047, 0.046, 19.541, 20.965, 23.958, 26.338, 29.633, 33.319,
37.401, 42.161, 47.465, 53.765, 61.887, 72.143, 83.478, 99.049) # cp
muw <- rep(0.25, length(p_pvt))
mug <- c(0.049, 0.047, 0.046, 0.046, 0.0449, 0.042, 0.0393, 0.0366, 0.0339,
0.0312, 0.0283, 0.0254, 0.0225, 0.0198, 0.0174, 0.0154, 0.014) # cp
pvt_table <- data.frame(p = p_pvt, Bo = Bo, Rs = Rs, Rv = Rv, Bg = Bg,
Bw = Bw, muo = muo, mug = mug, muw = muw)
p <- c(8000, 7500, 7280, 7250, 7000, 6500, 6000, 5500, 5000, 4500, 4000, 3500,
3000, 2500, 2000, 1500, 1000) # psia
Gi <- rep(0, length.out = length(p))
wf <- rep(1, length.out = length(p))
# generating a set of pseudo relative permeability curves using
# laboratory 'Kr' values
sg_lab <- c(0.05, 0.152, 0.248, 0.352, 0.448, 0.552, 0.65)
krg_lab <- c(0, 0.05, 0.09, 0.18, 0.3, 0.5, 1)
kro_lab <- c(1, 0.6, 0.35, 0.13, 0.04, 0.01, 0)
swcrit_lab <- 0.2
sgcrit_lab <- 0.05
sorgr_lab <- 0.15
fun_kr <- function(x, swcrit, sgcrit, sorg, sg, krg, kro) {
kr_est <- Rrelperm::kr2p_gl(SWCON = swcrit, SOIRG = sorg, SORG = sorg,
SGCON = sgcrit, SGCRIT = sgcrit, KRGCL = 1,
KROGCG = 1, NG = x[1], NOG = x[2], NP = 101)
krg_est_sub <- approx(x = kr_est[,1], y = kr_est[,3], xout = sg, rule = 2)$y
kro_est_sub <- approx(x = kr_est[,1], y = kr_est[,4], xout = sg, rule = 2)$y
error <- (krg - krg_est_sub) ^ 2 + (kro - kro_est_sub) ^ 2
return(error)
}
par <- c(2, 2)
opt_results <- minpack.lm::nls.lm(par = par, fn = fun_kr, swcrit = swcrit_lab,
sgcrit = sgcrit_lab, sorg = sorgr_lab, sg = sg_lab,
krg = krg_lab, kro = kro_lab,
lower = c(0.1,0.1), upper = c(10,10))
sol <- opt_results$par
sol
rel_perm <- as.data.frame(Rrelperm::kr2p_gl(SWCON = swcrit_lab, SOIRG = sorgr_lab,
SORG = sorgr_lab, SGCON = sgcrit_lab,
SGCRIT = sgcrit_lab, KRGCL = 1,
KROGCG = 1, NG = sol[1], NOG = sol[2],
NP = 101))
colnames(rel_perm) <- c("Sg", "Sl", "Krg", "Krog")
p_forecast <- p
Gi_forecast <- Gi
wf_forecast <- wf
time_lst_forecast <- mbal_time(c(1:length(p_forecast)), "year")
forecast_lst <- mbal_forecast_param_gas(input_unit = "Field", output_unit = "Field",
G = 1e6, phi = 0.1, swi = 0.2,
pd = 7255, p = p_forecast, pvt = pvt_table,
M = 0, cf = 2e-6, wf = wf_forecast,
rel_perm = rel_perm)
forecast_results <- mbal_forecast_gas(forecast_lst, time_lst_forecast)
forecast_results
reservoir_performance_table <- data.frame(p = p)
reservoir_performance_table$`RF_oil` <- c(0, 2.4, 3.5, 4.2, 5.8, 7.6, 10.4,
13.8, 16.6, 20, 23.1, 25.9, 28.5, 30.8,
32.7, 34.1, 35.2) / 100
reservoir_performance_table$`Sg` <- c(100, 100, 100, 99.99, 99.91, 99.29, 98.68,
98.15, 97.66, 97.21, 96.79, 96.45, 96.15,
95.99, 95.9, 95.9, 95.95) * 0.8 / 100
reservoir_performance_table$`GOR` <- c(22.53, 22.53, 22.53, 22.55, 22.78, 24.81,
27.4, 30.4, 34.25, 39.37, 46.73, 56.82, 72.99,
95.24, 126.58, 172.41, 227.27) * 1000 # SCF/STB
p1 <- forecast_results %>% ggplot(aes(`P (psia)`, SGg, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = reservoir_performance_table, aes(`p`, Sg, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("Gas Saturation Plot") +
theme_bw()
p1
p2 <- forecast_results %>% ggplot(aes(`P (psia)`, `GOR (SCF/STB)`, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = reservoir_performance_table, aes(`p`, GOR, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("GOR Plot") +
theme_bw()
p2
p3 <- forecast_results %>% ggplot(aes(`P (psia)`, `RF_oil`, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = reservoir_performance_table, aes(`p`, RF_oil, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("Oil Recovery Plot") +
theme_bw()
p3
p4 <- forecast_results %>% ggplot(aes(`P (psia)`, `Liq_volume`, color = "Forecast")) +
geom_point(size = 3) +
scale_color_manual(name="Data", values=c("Forecast" = "red")) +
ggtitle("CCE Liquid Volume Plot") +
theme_bw()
p4
p5 <- forecast_results %>%
tidyr::pivot_longer(cols = c("Igd", "Ifwd"), names_to = "Drive Mechanism",
values_to = "Fraction", values_drop_na = TRUE) %>%
ggplot(aes(`P (psia)`, Fraction, fill = `Drive Mechanism`)) +
geom_area() +
ggtitle("Energy Plot") +
theme_bw()
p5
## ---- fig.width = 6, fig.height= 4, fig.align = "center", warning = FALSE-----
library(Rmbal)
library(Rrelperm)
library(pracma)
library(ggplot2)
library(dplyr)
library(magrittr)
p_pvt <- c(5800, 5550, 5450, 5420, 5300, 4800, 4300, 3800, 3300, 2800, 2300,
1800, 1300, 800) # psia
Bo <- c(4.382, 4.441, 4.468, 2.378, 2.366, 2.032, 1.828, 1.674, 1.554,
1.448, 1.36, 1.279, 1.2, 1.131) # RB/STB
Rs <- c(6042, 6042, 6042, 2795, 2750, 2128, 1730, 1422, 1177, 960, 776, 607,
443, 293) #SCF/STB
Bg <- c(0.725, 0.735, 0.739, 0.74, 0.743, 0.758, 0.794, 0.854, 0.947,
1.09, 1.313, 1.677, 2.316, 3.695) / 1000 # RB/SCF
Rv <- c(165.5, 165.5, 165.5, 164.2, 156.6, 114, 89, 65.2, 48.3, 35, 25,
19, 15, 13.5) / 1e6 # STB/SCF
cw <- 2e-6
Bwi <- 1.0
Bw <- Bwi * exp(cw * (p_pvt[1] - p_pvt))
muo <- c(0.0612, 0.062, 0.0587, 0.135, 0.1338, 0.1826, 0.2354, 0.3001, 0.3764,
0.4781, 0.6041, 0.7746, 1.0295, 1.358) # cp
muw <- rep(0.25, length(p_pvt))
mug <- c(0.0612, 0.062, 0.0587, 0.0581, 0.0554, 0.0436, 0.0368, 0.0308, 0.0261,
0.0222, 0.0191, 0.0166, 0.0148, 0.0135) # cp
pvt_table <- data.frame(p = p_pvt, Bo = Bo, Rs = Rs, Rv = Rv, Bg = Bg,
Bw = Bw, muo = muo, mug = mug, muw = muw)
p <- c(5800, 5550, 5450, 5420, 5300, 4800, 4300, 3800, 3300, 2800, 2300,
1800, 1300, 800) # psia
Gi <- rep(0, length.out = length(p))
wf <- rep(1, length.out = length(p))
# generating a set of pseudo relative permeability curves using
# laboratory 'Kr' values
sg_lab <- c(0.05, 0.152, 0.248, 0.352, 0.448, 0.552, 0.65)
krg_lab <- c(0, 0.05, 0.09, 0.18, 0.3, 0.5, 1)
kro_lab <- c(1, 0.6, 0.35, 0.13, 0.04, 0.01, 0)
swcrit_lab <- 0.2
sgcrit_lab <- 0.05
sorgr_lab <- 0.15
fun_kr <- function(x, swcrit, sgcrit, sorg, sg, krg, kro) {
kr_est <- Rrelperm::kr2p_gl(SWCON = swcrit, SOIRG = sorg, SORG = sorg,
SGCON = sgcrit, SGCRIT = sgcrit, KRGCL = 1,
KROGCG = 1, NG = x[1], NOG = x[2], NP = 101)
krg_est_sub <- approx(x = kr_est[,1], y = kr_est[,3], xout = sg, rule = 2)$y
kro_est_sub <- approx(x = kr_est[,1], y = kr_est[,4], xout = sg, rule = 2)$y
error <- (krg - krg_est_sub) ^ 2 + (kro - kro_est_sub) ^ 2
return(error)
}
par <- c(2, 2)
opt_results <- minpack.lm::nls.lm(par = par, fn = fun_kr, swcrit = swcrit_lab,
sgcrit = sgcrit_lab, sorg = sorgr_lab, sg = sg_lab,
krg = krg_lab, kro = kro_lab,
lower = c(0.1,0.1), upper = c(10,10))
sol <- opt_results$par
sol
rel_perm <- as.data.frame(Rrelperm::kr2p_gl(SWCON = swcrit_lab, SOIRG = sorgr_lab,
SORG = sorgr_lab, SGCON = sgcrit_lab,
SGCRIT = sgcrit_lab, KRGCL = 1,
KROGCG = 1, NG = sol[1], NOG = sol[2],
NP = 101))
colnames(rel_perm) <- c("Sg", "Sl", "Krg", "Krog")
p_forecast <- p
Gi_forecast <- Gi
wf_forecast <- wf
time_lst_forecast <- mbal_time(c(1:length(p_forecast)), "year")
forecast_lst <- mbal_forecast_param_gas(input_unit = "Field", output_unit = "Field",
G = 1e6, phi = 0.1, swi = 0.2,
pd = 5430, p = p_forecast, pvt = pvt_table,
M = 0, cf = 2e-6, wf = wf_forecast,
rel_perm = rel_perm)
forecast_results <- mbal_forecast_gas(forecast_lst, time_lst_forecast)
forecast_results
reservoir_performance_table <- data.frame(p = p)
reservoir_performance_table$`RF_oil` <- c(0, 1.3, 1.9, 2.1, 2.6, 7, 10.1, 13.3,
16.2, 18.4, 20.2, 21.6, 22.8, 23.7) / 100
reservoir_performance_table$`Sg` <- c(100, 100, 100, 99.2, 95, 81.8, 78.7, 76.7, 76.3,
76.6, 77.2, 78.3, 79.5, 80.6) * 0.8 / 100
reservoir_performance_table$`GOR` <- c(6.04, 6.04, 6.04, 6.09, 6.39, 8.77, 11.22,
15.29, 20.62, 28.45, 39.82, 52.42, 66.48,
73.99) * 1000 # SCF/STB
p1 <- forecast_results %>% ggplot(aes(`P (psia)`, SGg, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = reservoir_performance_table, aes(`p`, Sg, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("Gas Saturation Plot") +
theme_bw()
p1
p2 <- forecast_results %>% ggplot(aes(`P (psia)`, `GOR (SCF/STB)`, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = reservoir_performance_table, aes(`p`, GOR, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("GOR Plot") +
theme_bw()
p2
p3 <- forecast_results %>% ggplot(aes(`P (psia)`, `RF_oil`, color = "Forecast")) +
geom_point(size = 3) +
geom_point(data = reservoir_performance_table, aes(`p`, RF_oil, color = "Field"))+
scale_color_manual(name="Data", values=c("Forecast" = "red", "Field" = "black")) +
ggtitle("Oil Recovery Plot") +
theme_bw()
p3
p4 <- forecast_results %>% ggplot(aes(`P (psia)`, `Liq_volume`, color = "Forecast")) +
geom_point(size = 3) +
scale_color_manual(name="Data", values=c("Forecast" = "red")) +
ggtitle("CCE Liquid Volume Plot") +
theme_bw()
p4
p5 <- forecast_results %>%
tidyr::pivot_longer(cols = c("Igd", "Ifwd"), names_to = "Drive Mechanism",
values_to = "Fraction", values_drop_na = TRUE) %>%
ggplot(aes(`P (psia)`, Fraction, fill = `Drive Mechanism`)) +
geom_area() +
ggtitle("Energy Plot") +
theme_bw()
p5
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.