Commit 88a879aa authored by Marco De Lucia's avatar Marco De Lucia
Browse files

Fixed demos with imports, out_bal & co in ReactTrans now assigned at the end to GlobalEnv

parent 4a6597de
## Functions for dealing with simulations with kinetics ## Functions for dealing with simulations with kinetics
### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2021 ### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2021
### Time-stamp: "Last modified 2021-04-28 16:06:07 delucia" ### Time-stamp: "Last modified 2021-04-28 17:23:36 delucia"
##' This function runs the generated input buffer (or a list thereof) ##' This function runs the generated input buffer (or a list thereof)
##' through \code{phreeqc}, which has been already loaded as ##' through \code{phreeqc}, which has been already loaded as
##' dependency. Also the database must be previously and explicitely ##' dependency. Also the database must be previously and explicitely
##' loaded by the user using, e.g., \code{db <- ##' loaded by the user using, e.g.,
##' RPhreeFile('/path/to/database.dat', is.db=TRUE)} followed by ##' \code{db <- RPhreeFile('/path/to/database.dat', is.db=TRUE)}
##' \code{phreeqc::phrLoadDatabaseString(db)} or using one of the ##' followed by
##' databases included in \code{phreeqc}, e.g., ##' \code{phreeqc::phrLoadDatabaseString(db)}
##' \code{phreeqc::phrLoadDatabaseString(phreeqc.dat)}. ##' or using one of the databases included in \code{phreeqc}, e.g.,
##' ##' \code{phreeqc::phrLoadDatabaseString(phreeqc.dat)}
##' Parallelization is achieved using \code{%dopar%} from package ##' Parallelization is achieved using \code{\%dopar\%} from package
##' \code{foreach}, while the backend should have been previously set ##' \code{foreach}, while the backend should have been previously set
##' with appropriate call to \code{doParallel}. At the moment the ##' with appropriate call to \code{doParallel}. At the moment the
##' function doesn't check that the SELCTED_OUTPUT block is ##' function doesn't check that the SELCTED_OUTPUT block is
...@@ -51,8 +51,6 @@ RunPQC <- function(input, procs=1, second=TRUE) { ...@@ -51,8 +51,6 @@ RunPQC <- function(input, procs=1, second=TRUE) {
} }
} else { } else {
## go parallel! ## go parallel!
## old one
## res <- parallel::mclapply(input, .runPQC, mc.silent=TRUE, mc.cores=procs) ## res <- parallel::mclapply(input, .runPQC, mc.silent=TRUE, mc.cores=procs)
res <- foreach::foreach(i=seq_along(input), .combine=rbind) %dopar% .runPQC(input[[i]], onlysecond=second) res <- foreach::foreach(i=seq_along(input), .combine=rbind) %dopar% .runPQC(input[[i]], onlysecond=second)
......
### Licence: LGPL version 2.1 ### Licence: LGPL version 2.1
## Time-stamp: "Last modified 2021-04-28 17:04:27 delucia" ## Time-stamp: "Last modified 2021-04-28 17:20:11 delucia"
##' @title Revert a string ##' @title Revert a string
##' @param x character vector. All element of the vector will be ##' @param x character vector. All element of the vector will be
...@@ -439,7 +439,7 @@ Pourbaix <- function(element, pe=seq(-5,8), pH=seq(6,12), Tc=25, Patm=1, base, f ...@@ -439,7 +439,7 @@ Pourbaix <- function(element, pe=seq(-5,8), pH=seq(6,12), Tc=25, Patm=1, base, f
if (missing(ann.title)) { if (missing(ann.title)) {
mymain <- paste("Pourbaix diagram for", mymain <- paste("Pourbaix diagram for",
ifelse(missing(element), "the specified system", ifelse(has_valence, valence, element)), ifelse(missing(element), "the specified system", ifelse(has_valence, valence, element)),
"\n T=", Tc, "\U3f43 C / P=", Patm, "atm") "\n T=", Tc, "\u00B0C / P=", Patm, "atm")
} else } else
mymain <- ann.title mymain <- ann.title
......
## Functions for dealing with surrogate simulations ## Functions for dealing with surrogate simulations
### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2021 ### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2021
### Time-stamp: "Last modified 2021-04-28 15:54:48 delucia" ### Time-stamp: "Last modified 2021-04-28 18:13:39 delucia"
##' This function takes the current state of a chemical system in form ##' This function takes the current state of a chemical system in form
...@@ -294,8 +294,6 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix ...@@ -294,8 +294,6 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix
spr_string <- paste0("%0", pad,"d") spr_string <- paste0("%0", pad,"d")
names(saved_complete) <- paste0("step_",sprintf(spr_string,seq(0,maxiter))) names(saved_complete) <- paste0("step_",sprintf(spr_string,seq(0,maxiter)))
utils::globalVariables(c("out_inp","out_res","out_bal"))
out_inp <- vector(mode="list",length=length(maxiter)) out_inp <- vector(mode="list",length=length(maxiter))
out_res <- vector(mode="list",length=length(maxiter)) out_res <- vector(mode="list",length=length(maxiter))
out_bal <- vector(mode="list",length=length(maxiter)) out_bal <- vector(mode="list",length=length(maxiter))
...@@ -406,7 +404,6 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix ...@@ -406,7 +404,6 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix
# so that the first run on Rphree is actually one and the next ones end # so that the first run on Rphree is actually one and the next ones end
time <- time + dt time <- time + dt
iter <- iter + 1 iter <- iter + 1
#cat("ENTERING ADV","\n")
## transform the pH/pe back into activities ## transform the pH/pe back into activities
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U) state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U)
...@@ -485,12 +482,16 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix ...@@ -485,12 +482,16 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix
} }
on.exit() on.exit()
if (procs > 1) { if (procs > 1) {
parallel::stopCluster(ThisRunCluster) parallel::stopCluster(ThisRunCluster)
msg("Parallelization cluster stopped") msg("Parallelization cluster stopped")
} }
msg("assign out_balance, out_input and out_results to global environment for possible later use")
assign("out_balance", out_bal, envir = .GlobalEnv)
assign("out_results", out_res, envir = .GlobalEnv)
assign("out_input", out_inp, envir = .GlobalEnv)
msg("total time of chemistry: ",round(sum(timing[,3]),3),"seconds") msg("total time of chemistry: ",round(sum(timing[,3]),3),"seconds")
msg("total number of simulations: ",sum(timing[,1])) msg("total number of simulations: ",sum(timing[,1]))
...@@ -604,8 +605,6 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi ...@@ -604,8 +605,6 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi
spr_string <- paste0("%0", pad,"d") spr_string <- paste0("%0", pad,"d")
names(saved_complete) <- paste0("step_",sprintf(spr_string,seq(0,maxiter))) names(saved_complete) <- paste0("step_",sprintf(spr_string,seq(0,maxiter)))
utils::globalVariables(c("out_inp","out_res","out_bal"))
out_inp <- vector(mode="list",length=length(maxiter)) out_inp <- vector(mode="list",length=length(maxiter))
out_res <- vector(mode="list",length=length(maxiter)) out_res <- vector(mode="list",length=length(maxiter))
out_bal <- vector(mode="list",length=length(maxiter)) out_bal <- vector(mode="list",length=length(maxiter))
...@@ -711,7 +710,6 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi ...@@ -711,7 +710,6 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi
# so that the first run on Rphree is actually one and the next ones end # so that the first run on Rphree is actually one and the next ones end
time <- time + dt time <- time + dt
iter <- iter + 1 iter <- iter + 1
#cat("ENTERING ADV","\n")
## transform the pH/pe back into activities ## transform the pH/pe back into activities
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U) state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U)
...@@ -737,10 +735,8 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi ...@@ -737,10 +735,8 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi
out_bal[[iter]] <- bal out_bal[[iter]] <- bal
totbal <- apply(bal, 1, mae) totbal <- apply(bal, 1, mae)
## cat(":: Surrogate Balance Summary:\n")
## print(summary(totbal)) ## print(summary(totbal))
## tol2 <- unname(quantile(totbal, .85)) ## tol2 <- unname(quantile(totbal, .85))
## cat(" :RT: Balance Checked\n")
nonok <- which(totbal>tol) nonok <- which(totbal>tol)
if (length(nonok)>0){ if (length(nonok)>0){
...@@ -804,6 +800,11 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi ...@@ -804,6 +800,11 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi
msg("Parallelization cluster stopped") msg("Parallelization cluster stopped")
} }
msg("assign out_balance, out_input and out_results to global environment for possible later use")
assign("out_balance", out_bal, envir = .GlobalEnv)
assign("out_results", out_res, envir = .GlobalEnv)
assign("out_input", out_inp, envir = .GlobalEnv)
msg("total time of chemistry: ",round(sum(timing[,3]),3),"seconds") msg("total time of chemistry: ",round(sum(timing[,3]),3),"seconds")
msg("total number of simulations: ",sum(timing[,1])) msg("total number of simulations: ",sum(timing[,1]))
msg(" Bye.") msg(" Bye.")
......
### Licence: LGPL version 2.1 ### Licence: LGPL version 2.1
## Time-stamp: "Last modified 2018-08-29 15:15:50 delucia" ## Time-stamp: "Last modified 2021-04-28 18:23:19 delucia"
require(e1071)
require(import)
db <- RPhreeFile(system.file("extdata", "phreeqc_kin.dat", package="RedModRphree"), is.db=TRUE) db <- RPhreeFile(system.file("extdata", "phreeqc_kin.dat", package="RedModRphree"), is.db=TRUE)
phreeqc::phrLoadDatabaseString(db) phreeqc::phrLoadDatabaseString(db)
......
## Licence: LGPL version 2.1 ## Licence: LGPL version 2.1
## Time-stamp: "Last modified 2018-05-06 21:29:53 delucia" ## Time-stamp: "Last modified 2021-04-28 18:27:32 delucia"
library(RedModRphree) library(RedModRphree)
require(e1071)
require(import)
db <- RPhreeFile(system.file("extdata", "phreeqc_kin.dat", package="RedModRphree"), is.db=TRUE) db <- RPhreeFile(system.file("extdata", "phreeqc_kin.dat", package="RedModRphree"), is.db=TRUE)
base <- c("SOLUTION 1", "units mol/kgw", "temp 25.0", "water 1", base <- c("SOLUTION 1", "units mol/kgw", "temp 25.0", "water 1",
......
...@@ -22,14 +22,17 @@ a data.frame containing the SELECTED_OUTPUT, ordered if ...@@ -22,14 +22,17 @@ a data.frame containing the SELECTED_OUTPUT, ordered if
This function runs the generated input buffer (or a list thereof) This function runs the generated input buffer (or a list thereof)
through \code{phreeqc}, which has been already loaded as through \code{phreeqc}, which has been already loaded as
dependency. Also the database must be previously and explicitely dependency. Also the database must be previously and explicitely
loaded by the user using, e.g., \code{db <- loaded by the user using, e.g.,
RPhreeFile('/path/to/database.dat', is.db=TRUE)} followed by \code{db <- RPhreeFile('/path/to/database.dat', is.db=TRUE)}
\code{phreeqc::phrLoadDatabaseString(db)} or using one of the followed by
databases included in \code{phreeqc}, e.g., \code{phreeqc::phrLoadDatabaseString(db)}
\code{phreeqc::phrLoadDatabaseString(phreeqc.dat)}. or using one of the databases included in \code{phreeqc}, e.g.,
} \code{phreeqc::phrLoadDatabaseString(phreeqc.dat)}
\details{ Parallelization is achieved using \code{\%dopar\%} from package
\code{foreach}, while the backend should have been previously set
with appropriate call to \code{doParallel}. At the moment the
function doesn't check that the SELCTED_OUTPUT block is
well-formed, so this is responsibility of the user.
} }
\author{ \author{
MDL MDL
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment