max/poet
1
0
Fork 0
mirror of https://git.gfz-potsdam.de/naaice/poet.git synced 2025-12-15 20:38:23 +01:00
Code Issues Packages Projects Releases Wiki Activity

fix: update grid dimensions and refine species definitions in barite_fgcs_2.R

Browse Source
This commit is contained in:
Max Lübke 2024-12-20 10:05:37 +01:00
parent b93cdeed83
commit cefea926cf
1 changed files with 19 additions and 19 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
bench/fgcs/barite_fgcs_2.R
View File

@ -1,28 +1,28 @@
## Time-stamp: "Last modified 2024-12-11 16:08:11 delucia" ## Time-stamp: "Last modified 2024-12-11 16:08:11 delucia"
cols <- 200 cols <- 1000
rows <- 200 rows <- 1000
dim_cols <- 10 dim_cols <- 50
dim_rows <- 10 dim_rows <- 50
ncirc <- 20 ## number of crystals ncirc <- 20 ## number of crystals
rmax <- cols/10 ## max radius (in nodes) rmax <- cols / 10 ## max radius (in nodes)
set.seed(22933) set.seed(22933)
centers <- cbind(sample(seq_len(cols), ncirc), sample(seq_len(rows), ncirc)) centers <- cbind(sample(seq_len(cols), ncirc), sample(seq_len(rows), ncirc))
radii <- sample(seq_len(rmax), ncirc, replace=TRUE) radii <- sample(seq_len(rmax), ncirc, replace = TRUE)
mi <- matrix(rep(seq_len(cols), rows), byrow = TRUE, nrow = rows) mi <- matrix(rep(seq_len(cols), rows), byrow = TRUE, nrow = rows)
mj <- matrix(rep(seq_len(cols), each = rows), byrow = TRUE, nrow = rows) mj <- matrix(rep(seq_len(cols), each = rows), byrow = TRUE, nrow = rows)
tmpl <- lapply(seq_len(ncirc), function(x) which((mi-centers[x, 1])^2 + (mj-centers[x, 2])^2 < radii[x]^2, arr.ind = TRUE)) tmpl <- lapply(seq_len(ncirc), function(x) which((mi - centers[x, 1])^2 + (mj - centers[x, 2])^2 < radii[x]^2, arr.ind = TRUE))
inds <- do.call(rbind, tmpl) inds <- do.call(rbind, tmpl)
grid <- matrix(1, nrow = rows, ncol = cols) grid <- matrix(1, nrow = rows, ncol = cols)
grid[inds] <- 2 grid[inds] <- 2
alpha <- matrix( 1e-5, ncol = cols, nrow = rows) alpha <- matrix(1e-5, ncol = cols, nrow = rows)
alpha[inds] <- 1e-7 alpha[inds] <- 1e-7
## image(grid, asp=1) ## image(grid, asp=1)
@ -30,13 +30,13 @@ alpha[inds] <- 1e-7
## Define grid configuration for POET model ## Define grid configuration for POET model
grid_setup <- list( grid_setup <- list(
pqc_in_file = "./barite_fgcs_2.pqi", pqc_in_file = "./barite_fgcs_2.pqi",
pqc_db_file = "./db_barite.dat", ## database file pqc_db_file = "../barite/db_barite.dat", ## database file
grid_def = grid, ## grid definition, IDs according to the Phreeqc input grid_def = grid, ## grid definition, IDs according to the Phreeqc input
grid_size = c(dim_cols, dim_rows), ## grid size in meters grid_size = c(dim_cols, dim_rows), ## grid size in meters
constant_cells = c() ## IDs of cells with constant concentration constant_cells = c() ## IDs of cells with constant concentration
) )
bound_length <- cols/10 bound_length <- cols / 10
bound_N <- list( bound_N <- list(
"type" = rep("constant", bound_length), "type" = rep("constant", bound_length),
@ -52,13 +52,13 @@ bound_W <- list(
bound_E <- list( bound_E <- list(
"type" = rep("constant", bound_length), "type" = rep("constant", bound_length),
"sol_id" = rep(4, bound_length), "sol_id" = rep(4, bound_length),
"cell" = seq(rows-bound_length+1, rows) "cell" = seq(rows - bound_length + 1, rows)
) )
bound_S <- list( bound_S <- list(
"type" = rep("constant", bound_length), "type" = rep("constant", bound_length),
"sol_id" = rep(4, bound_length), "sol_id" = rep(4, bound_length),
"cell" = seq(cols-bound_length+1, cols) "cell" = seq(cols - bound_length + 1, cols)
) )
diffusion_setup <- list( diffusion_setup <- list(
@ -77,7 +77,7 @@ dht_species <- c(
"O" = 7, "O" = 7,
"Ba" = 7, "Ba" = 7,
"Cl" = 7, "Cl" = 7,
"S(6)" = 7, "S" = 7,
"Sr" = 7, "Sr" = 7,
"Barite" = 4, "Barite" = 4,
"Celestite" = 4 "Celestite" = 4
@ -86,10 +86,10 @@ dht_species <- c(
pht_species <- c( pht_species <- c(
"Ba" = 4, "Ba" = 4,
"Cl" = 3, "Cl" = 3,
"S(6)" = 3, "S" = 3,
"Sr" = 3, "Sr" = 3,
"Barite" = 2, "Barite" = 0,
"Celestite" = 2 "Celestite" = 0
) )
chemistry_setup <- list( chemistry_setup <- list(