Commit 7b53474e authored by Marco De Lucia's avatar Marco De Lucia
Browse files

Fixing stuff here and there, too much to track at this point. It works.

parent 7a4acb99
......@@ -4,16 +4,8 @@ Version: 0.0.0.1
Authors@R: c(person("Marco", "De Lucia", email = "delucia@gfz-potsdam.de", role = c("aut", "cre")),
person("Janis", "Jatnieks", email = "deltaxzz@gmail.com", role = c("ctb")))
Description: In the framework of RedMod project, we develop an R infrastructure for training and applying geochemical surrogate models to reactive transport.
Depends: R (>= 3.2.0)
Depends: R (>= 3.2.0), parallel, phreeqc, caret, mgcv, randomForest, graphics, methods
License: LGPL-2.1
Encoding: UTF-8
Imports:
parallel,
methods,
phreeqc,
mgcv,
caret,
randomForest,
graphics
LLazyData: true
RoxygenNote: 6.0.1
......@@ -9,6 +9,7 @@ export(CheckBalance)
export(Distribute)
export(DistributeKin)
export(DistributeKinMatrix)
export(DistributeMatrix)
export(ElementalBalanceMin)
export(ExtractPphases)
export(ExtractSpecies)
......@@ -26,6 +27,7 @@ export(RPhreeFile)
export(RPhreeWriteInp)
export(RReadOut)
export(RReadOutKin)
export(ReactTranspBalanceEq)
export(ReactTranspBalanceKin)
export(RecomposeState)
export(ReduceState)
......
## Functions for dealing with simulations with kinetics
### Marco De Lucia, delucia@gfz-potsdam.de, 2009-2018
### Time-stamp: "Last modified 2018-05-03 15:18:50 delucia"
### Time-stamp: "Last modified 2018-05-04 20:47:46 delucia"
##' This function just runs the generated input buffer - or a list
##' thereof - through \code{phreeqc}. Obviously it requires the
......@@ -33,10 +33,9 @@ RunPQC <- function(input, procs=1) {
nso <- sub('.mol.kgw.','',nso, fixed=TRUE)
nso <- sub('^k_','',nso)
nso <- sub('.g.','(g)',nso, fixed=TRUE)
out <- out[seq(2, nrow(out), by=2),]
colnames(out) <- nso
if ("time" %in% nso)
out <- out[seq(2, nrow(out), by=2),]
return(out)
return(data.matrix(out))
}
## is lin a list of inputs?
......@@ -44,7 +43,7 @@ RunPQC <- function(input, procs=1) {
if(is.character(input)) { ## normal sequential run
res <- .runPQC(input)
} else {
stop(" :RunPQC: something wrong with the input, dying.\n")
stopmsg("something wrong with the input, dying!")
}
} else {
## go parallel!
......@@ -107,7 +106,7 @@ DistributeKin <- function(input, prop=NULL, values, ident="-m0")
ids <- idline + initline -1
if (length(ids)!=n)
stop(" :DistributeKin: Length must match the number of values\n")
stopmsg("Length must match the number of values")
newinp[ids] <- paste(ident,values,sep=" ")
return(newinp)
......@@ -153,7 +152,6 @@ splitMultiKin <- function(data, procs, base, first, prop, minerals, kin, dt, ann
newinp <- base
for (i in indAq) {
## cat(paste("distributing", prop[i],", i=",i, "values",paste(values[,i], collapse=" "),"\n"))
newinp <- Distribute(newinp, prop[i], tmpdata[i])
}
......@@ -169,13 +167,12 @@ splitMultiKin <- function(data, procs, base, first, prop, minerals, kin, dt, ann
}
for (i in kin) {
## cat(paste("distributing Kinetic", prop[i],", i=",i, "values",paste(values[,i], collapse=" "),"\n"))
newinp <- DistributeKin(newinp, prop[i], tmpdata[i], ident="-m")
}
## And finally the time
newinp <- Distribute(newinp, "-steps", dt)
return(newinp)
return(c(first,newinp))
}
......@@ -214,7 +211,56 @@ splitMultiKin <- function(data, procs, base, first, prop, minerals, kin, dt, ann
return(input)
}
## if we arrive here must be gone something wrong
stop("\n\n BUMMER!!! \n Something wrong in splitMultiKin, shutting down\n\n")
stopmsg("BUMMER!!! Something wrong in here, shutting down")
}
##' @title Distribute all values from a matrix into a buffer, EQUILIBRIUM case
##' @param input The input buffer to manipulate
##' @param prop The names of properties to distribute
##' @param values matrix or data.frame containing the values to distribute
##' @param minerals indeces of the EQUILIBRIUM minerals
##' @param SI logical, defaults to FALSE. Should we paste SI before value of a mineral?
##' @return input buffer with the values in correct places
##' @author MDL
##' @export
DistributeMatrix <- function(input, prop, values, minerals, SI=FALSE)
### Distribute properties in an input buffer
### prop is a matrix with colnames
{
if (!is.matrix(values)&&!is.data.frame(values))
stopmsg("I need a matrix to distribute")
## correct quoting
qprop <- gsub("\\(","\\\\(",prop)
qprop <- gsub("\\)","\\\\)",qprop)
n <- length(qprop)
if (n!=ncol(values)) {
stopmsg("prop names and values matrix of different dimension.")
}
if (is.character(minerals))
{
qmin <- gsub("\\(","\\\\(",minerals)
qmin <- gsub("\\)","\\\\)",qmin)
mins <- match(qmin, qprop)
}
else
mins <- as.numeric(minerals)
indAq <- seq(1:n)[-mins]
newinp <- input
for (i in indAq) {
newinp <- Distribute(newinp, prop[i], values[,i])
}
for (i in mins) {
if (SI)
newinp <- Distribute(newinp, prop[i], values[,i])
else
newinp <- Distribute(newinp, prop[i], paste("0.0 ",values[,i]))
}
return(newinp)
}
##' This function \link{DistributeKin} a matrix of properties values into a buffer
......@@ -241,12 +287,6 @@ DistributeKinMatrix <- function(input, prop, values, minerals, kin, dt, ann)
n <- length(qprop)
if (n!=ncol(values)) {
msg(" BUMMER!")
cat("prop:", paste(qprop,collapse=" "),"\n")
cat("colnames(values):", paste(qprop,collapse=", "), "\n")
cat("values:\n")
print(values)
cat("\ninput:\n")
print(input)
stopmsg("prop names and values matrix of different dimensions?")
}
......@@ -277,7 +317,6 @@ DistributeKinMatrix <- function(input, prop, values, minerals, kin, dt, ann)
}
for (i in indAq) {
## cat(paste("distributing", prop[i],", i=",i, "values",paste(values[,i], collapse=" "),"\n"))
newinp <- Distribute(newinp, prop[i], values[,i])
}
for (i in setdiff(mins, kin)) {
......@@ -291,7 +330,6 @@ DistributeKinMatrix <- function(input, prop, values, minerals, kin, dt, ann)
}
}
for (i in kin) {
## cat(paste("distributing Kinetic", prop[i],", i=",i, "values",paste(values[,i], collapse=" "),"\n"))
newinp <- DistributeKin(newinp, prop[i], values[,i], ident="-m")
}
......@@ -326,14 +364,12 @@ splitMultiFix <- function(data, procs, base, first, prop, minerals, nmax=200)
if ( ntot ==1) ## try and call normal Distribute (i.e., ##!is.matrix(data) ||
## just 1 simulation!)
{
msg("just 1 simulation, reverting back to normal Distribute\n")
msg("just 1 simulation, reverting back to normal Distribute")
BigInp <- RepSol(base, 1, first=first)
namdis <- colnames(data)
indminer <- minerals
nammin <- namdis[indminer]
## cat("nammin", nammin,"\n")
namtot <- namdis[-indminer]
## cat("namtot", namtot,"\n")
for (dis in seq_along(nammin)) {
BigInp <- Distribute(BigInp, nammin[dis], paste("0.0 ", data[1,nammin[dis] ]))
}
......@@ -342,8 +378,6 @@ splitMultiFix <- function(data, procs, base, first, prop, minerals, nmax=200)
}
return(BigInp)
}
## stop(" ERR - splitMulti: currently only Distribute.matrix supported!")
if (ntot > nmax && procs > 1) {
msg("Splitting up ", ntot, "simulations to", procs, "processes...")
......@@ -364,15 +398,14 @@ splitMultiFix <- function(data, procs, base, first, prop, minerals, nmax=200)
input <- vector(mode="list", length=procs)
attr(input,"sims") <- breaks
input[[1]] <- Distribute.matrix(BigInplong, prop, newlist[[1]], minerals=minerals)
input[[1]] <- DistributeMatrix(BigInplong, prop, newlist[[1]], minerals=minerals)
input[c(2:procs)] <- lapply(newlist[c(2:procs)],
function(x)
return(Distribute.matrix(BigInp, prop, x, minerals=minerals))
return(DistributeMatrix(BigInp, prop, x, minerals=minerals))
)
} else {
## msg("NO NEED TO PARALLELIZE...\n"))
BigInp <- RepSol(base, ntot, first=first)
input <- Distribute.matrix(input=BigInp, prop=prop, values=data, minerals=minerals)
input <- DistributeMatrix(input=BigInp, prop=prop, values=data, minerals=minerals)
}
return(input)
}
......
This diff is collapsed.
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