Commit 2f2350b8 authored by Marco De Lucia's avatar Marco De Lucia
Browse files

Added RPinfo/related functions and corresponding man pages

parent 987ef545
Package: RedModRphree
Title: Programmable interface to the phreeqc geochemical solver adding features such as surrogate models, reactive transport and Pourbaix diagrams
Version: 0.2
Authors@R: c(person("Marco", "De Lucia", email = "delucia@gfz-potsdam.de", role = c("aut", "cre")),
Authors@R: c(person("Marco", "De Lucia", email = "delucia@gfz-potsdam.de", role = c("aut", "cre"), comment=, comment = c(ORCID = "0000-0003-0918-3766")),
person("Janis", "Jatnieks", email = "deltaxzz@gmail.com", role = c("ctb")))
Author: Marco De Lucia [aut, cre],
Janis Jatnieks [ctb]
......
......@@ -25,10 +25,14 @@ export(Matplot)
export(MatplotSingle)
export(ParseFormula)
export(PlotModsInSample)
export(RGetPhases)
export(RPhreeExt)
export(RPhreeFile)
export(RPhreeWriteInp)
export(RPinfo)
export(RReadOut)
export(RReadOutKin)
export(RReadPhases)
export(ReactTranspBalanceEq)
export(ReactTranspBalanceKin)
export(RecomposeState)
......
## Functions for dealing with surrogate simulations
### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2018
### Time-stamp: "Last modified 2021-04-12 12:37:10 delucia"
### Time-stamp: "Last modified 2021-04-12 13:05:02 delucia"
##' Todo
......@@ -9,12 +9,13 @@
##' @title Simple explicit 1D forward-Euler advection, multispecies
##' @param conc matrix containing the chemical state
##' @param inflow vector of concentrations entering the inlet
##' @param dx the grid spacing
##' @param dt the required time step
##' @param dx the grid spacing in m
##' @param dt the required time step in s
##' @param U Darcy velocity in m^3/s
##' @return a matrix containing the state after the advection
##' @author MDL
##' @export
AdvectionPQC <- function(conc, inflow=rep(0,ncol(conc)), dx, dt)
AdvectionPQC <- function(conc, inflow=rep(0,ncol(conc)), dx, dt, U)
{
if (is.matrix(conc)) {
## msg("conc is matrix")
......@@ -333,7 +334,7 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix
## transform the pH/pe back into activities for Advection
msg("first step, Advection()...")
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt)
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U)
msg("Advection done...")
Tr <- Ts <- 0
......@@ -401,7 +402,7 @@ ReactTranspBalanceEq <- function(setup, init, maxtime, step=c("time","iter","fix
#cat("ENTERING ADV","\n")
## transform the pH/pe back into activities
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt)
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U)
## treat special H+/pH, e-/pe cases
state_T <- Act2pH(state_T)
......@@ -636,7 +637,7 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi
## transform the pH/pe back into activities for Advection
msg("first step, Advection()...")
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt)
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U)
Tr <- Ts <- 0
## treat special H+/pH, e-/pe cases
......@@ -704,7 +705,7 @@ ReactTranspBalanceKin <- function(setup, init, maxtime, step=c("time","iter","fi
#cat("ENTERING ADV","\n")
## transform the pH/pe back into activities
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt)
state_T <- AdvectionPQC(pH2Act(state_C), inflow=bound, dx=dx, dt=dt, U=U)
## reduction of the problem
if(reduce)
......
### Utility functions for RedModRphree
### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2018
### Time-stamp: "Last modified 2018-05-03 11:03:54 delucia"
### Time-stamp: "Last modified 2021-04-12 12:58:36 delucia"
##' Replicates an input buffer containing only one SOLUTION, taking
##' care of SOLUTION/KINETICS/PURE identifiers. Eventually insert a
......@@ -216,3 +216,145 @@ SuppressSim <- function(biginp, n=1L)
newinp <- biginp[ -seq(startsols[n], startsols[n+1]-1) ]
return(newinp)
}
##' This functions extracts informations from calculated Rphree
##' formatted solutions.
##'
##' @title Extract specific values from an Rphree output
##' @param lin The list containg the Rphree solution(s)
##' @param cat A string indicating the category (or output block) to
##' search in, must be one of: "tot, desc, pphases, master, species,
##' kin, SI"
##' @param prop The name of the inquired property (element, species,
##' mineral phase,\dots) whose value(s)
##' @param force Logical. If TRUE (default if left unspecified), a valid numeric value (0) is returned even if the property is not found
##' value is returned if the inquired property is not found in the
##' solution
##' @param flex Logical. If TRUE, expects no "ListInfo" in the
##' formatted solution list and performs heuristics to circumvent this
##' absence
##' @return A numeric vector containing the inquired properties
##' @author MDL
##' @examples
##' \dontrun{
##' RPinfo(SOL,"tot","Fe")
##' RPinfo(SOL,"pphases","Anhydrite")
##' }
##' @export
RPinfo <- function(lin, cat=c("tot","desc","pphases","master","species","kin","SI"),
prop, force=TRUE, flex=FALSE)
{
cat <- match.arg(cat)
if ("ListInfo" %in% names(lin))
nsim <- lin$ListInfo$n
else {
if (!flex)
stop("RPinfo:: no ListInfo! Perhaps not a valid list produced by Rphree?
Specify flex=TRUE if you still want to try")
else
nsim <- length(lin)
}
if ( nsim > 1 )
## we have a list of solutions, check if the name exists
{
if ( cat %in% names(lin[[1]])) {
## everything seems fine, let's (s)apply
lout <- sapply(lin[1:nsim], RPhreeExt, cat=cat, prop=prop)
} else {
if (force)
{
lout <- rep(0,length(lin))
} else {
stop("RPinfo:: Sure that the right output was selected in the simulation?")
}
}
} else {
## just 1 solution
if ( cat %in% names(lin)) {
lout <- RPhreeExt(lin, cat, prop)
} else {
if (force)
{
lout <- 0
} else {
stop("RPinfo:: Sure that the right output was selected in the simulation?")
}
}
}
return(lout)
}
##' Workhorse function for extraction of results from a solution,
##' primarily intended to be used by \code{\link{RPinfo}}, where it is
##' called inside a \code{sapply} statement.
##'
##' @title Workhorse function for extraction of informations from a
##' solution
##' @param sol The solution outputted by Rphree
##' @param cat String: the category or output block where to look for
##' the property \code{prop}
##' @param prop The name of the property to look for
##' @return Numeric of length 1, the inquired value for the property
##' @author MDL
##' @export
RPhreeExt <- function(sol, cat, prop)
{
if (is.data.frame(sol[[cat]]))
return(sol[[cat]][prop,1])
else
return(0)
}
##' Workhorse function for extraction of results from a solution,
##' specifically equilibrium phases.
##'
##' @title Workhorse function for extraction of
##' @param lin The solution outputted by Rphree
##' @param which String containing the name of the specific pphase to
##' look for. If omitted, all pphases are returned.
##' @param delta Logical. If TRUE, delta min are returned, if FALSE,
##' total amount
##' @return A data.frame containing the results of equilibrium phases
##' @author MDL
##' @export
RGetPhases <- function(lin, which=NULL, delta=TRUE)
{
i <- 1
if (delta) i <- 2
tab <- lin$pphases
if (!is.character(which)) {
return(tab[,i])
} else {
inde <- match(which,rownames(tab))
return(tab[inde,i])
}
}
##' Workhorse function for extraction of results from a solution,
##' specifically equilibrium phases. Variant.
##'
##' @title extract the phases from a Rphree solution
##' @param lin The solution outputted by Rphree
##' @return A data.frame containing the results of equilibrium phases
##' @author MDL
##' @export
RReadPhases <- function(lin)
{
ns <- length( phas <- unique(unlist(lapply(lin,function (x) return(rownames(x$pphases))))) )
if (ns==0)
return(NULL)
mat <- matrix(0,ncol=ns,nrow=n)
for (i in seq_along(phas))
{
mat[,i] <- RPinfo(lin,"pphases",phas[i])
}
colnames(mat) <- phas
mat[is.na(mat)] <- 0
return(mat)
}
### RedModRphree, demo of kinetics
### Licence: LGPL version 2.1
### Time-stamp: "Last modified 2018-05-06 21:26:21 delucia"
### Time-stamp: "Last modified 2021-04-12 17:59:02 delucia"
library(RedModRphree)
db <- RPhreeFile(system.file("extdata", "phreeqc_kin.dat", package="RedModRphree"), is.db=TRUE)
......
......@@ -4,16 +4,18 @@
\alias{AdvectionPQC}
\title{Simple explicit 1D forward-Euler advection, multispecies}
\usage{
AdvectionPQC(conc, inflow = rep(0, ncol(conc)), dx, dt)
AdvectionPQC(conc, inflow = rep(0, ncol(conc)), dx, dt, U)
}
\arguments{
\item{conc}{matrix containing the chemical state}
\item{inflow}{vector of concentrations entering the inlet}
\item{dx}{the grid spacing}
\item{dx}{the grid spacing in m}
\item{dt}{the required time step}
\item{dt}{the required time step in s}
\item{U}{Darcy velocity in m^3/s}
}
\value{
a matrix containing the state after the advection
......
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/Rphree_Utils.R
\name{RGetPhases}
\alias{RGetPhases}
\title{Workhorse function for extraction of}
\usage{
RGetPhases(lin, which = NULL, delta = TRUE)
}
\arguments{
\item{lin}{The solution outputted by Rphree}
\item{which}{String containing the name of the specific pphase to
look for. If omitted, all pphases are returned.}
\item{delta}{Logical. If TRUE, delta min are returned, if FALSE,
total amount}
}
\value{
A data.frame containing the results of equilibrium phases
}
\description{
Workhorse function for extraction of results from a solution,
specifically equilibrium phases.
}
\author{
MDL
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/Rphree_Utils.R
\name{RPhreeExt}
\alias{RPhreeExt}
\title{Workhorse function for extraction of informations from a
solution}
\usage{
RPhreeExt(sol, cat, prop)
}
\arguments{
\item{sol}{The solution outputted by Rphree}
\item{cat}{String: the category or output block where to look for
the property \code{prop}}
\item{prop}{The name of the property to look for}
}
\value{
Numeric of length 1, the inquired value for the property
}
\description{
Workhorse function for extraction of results from a solution,
primarily intended to be used by \code{\link{RPinfo}}, where it is
called inside a \code{sapply} statement.
}
\author{
MDL
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/Rphree_Utils.R
\name{RPinfo}
\alias{RPinfo}
\title{Extract specific values from an Rphree output}
\usage{
RPinfo(
lin,
cat = c("tot", "desc", "pphases", "master", "species", "kin", "SI"),
prop,
force = TRUE,
flex = FALSE
)
}
\arguments{
\item{lin}{The list containg the Rphree solution(s)}
\item{cat}{A string indicating the category (or output block) to
search in, must be one of: "tot, desc, pphases, master, species,
kin, SI"}
\item{prop}{The name of the inquired property (element, species,
mineral phase,\dots) whose value(s)}
\item{force}{Logical. If TRUE (default if left unspecified), a valid numeric value (0) is returned even if the property is not found
value is returned if the inquired property is not found in the
solution}
\item{flex}{Logical. If TRUE, expects no "ListInfo" in the
formatted solution list and performs heuristics to circumvent this
absence}
}
\value{
A numeric vector containing the inquired properties
}
\description{
This functions extracts informations from calculated Rphree
formatted solutions.
}
\examples{
\dontrun{
RPinfo(SOL,"tot","Fe")
RPinfo(SOL,"pphases","Anhydrite")
}
}
\author{
MDL
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/Rphree_Utils.R
\name{RReadPhases}
\alias{RReadPhases}
\title{extract the phases from a Rphree solution}
\usage{
RReadPhases(lin)
}
\arguments{
\item{lin}{The solution outputted by Rphree}
}
\value{
A data.frame containing the results of equilibrium phases
}
\description{
Workhorse function for extraction of results from a solution,
specifically equilibrium phases. Variant.
}
\author{
MDL
}
Markdown is supported
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