From 0d6f3cf3407d19085de8a9481f2f1e8953474e0d Mon Sep 17 00:00:00 2001
From: farid <farid.mohammadi@iws.uni-stuttgart.de>
Date: Fri, 7 Jan 2022 12:23:20 +0100
Subject: [PATCH] [PA-A][PNM] update script with only one Transmissibility
param.
---
.../tests/PA-A/ffpm_validation_stokespnm.py | 48 ++++++++++---------
1 file changed, 25 insertions(+), 23 deletions(-)
diff --git a/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py b/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py
index dfb0e0c3d..3d09b8ddc 100755
--- a/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py
+++ b/BayesValidRox/tests/PA-A/ffpm_validation_stokespnm.py
@@ -115,18 +115,22 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
Inputs.addMarginals() #VyMaxTop
Inputs.Marginals[0].Name = '$V^{top}$' # [5e-4, 1.5e-3]
Inputs.Marginals[0].InputValues = stats.uniform(loc=5.0e-04, scale=1.5e-03-5.0e-04).rvs(size=MCSize)
-
- Inputs.addMarginals() #TransmissibilityThroat '$g_{t,ij}$'
- Inputs.Marginals[1].Name = '$g_{t,ij}$' # [1.0e-7, 1.0e-04]
- Inputs.Marginals[1].InputValues = stats.uniform(loc=1.0e-7, scale=1.0e-04-1.0e-7).rvs(size=MCSize)
-
- Inputs.addMarginals() #TransmissibilityHalfPore '$g_{p,i}$'
- Inputs.Marginals[2].Name = '$g_{p,i}$' # [1.0e-7, 1.0e-04]
- Inputs.Marginals[2].InputValues = stats.uniform(loc=1.0e-7, scale=1.0e-04-1.0e-7).rvs(size=MCSize)
+
+ Inputs.addMarginals() #TransmissibilityTotal '$g_{t,}$'
+ Inputs.Marginals[1].Name = '$g_{t}$' # [1.0e-7, 1.0e-05]
+ Inputs.Marginals[1].InputValues = stats.uniform(loc=1e-7, scale=1e-5-1e-7).rvs(size=MCSize)
+
+ # Inputs.addMarginals() #TransmissibilityThroat '$g_{t,ij}$'
+ # Inputs.Marginals[1].Name = '$g_{t,ij}$' # [1.0e-7, 1.0e-04]
+ # Inputs.Marginals[1].InputValues = stats.uniform(loc=1e-7, scale=1e-4-1e-7).rvs(size=MCSize)
+
+ # Inputs.addMarginals() #TransmissibilityHalfPore '$g_{p,i}$'
+ # Inputs.Marginals[2].Name = '$g_{p,i}$' # [1.0e-7, 1.0e-04]
+ # Inputs.Marginals[2].InputValues = stats.uniform(loc=1e-7, scale=1e-4-1e-7).rvs(size=MCSize)
Inputs.addMarginals() #'$\\beta_{pore}$'
- Inputs.Marginals[3].Name = '$\\beta_{pore}$' # [1000.0, 100000.0]
- Inputs.Marginals[3].InputValues = stats.uniform(loc=1e3, scale=1e5-1e3).rvs(size=MCSize)
+ Inputs.Marginals[2].Name = '$\\beta_{pore}$' # [1000.0, 100000.0]
+ Inputs.Marginals[2].InputValues = stats.uniform(loc=1e3, scale=1e5-1e3).rvs(size=MCSize)
#=====================================================
#========== DEFINITION OF THE METAMODEL ============
@@ -181,10 +185,10 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
# Sampling methods
# 1) random 2) latin_hypercube 3) sobol 4) halton 5) hammersley 6) chebyshev(FT)
# 7) korobov 8) grid(FT) 9) nested_grid(FT) 10)user
- MetaModelOpts.ExpDesign.SamplingMethod = 'user'
+ MetaModelOpts.ExpDesign.SamplingMethod = 'latin_hypercube'
# Provide the experimental design object with a hdf5 file
- MetaModelOpts.ExpDesign.hdf5File = "ExpDesign_{}.hdf5".format(Model.Name)
+ # MetaModelOpts.ExpDesign.hdf5File = "ExpDesign_{}.hdf5".format(Model.Name)
# For calculation of validation error for SeqDesign
import h5py
@@ -215,7 +219,7 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
# Remove zip file
if os.path.isfile(Model.Name + '.zip'): os.remove(Model.Name + '.zip')
- # # Save PCE models
+ # Save PCE models
with open('PCEModel_'+Model.Name+'.pkl', 'wb') as output:
joblib.dump(PCEModel, output, 2)
@@ -253,12 +257,10 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
# Select the inference method
BayesOptsCalib.SamplingMethod = 'MCMC'
BayesOptsCalib.MCMCnwalkers = 30
- # BayesOptsCalib.selectedIndices = {'velocity [m/s]':[0,1,2,3,4,5,6,7],
+ # BayesOptsCalib.selectedIndices = {'velocity [m/s]':[0,2,3,4,6,7,8,9],
# 'p':[0,1,2]}
# BayesOptsCalib.ReqOutputType = ['velocity [m/s]']
- # import emcee
- # BayesOptsCalib.MCMCmoves = emcee.moves.KDE()
- BayesOptsCalib.MCMCnSteps = 200000
+ # BayesOptsCalib.MCMCnSteps = 50000#00
# BayesOptsCalib.MCMCverbose = True
# Maximum a Posteriori based on primary analysis
# theta_MAP = np.array([1e-3,5.26480244078e-06,1.7642875191945115e-05,
@@ -293,11 +295,11 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
# Mean absolute error
# delta_dict[out] = np.mean(abs(data - ED_Y[out]),axis=0)
# Mean squared error
- delta_dict[out] = np.mean((data - ED_Y[out])**2,axis=0)
- # delta_dict[out] = np.ones_like(data)
+ # delta_dict[out] = np.mean((data - ED_Y[out])**2,axis=0)
+ delta_dict[out] = np.ones_like(data)
errormodels = pd.DataFrame.from_dict(delta_dict,orient='index').T
- DiscrepancyOpts.Parameters = errormodels/errormodels.sum()
+ DiscrepancyOpts.Parameters = errormodels#/errormodels.sum()
BayesOptsCalib.Discrepancy = DiscrepancyOpts
@@ -375,7 +377,7 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
# Plot posterior predictive
postPredictiveplot(PCEModel.ModelObj.Name, errorPerc, averaging,
inletLoc=inletLoc, case='Calib', bins=20)
-
+ return BayesCalib
#=====================================================
#================== VALIDATION =====================
#=====================================================
@@ -591,11 +593,11 @@ def run(params, averaging=True,errorPerc=0.05, inletLoc='top', PCEEDMethod='norm
#==========================================================================
# ---------- set BayesValidRox params ----------
# PCEExpDesignMethod = 'normal' # 'normal' or 'sequential'
-# nInitSamples = 10 #50 Initial No. of orig. Model runs for surrogate training
+# nInitSamples = 100 #50 Initial No. of orig. Model runs for surrogate training
# nTotalSamples = 200 #100 Total No. of orig. Model runs for surrogate training
# nBootstrapItr = 1000 # No. of bootstraping iterations for Bayesian analysis
# BootstrapNoise = 0.005 # Noise amount for bootstraping in Bayesian analysis
-# perturbedData = np.loadtxt('./data/perturbedValidData_squared_inclusion_topInflow.csv',delimiter=',')
+perturbedData = np.loadtxt('./data/perturbedValidData_squared_inclusion_topInflow.csv',delimiter=',')
# perturbedData = np.loadtxt('./data/perturbedValidData_squared_inclusion_leftInflow.csv',delimiter=',')
# perturbedData = np.loadtxt('./data/perturbedValidDataAvg_squared_inclusion_topInflow.csv',delimiter=',')
# perturbedData = np.loadtxt('./data/perturbedValidDataAvg_squared_inclusion_leftInflow.csv',delimiter=',')
--
GitLab