Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
Marco De Lucia
RedModRphree
Commits
0283c6d1
Commit
0283c6d1
authored
May 07, 2018
by
Marco De Lucia
Browse files
renamed demos for coherence
parent
b6ea51da
Changes
2
Hide whitespace changes
Inline
Side-by-side
demo/demo-eq-TrainMultiple.R
0 → 100644
View file @
0283c6d1
des_path
<-
system.file
(
"extdata"
,
"demo_equilibrium_design.csv"
,
package
=
"RedModRphree"
)
res_path
<-
system.file
(
"extdata"
,
"demo_equilibrium_result.csv"
,
package
=
"RedModRphree"
)
des
<-
data.table
::
fread
(
des_path
)
res
<-
data.table
::
fread
(
res_path
)
## We create a directory where to dump the created data structures
path_dumps
<-
"TrainedModels/"
## if it does not exist we create it
if
(
!
dir.exists
(
path_dumps
))
dir.create
(
path_dumps
,
showWarnings
=
FALSE
)
Train
(
c
(
"PolyMARS"
,
"glmStepAIC"
,
"MARS"
),
preprocessing_ind
=
c
(
17
),
input_data
=
des
,
output_data
=
res
,
seed
=
2305
,
train_para
=
TRUE
,
run_para
=
TRUE
,
preproc_para
=
FALSE
,
selection_criteria
=
"MAD"
,
write_full_residuals
=
FALSE
,
tuner
=
FALSE
,
use_cores
=
8
,
training_samples
=
0.95
,
write_filepath_prefix
=
path_dumps
)
a
<-
SelectedSurrogate
(
des
)
SelectedMethods
<-
GetModelNames
(
model_list
)
names
(
SelectedMethods
)
<-
colnames
(
res
)
par
(
mfrow
=
c
(
4
,
2
))
for
(
i
in
colnames
(
des
))
{
plot
(
res
[,
..i
][[
1
]],
a
[,
i
],
"p"
,
pch
=
3
,
main
=
paste
(
i
,
SelectedMethods
[
i
],
sep
=
", "
),
xlab
=
"Full physics"
,
ylab
=
"Surrogate"
)
abline
(
0
,
1
,
lty
=
"dashed"
,
col
=
"grey"
)
}
demo/demo-kin-surr-RF.R
0 → 100644
View file @
0283c6d1
## Licence: LGPL version 2.1
## Time-stamp: "Last modified 2018-05-06 21:29:53 delucia"
library
(
RedModRphree
)
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"
,
"pH 9.91 charge"
,
"pe 4.0"
,
"C 1.2279E-04"
,
"Ca 1.2279E-04"
,
"Mg 0.001"
,
"Cl 0.002"
,
"PURE 1"
,
"O2g -0.1675 10"
,
"KINETICS 1"
,
"-steps 100"
,
"-bad_step_max 1000"
,
"Calcite"
,
"-m 0.000207"
,
"-parms 3.20"
,
"Dolomite"
,
"-m 0.0"
,
"-parms 0.32"
,
"END"
)
initsim
<-
c
(
"SELECTED_OUTPUT"
,
"-high_precision true"
,
"-reset false"
,
"USER_PUNCH"
,
"-head C Ca Cl Mg H e O2g Calcite Dolomite"
,
"10 PUNCH TOT(\"C\"), TOT(\"Ca\"), TOT(\"Cl\"), TOT(\"Mg\"), ACT(\"H+\"), ACT(\"e-\"), EQUI(\"O2g\"), EQUI(\"Calcite\"), EQUI(\"Dolomite\")"
,
"SOLUTION 1"
,
"units mol/kgw"
,
"temp 25.0"
,
"water 1"
,
"pH 9.91 charge"
,
"pe 4.0"
,
"C 1.2279E-04"
,
"Ca 1.2279E-04"
,
"Cl 1E-12"
,
"Mg 1E-12"
,
"PURE 1"
,
"O2g -0.6788 10.0"
,
"Calcite 0.0 2.07E-3"
,
"Dolomite 0.0 0.0"
,
"END"
)
selout
<-
c
(
"SELECTED_OUTPUT"
,
"-high_precision true"
,
"-reset false"
,
"-time"
,
"-soln"
,
"-temperature true"
,
"-water true"
,
"-pH"
,
"-pe"
,
"-totals C Ca Cl Mg"
,
"-kinetic_reactants Calcite Dolomite"
,
"-equilibrium O2g"
)
prop
<-
c
(
"C"
,
"Ca"
,
"Cl"
,
"Mg"
,
"pH"
,
"pe"
,
"O2g"
,
"Calcite"
,
"Dolomite"
)
maxsim
<-
20
## minimum number of single chemical simulations in each time step to spawn parallelization
grid
<-
50
## grid discretization (number of elements)
L
<-
0.5
## total system length [m]
U
<-
9.375e-6
## imposed (constant) Darcy velocity [m3/s]
Cf
<-
0.9
## Courant number to add some numerical dispersion
## boundary conditions (concentration of the injected fluid entering into the 1D column)
bound
<-
c
(
H
=
1E-7
,
e
=
1E-4
,
Mg
=
0.001
,
Cl
=
0.002
)
exper
<-
c
(
50
,
100
)
expsim
<-
vector
(
mode
=
"list"
,
length
=
length
(
exper
))
timexpsim
<-
vector
(
mode
=
"list"
,
length
=
length
(
exper
))
for
(
i
in
seq_along
(
exper
))
{
setup
<-
list
(
n
=
exper
[
i
],
len
=
L
,
U
=
U
,
base
=
base
,
first
=
selout
,
Cf
=
Cf
,
bound
=
bound
,
dt
=
210
,
prop
=
prop
,
immobile
=
c
(
7
,
8
,
9
),
kin
=
c
(
8
,
9
),
ann
=
list
(
O2g
=
-0.1675
),
initsim
=
initsim
)
timexpsim
[[
i
]]
<-
system.time
({
expsim
[[
i
]]
<-
ReactTranspBalanceKin
(
setup
=
setup
,
db
=
db
,
step
=
"fix_dt"
,
maxtime
=
21000
,
procs
=
4
,
maxsim
=
maxsim
,
writeout
=
FALSE
,
reduce
=
FALSE
,
surrogate
=
FALSE
)
})
}
par
(
mfrow
=
c
(
2
,
1
))
MatplotSingle
(
sim
=
expsim
[[
1
]],
xlab
=
"Domain element"
,
ylab
=
"Ca, Mg [molal], Calcite, Dolomite [mol]"
,
main
=
"Kinetic simulation, grid 50"
)
MatplotSingle
(
sim
=
expsim
[[
2
]],
xlab
=
"Domain element"
,
ylab
=
"Ca, Mg [molal], Calcite, Dolomite [mol]"
,
main
=
"Kinetic simulation, grid 100"
)
sam
<-
ExtractSamples
(
expsim
,
var_design
=
c
(
"C"
,
"Ca"
,
"Cl"
,
"Mg"
,
"pH"
,
"pe"
,
"Calcite"
,
"Dolomite"
),
var_result
=
c
(
"C"
,
"Ca"
,
"Cl"
,
"Mg"
,
"pH"
,
"pe"
,
"Calcite"
,
"Dolomite"
))
## library(data.table)
## data.table::fwrite(sam$design, file = "RedModRphree/inst/extdata/demo_kinetic_design.csv", append = FALSE, quote = "auto",
## sep = ",", eol = "\n", na = "NA", dec = ".", row.names = FALSE, col.names = TRUE)
## data.table::fwrite(sam$result, file = "RedModRphree/inst/extdata/demo_kinetic_result.csv", append = FALSE, quote = "auto",
## sep = ",", eol = "\n", na = "NA", dec = ".", row.names = FALSE, col.names = TRUE)
set.seed
(
3560
)
library
(
caret
)
fitControl
<-
trainControl
(
method
=
"none"
,
allowParallel
=
TRUE
,
verboseIter
=
TRUE
,
returnData
=
FALSE
,
returnResamp
=
"none"
,
savePredictions
=
"none"
,
predictionBounds
=
rep
(
TRUE
,
2
))
modk
<-
parallel
::
mclapply
(
colnames
(
sam
$
result
),
function
(
x
)
return
(
train
(
x
=
sam
$
design
,
y
=
sam
$
result
[,
x
],
trControl
=
fitControl
,
##tuneGrid = tg,
method
=
"parRF"
)),
mc.cores
=
8
)
names
(
modk
)
<-
colnames
(
sam
$
result
)
saveRDS
(
modk
,
"models_kinetics_RF.rds"
)
PlotModsInSample
(
modk
,
sam
$
result
,
sam
$
design
)
baleq
<-
list
(
Calcite
=
c
(
C
=
1
,
Ca
=
1
),
Dolomite
=
c
(
C
=
2
,
Ca
=
1
,
Mg
=
1
))
balance
<-
CheckBalance
(
baleq
,
sam
$
design
,
sam
$
result
)
MySurrKin
<-
function
(
state
,
model
){
excl
<-
which
(
colnames
(
state
)
==
"O2g"
)
remember_state_names
<-
colnames
(
state
)
des
<-
state
[,
-
excl
]
if
(
NA
%in%
des
)
{
print
(
des
)
stop
(
"Mysurr: NA in data"
)
}
if
(
NaN
%in%
des
)
{
print
(
des
)
stop
(
"Mysurr: NaN in data"
)
}
order_for_surrogate
<-
c
(
"C"
,
"Ca"
,
"Cl"
,
"Mg"
,
"pH"
,
"pe"
,
"Calcite"
,
"Dolomite"
)
rem_attr
<-
attr
(
state
,
"immobile"
)
ss
<-
des
[,
order_for_surrogate
]
pred
<-
lapply
(
names
(
model
),
function
(
x
)
as.numeric
(
predict.train
(
model
[[
x
]],
ss
)))
out
<-
cbind
(
as.data.frame
(
pred
),
state
[,
excl
])
colnames
(
out
)
<-
c
(
names
(
model
),
"O2g"
)
## order back to how it is expected in "state"
out
<-
out
[,
remember_state_names
]
attr
(
out
,
"immobile"
)
<-
rem_attr
return
(
out
)
}
exper
<-
expand.grid
(
grid
=
c
(
50
,
100
),
tol
=
c
(
1E-6
,
1E-5
))
simsur
<-
vector
(
mode
=
"list"
,
length
=
nrow
(
exper
))
timesur
<-
vector
(
mode
=
"list"
,
length
=
nrow
(
exper
))
outbalsur
<-
vector
(
mode
=
"list"
,
length
=
nrow
(
exper
))
for
(
i
in
seq
(
nrow
(
exper
)))
{
setup
<-
list
(
n
=
exper
$
grid
[
i
],
len
=
L
,
U
=
U
,
base
=
base
,
first
=
selout
,
Cf
=
Cf
,
bound
=
bound
,
dt
=
210
,
prop
=
prop
,
immobile
=
c
(
7
,
8
,
9
),
kin
=
c
(
8
,
9
),
ann
=
list
(
O2g
=
-0.1675
),
initsim
=
initsim
)
timesur
[[
i
]]
<-
system.time
(
simsur
[[
i
]]
<-
ReactTranspBalanceKin
(
setup
=
setup
,
db
=
db
,
step
=
"fix_dt"
,
maxtime
=
21000
,
procs
=
4
,
maxsim
=
maxsim
,
writeout
=
FALSE
,
reduce
=
FALSE
,
surrogate
=
TRUE
,
surrogate.FUN
=
MySurrKin
,
model
=
modk
,
baleq
=
baleq
,
tol
=
exper
$
tol
[
i
],
call_pqc
=
TRUE
)
)
outbalsur
[[
i
]]
<-
out_bal
}
par
(
mfrow
=
c
(
2
,
2
))
inds
<-
rep
(
c
(
1
,
2
),
2
)
for
(
i
in
seq
(
nrow
(
exper
)))
{
Matplot
(
sur
=
simsur
[[
i
]],
sim
=
expsim
[[
inds
[
i
]
]],
xlab
=
"Domain element"
,
ylab
=
"Ca, Mg [molal], Calcite, Dolomite [mol]"
,
main
=
paste
(
"RandomForest surrogate, tol="
,
exper
$
tol
[
i
]))
}
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment