-
Notifications
You must be signed in to change notification settings - Fork 0
/
fullModelConcentrationEstimatorWrapper.py
178 lines (157 loc) · 7.45 KB
/
fullModelConcentrationEstimatorWrapper.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
############################################################################
#
# Imperial College London, United Kingdom
# Multifunctional Nanomaterials Laboratory
#
# Project: ERASE
# Year: 2021
# Python: Python 3.7
# Authors: Ashwin Kumar Rajagopalan (AK)
#
# Purpose:
# Wrapper function to estimate the concentration given a time-resolved
# measurement from the full model
#
# Last modified:
# - 2021-02-03, AK: Change inputParameters (add adsorbent density)
# - 2021-01-26, AK: Add noise to true measurement
# - 2021-01-25, AK: Integrate full model concentration estimator
# - 2021-01-21, AK: Initial creation
#
# Input arguments:
#
#
# Output arguments:
#
#
############################################################################
import auxiliaryFunctions
from simulateFullModel import simulateFullModel
from estimateConcentration import estimateConcentration
from estimateConcentrationFullModel import estimateConcentrationFullModel
from tqdm import tqdm # To track progress of the loop
import numpy as np
from numpy import savez
from joblib import Parallel, delayed # For parallel processing
import multiprocessing
import os
import time
import socket
# Find out the total number of cores available for parallel processing
num_cores = multiprocessing.cpu_count()
# Get the commit ID of the current repository
gitCommitID = auxiliaryFunctions.getCommitID()
# Get the current date and time for saving purposes
currentDT = auxiliaryFunctions.getCurrentDateTime()
# Flag to determine whether concentration estimated accounted for kinetics
# (False) or not (True)
flagIndTime = True
# Flag to determine whether constant pressure (False) or constant flow rate
# model to be used (True)
modelConstF = False
# Sensor ID
sensorID = [6,2] # [-]
# Number of Gases
numberOfGases = 2 # [-]
# Rate Constant
rateConstant = ([.01,.01,10000.0],
[.01,.01,10000.0]) # [1/s]
# Feed mole fraction
feedMoleFrac = [0.1,0.9,0.0] # [-]
# Time span for integration [tuple with t0 and tf] [s]
timeInt = (0,1000) # [s]
# Volumetric flow rate [m3/s]
flowIn = 5e-7 # [s]
# Measurement noise characteristics
meanNoise = 0.0 # [g/kg]
stdNoise = 0.0 # [g/kg]
# Loop over all rate constants
outputStruct = {}
for ii in tqdm(range(len(sensorID))):
# Call the full model with a given rate constant
timeSim, _ , sensorFingerPrint, inputParameters = simulateFullModel(sensorID = sensorID[ii],
rateConstant = rateConstant[ii],
feedMoleFrac = feedMoleFrac,
timeInt = timeInt,
flowIn = flowIn)
# Generate the noise data for all time instants and all materials
measurementNoise = np.random.normal(meanNoise, stdNoise, len(timeSim))
sensorFingerPrintRaw = sensorFingerPrint # Raw sensor finger print
# Add measurement noise to the sensor finger print
sensorFingerPrint = sensorFingerPrint + measurementNoise
outputStruct[ii] = {'timeSim':timeSim,
'sensorFingerPrint':sensorFingerPrint,
'sensorFingerPrintRaw':sensorFingerPrintRaw, # Without noise
'measurementNoise':measurementNoise,
'inputParameters':inputParameters} # Check simulateFullModel.py for entries
# Prepare time-resolved sendor finger print
timeSim = []
timeSim = outputStruct[0]["timeSim"]
sensorFingerPrint = np.zeros([len(timeSim),len(sensorID)])
for ii in range(len(sensorID)):
sensorFingerPrint[:,ii] = outputStruct[ii]["sensorFingerPrint"]
# flagIndTime - If true, each time instant evaluated without knowledge of
# actual kinetics in the estimate of the concentration (accounted only in the
# generation of the true sensor response above) - This would lead to a
# scenario of concentrations evaluated at each time instant
if flagIndTime:
# Start time for time elapsed
startTime = time.time()
# Initialize output matrix
arrayConcentration = np.zeros([len(timeSim),numberOfGases+len(sensorID)])
# Loop over all time instants and estimate the gas composition
arrayConcentrationTemp = Parallel(n_jobs=num_cores)(delayed(estimateConcentration)
(None,numberOfGases,None,sensorID,
fullModel = True,
fullModelResponse = sensorFingerPrint[ii,:])
for ii in tqdm(range(len(timeSim))))
# Convert the list to array
arrayConcentration = np.array(arrayConcentrationTemp)
# Stop time for time elapsed
stopTime = time.time()
# Time elapsed [s]
timeElapsed = stopTime - startTime
# flagIndTime - If false, only one concentration estimate obtained. The full
# model is simulated in concentration estimator accounting for the kinetics.
else:
# Start time for time elapsed
startTime = time.time()
# Initialize output matrix
arrayConcentration = np.zeros([len(timeSim)-1,numberOfGases+1+len(sensorID)])
# Loop over all time instants and estimate the gas composition
# The concentration is estimated for all the time instants in timeSim.
# This simulates using the full model as and when a new measurement is
# available. This assumes that the estimator has knowledge of the past
# measurements
# The first time instant is not simulated. This is fixed by adding a dummy
# row after the simulations are over
arrayConcentrationTemp = Parallel(n_jobs=num_cores)(delayed(estimateConcentrationFullModel)
(numberOfGases,sensorID,
fullModelResponse = sensorFingerPrint[0:ii+1,:],
rateConstant = rateConstant,
timeInt = (0,timeSim[ii+1]),flowIn = flowIn)
for ii in tqdm(range(len(timeSim)-1)))
# Convert the list to array
arrayConcentration = np.array(arrayConcentrationTemp)
# Add dummy row to be consistent (intialized to init condition)
firstRow = np.concatenate((np.array(sensorID), inputParameters[6]))
arrayConcentration = np.vstack([firstRow,arrayConcentration])
# Stop time for time elapsed
stopTime = time.time()
# Time elapsed [s]
timeElapsed = stopTime - startTime
# Check if simulationResults directory exists or not. If not, create the folder
if not os.path.exists('simulationResults'):
os.mkdir('simulationResults')
# Save the array concentration into a native numpy file
# The .npz file is saved in a folder called simulationResults (hardcoded)
filePrefix = "fullModelConcentrationEstimate"
sensorText = str(sensorID).replace('[','').replace(']','').replace(' ','-').replace(',','')
saveFileName = filePrefix + "_" + sensorText + "_" + currentDT + "_" + gitCommitID;
savePath = os.path.join('simulationResults',saveFileName)
savez (savePath, outputStruct = outputStruct, # True response
arrayConcentration = arrayConcentration, # Estimated response
eqbmModelFlag = flagIndTime, # Flag to tell whether eqbm. or full model used
noiseCharacteristics = [meanNoise, stdNoise], # Noise characteristics
timeElapsed = timeElapsed, # Time elapsed for the simulation
hostName = socket.gethostname()) # Hostname of the computer