Skip to content

Commit 45da2f3

Browse files
bef-18bef-18
authored
Merge pull request #153 from ArnauMiro/develop
Develop to main
2 parents e7d4658 + 5fe0a8e commit 45da2f3

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

43 files changed

+2522
-201
lines changed

Diff for: ‎.gitignore

+6-2
Original file line numberDiff line numberDiff line change
@@ -184,6 +184,10 @@ Examples/Data/*
184184
flow/*
185185
events*
186186
model_state
187-
*.out
188-
*.err
189187
submit*.sh
188+
*.err
189+
*.out
190+
*.json
191+
*.png
192+
*.pth
193+
*.pdf

Diff for: ‎Examples/NN/SHRED/example_01_SVD_and_sensors_cylinder.py

+68
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,68 @@
1+
#!/usr/bin/env python
2+
#
3+
# Example of POD to reduce dimensionality for SHRED architecture.
4+
#
5+
# Last revision: 19/07/2021
6+
from __future__ import print_function, division
7+
8+
import mpi4py
9+
mpi4py.rc.recv_mprobe = False
10+
11+
import numpy as np
12+
import pyLOM
13+
14+
15+
## Parameters
16+
DATAFILE = '/gpfs/scratch/bsc21/bsc021828/DATA_PYLOM/CYLINDER.h5'
17+
DATAFIL2 = './CYLINDER.h5'
18+
VARLIST = ['VELOX', 'VORTI']
19+
20+
21+
## Data loading
22+
m = pyLOM.Mesh.load(DATAFILE)
23+
d = pyLOM.Dataset.load(DATAFILE,ptable=m.partition_table)
24+
t = d.get_variable('time')
25+
N = t.shape[0]
26+
27+
28+
## Divide in training, validation and test and append mask to current dataset
29+
tridx, vaidx, teidx = d.split_data('time', mode='reconstruct')
30+
pyLOM.pprint(0, tridx.shape, vaidx.shape, teidx.shape, t.shape)
31+
# Regenerate output datafile with the variable masks
32+
m.save(DATAFIL2, nopartition=True, mode='w')
33+
d.save(DATAFIL2, nopartition=True)
34+
35+
36+
## Extract sensors
37+
# Generate random sensors
38+
nsens = 20 # Number of sensors
39+
x0, x1 = 0.5, 8 # Bounds at the X axis of the region where the sensor will be located
40+
y0, y1 = -1, 1 # Bounds at the Y axis of the region where the sensor will be located
41+
bounds = np.array([x0,x1,y0,y1])
42+
dsens = d.select_random_sensors(nsens, bounds, VARLIST)
43+
# Save the sensor dataset
44+
dsens.save('sensors.h5', nopartition=True, mode='w')
45+
46+
47+
## Compute POD separately for each variable in order to reduce memory usage during the SVD. POD is computed only for the training dataset. Validation and test are projected to the POD modes
48+
for var in VARLIST:
49+
## Fetch training dataset and compute POD
50+
Xtrai = d.mask_field(var, tridx)
51+
PSI,S,V = pyLOM.POD.run(Xtrai,remove_mean=True,randomized=True,r=8,q=3)
52+
## Save POD modes for training of each variable
53+
pyLOM.POD.save('POD_trai_%s.h5'%var,PSI,S,V,d.partition_table,nvars=1,pointData=d.point)
54+
## Fetch validation dataset and project POD modes
55+
Xvali = d.mask_field(var, vaidx)
56+
proj = pyLOM.math.matmulp(PSI.T, Xvali)
57+
Vvali = pyLOM.math.matmul(pyLOM.math.diag(1/S), proj)
58+
## Save POD projection of validation data of each variable
59+
pyLOM.POD.save('POD_vali_%s.h5'%var,None,None,Vvali,d.partition_table,nvars=1,pointData=d.point)
60+
## Fetch test dataset and project POD modes
61+
Xtest = d.mask_field(var, teidx)
62+
proj = pyLOM.math.matmulp(PSI.T, Xtest)
63+
Vtest = pyLOM.math.matmul(pyLOM.math.diag(1/S), proj)
64+
## Save POD projection of validation data of each variable
65+
pyLOM.POD.save('POD_test_%s.h5'%var,None,None,Vtest,d.partition_table,nvars=1,pointData=d.point)
66+
67+
68+
pyLOM.cr_info()

Diff for: ‎Examples/NN/SHRED/example_02_fit_shred_cylinder.py

+98
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,98 @@
1+
import numpy as np
2+
import torch
3+
import pyLOM, pyLOM.NN
4+
5+
pyLOM.style_plots()
6+
7+
8+
## Set device
9+
device = pyLOM.NN.select_device() # Force CPU for this example, if left in blank it will automatically select the device
10+
11+
12+
## Input parameters
13+
# Data paths
14+
sensvar = 'VELOX' # Variable from the sensor measurements we'll be working with
15+
podvar = 'VELOX' # Variable from the POD we'll be working with
16+
17+
# Output paths
18+
inscaler = 'out/scalers/config_'
19+
outscale = 'out/scalers/pod.json'
20+
shreds = 'out/shreds/config_'
21+
22+
# SHRED sensor configurations for uncertainty quantification
23+
nconfigs = 20
24+
25+
26+
## Import sensor measurements
27+
sensors = pyLOM.Dataset.load('sensors.h5')
28+
mask_trai = sensors.get_variable('training_time')
29+
mask_vali = sensors.get_variable('validation_time')
30+
mask_test = sensors.get_variable('test_time')
31+
# Training
32+
sens_trai = sensors.mask_field(sensvar, mask_trai)
33+
nsens = sens_trai.shape[0]
34+
# Validation
35+
sens_vali = sensors.mask_field(sensvar, mask_vali)
36+
# Test
37+
sens_test = sensors.mask_field(sensvar, mask_test)
38+
# Compute total timesteps
39+
time = sensors.get_variable('time')
40+
ntimeG = time.shape[0]
41+
42+
43+
## Import POD coefficients and rescale them.
44+
# Training
45+
S, pod_trai = pyLOM.POD.load('POD_trai_%s.h5' % podvar, vars=['S','V'])
46+
pod_trai = pod_trai.astype(np.float32)
47+
pod_scaler = pyLOM.NN.MinMaxScaler(column=True)
48+
pod_scaler.fit(pod_trai)
49+
pod_scaler.save(outscale)
50+
trai_out = pod_scaler.transform(pod_trai)
51+
output_size = trai_out.shape[0]
52+
Sscale = S/np.sum(S)
53+
# Validation
54+
pod_vali = pyLOM.POD.load('POD_vali_%s.h5' % podvar, vars='V')[0].astype(np.float32)
55+
vali_out = pod_scaler.transform(pod_vali)
56+
# Test
57+
pod_test = pyLOM.POD.load('POD_test_%s.h5' % podvar, vars='V')[0].astype(np.float32)
58+
test_out = pod_scaler.transform(pod_test)
59+
# Full POD
60+
full_pod = np.zeros((output_size,ntimeG), dtype=pod_trai.dtype)
61+
full_pod[:,mask_trai] = pod_trai
62+
full_pod[:,mask_vali] = pod_vali
63+
full_pod[:,mask_test] = pod_test
64+
65+
66+
## Build SHRED architecture
67+
shred = pyLOM.NN.SHRED(output_size, device, nsens, nconfigs=nconfigs)
68+
69+
70+
## Fit all SHRED configurations using the data from the sensors
71+
for kk, mysensors in enumerate(shred.configs):
72+
# Get the sensor values for the training, validation and test
73+
mytrai = sens_trai[mysensors,:]
74+
myvali = sens_vali[mysensors,:]
75+
mytest = sens_test[mysensors,:]
76+
# Scale the data
77+
myscaler = pyLOM.NN.MinMaxScaler()
78+
scalpath = '%s%i.json' % (inscaler, kk)
79+
myscaler.fit(mytrai)
80+
myscaler.save(scalpath)
81+
trai_sca = myscaler.transform(mytrai)
82+
vali_sca = myscaler.transform(myvali)
83+
test_sca = myscaler.transform(mytest)
84+
# Concatenate train, test and validation data to generate the embeddings correctly
85+
embedded = np.zeros((trai_sca.shape[0],ntimeG), dtype=trai_sca.dtype)
86+
embedded[:,mask_trai] = trai_sca
87+
embedded[:,mask_vali] = vali_sca
88+
embedded[:,mask_test] = test_sca
89+
delayed = pyLOM.math.time_delay_embedding(embedded, dimension=50)
90+
# Generate training validation and test datasets both for reconstruction of states
91+
train_dataset = pyLOM.NN.Dataset(trai_out, variables_in=delayed[:,mask_trai,:])
92+
valid_dataset = pyLOM.NN.Dataset(vali_out, variables_in=delayed[:,mask_vali,:])
93+
# Fit SHRED
94+
shred.fit(train_dataset, valid_dataset, epochs=1500, patience=100, verbose=False, mod_scale=torch.tensor(Sscale))
95+
shred.save('%s%i' % (shreds,kk), scalpath, outscale, mysensors)
96+
97+
98+
pyLOM.cr_info()

Diff for: ‎Examples/NN/SHRED/example_03_inference_shred_cylinder.py

+120
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,120 @@
1+
import numpy as np
2+
import torch
3+
import pyLOM, pyLOM.NN
4+
5+
pyLOM.style_plots()
6+
7+
8+
## Set device
9+
device = pyLOM.NN.select_device() # Force CPU for this example, if left in blank it will automatically select the device
10+
11+
12+
## Input parameters
13+
# Data paths
14+
sensvar = 'VELOX' # Variable from the sensor measurements we'll be working with
15+
podvar = 'VELOX' # Variable from the POD we'll be working with
16+
shreds = 'out/shreds/config_' # Where the SHREDS models are saved
17+
18+
# Select which sensor configurations have to be postprocessed
19+
nconfigs = 20
20+
configIDs = np.arange(nconfigs, dtype=int)
21+
22+
23+
## Import sensor measurements
24+
sensors = pyLOM.Dataset.load('sensors.h5')
25+
mask_trai = sensors.get_variable('training_time')
26+
mask_vali = sensors.get_variable('validation_time')
27+
mask_test = sensors.get_variable('test_time')
28+
# Training
29+
sens_trai = sensors.mask_field(sensvar, mask_trai)
30+
nsens = sens_trai.shape[0]
31+
# Validation
32+
sens_vali = sensors.mask_field(sensvar, mask_vali)
33+
# Test
34+
sens_test = sensors.mask_field(sensvar, mask_test)
35+
# Compute total timesteps
36+
time = sensors.get_variable('time')
37+
ntimeG = time.shape[0]
38+
39+
40+
## Import POD coefficients and rescale them.
41+
# Training
42+
S, pod_trai = pyLOM.POD.load('POD_trai_%s.h5' % podvar, vars=['S','V'])
43+
pod_trai = pod_trai.astype(np.float32)
44+
Sscale = S/np.sum(S)
45+
output_size = pod_trai.shape[0]
46+
# Validation
47+
pod_vali = pyLOM.POD.load('POD_vali_%s.h5' % podvar, vars='V')[0].astype(np.float32)
48+
# Test
49+
pod_test = pyLOM.POD.load('POD_test_%s.h5' % podvar, vars='V')[0].astype(np.float32)
50+
# Full POD
51+
full_pod = np.zeros((output_size,ntimeG), dtype=np.float32)
52+
full_pod[:,mask_trai] = pod_trai
53+
full_pod[:,mask_vali] = pod_vali
54+
full_pod[:,mask_test] = pod_test
55+
56+
57+
## Build SHRED architecture
58+
shred = pyLOM.NN.SHRED(output_size, device, nsens)
59+
60+
61+
## Load all SHRED configurations and evaluate them using the data from the corresponding sensors
62+
outres = np.zeros((output_size, ntimeG, nconfigs), dtype=full_pod.dtype)
63+
MRE = np.zeros((output_size, nconfigs), dtype=full_pod.dtype)
64+
for kk in range(nconfigs):
65+
# Load the SHRED model for this configuration
66+
load_dict = torch.load('%s%i.pth' % (shreds, kk))
67+
mysensors = load_dict['sensors']
68+
shred.load_state_dict(load_dict['model_state_dict'])
69+
# Get the sensor values for the training, validation and test
70+
mytrai = sens_trai[mysensors,:]
71+
myvali = sens_vali[mysensors,:]
72+
mytest = sens_test[mysensors,:]
73+
# Load the appropiate scaler
74+
myscaler = pyLOM.NN.MinMaxScaler(column=True)
75+
myscaler = myscaler.load(load_dict['scaler_path'])
76+
trai_sca = myscaler.transform(mytrai)
77+
vali_sca = myscaler.transform(myvali)
78+
test_sca = myscaler.transform(mytest)
79+
# Concatenate train, test and validation data to generate the embeddings correctly
80+
embedded = np.zeros((trai_sca.shape[0],ntimeG), dtype=np.float32)
81+
embedded[:,mask_trai] = trai_sca
82+
embedded[:,mask_vali] = vali_sca
83+
embedded[:,mask_test] = test_sca
84+
delayed = pyLOM.math.time_delay_embedding(embedded, dimension=50)
85+
# Generate training validation and test datasets both for reconstruction of states
86+
output = shred(torch.from_numpy(delayed).permute(1,2,0).to(device)).cpu().detach().numpy().T
87+
# Load the POD scaler of that SHRED to transform POD data
88+
podscale = pyLOM.NN.MinMaxScaler(column=True)
89+
podscale = podscale.load(load_dict['podscale_path'])
90+
outres[:,:,kk] = podscale.inverse_transform(output)
91+
# Compute mean relative error
92+
MRE[:,kk] = pyLOM.math.MRE_array(full_pod, outres[:,:,kk])
93+
94+
95+
## Compute reconstruction statistics
96+
meanout = np.mean(outres, axis=2)
97+
stdout = np.std(outres, axis=2)
98+
meanMRE = np.mean(MRE, axis=1)
99+
stdMRE = np.std(MRE, axis=1)
100+
print('The mean MRE of the %i configurations is %.2f' % (nconfigs, np.mean(meanMRE)))
101+
102+
103+
## Save the mean output between configurations as the POD coefficients for reconstruction
104+
pyLOM.POD.save('POD_predicted_%s.h5'% podvar,None,None,meanout,sensors.partition_table,nvars=1)
105+
106+
107+
## Plot error bars
108+
#Non-scaled
109+
fig, _ = pyLOM.NN.plotModalErrorBars(meanMRE)
110+
fig.savefig('errorbars.pdf', dpi=300, bbox_inches='tight')
111+
#Scaled
112+
fig, _ = pyLOM.NN.plotModalErrorBars(Sscale*meanMRE)
113+
fig.savefig('errorbars_scaled.pdf', dpi=300, bbox_inches='tight')
114+
115+
## Plot POD modes reconstruction
116+
fig, _ = pyLOM.NN.plotTimeSeries(time, full_pod, meanout, stdout)
117+
fig.savefig('output_modes.pdf', dpi=600)
118+
119+
120+
pyLOM.cr_info()

Diff for: ‎Examples/NN/SHRED/example_04_reconstruct_cylinder.py

+51
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,51 @@
1+
#!/usr/bin/env python
2+
#
3+
# Example of POD.
4+
#
5+
# Last revision: 19/07/2021
6+
from __future__ import print_function, division
7+
8+
import mpi4py
9+
mpi4py.rc.recv_mprobe = False
10+
11+
import numpy as np
12+
import pyLOM
13+
14+
15+
## Parameters
16+
DATAFILE = '/gpfs/scratch/bsc21/bsc021828/DATA_PYLOM/CYLINDER.h5'
17+
VARIABLE = 'VELOX'
18+
19+
20+
## Data loading
21+
m = pyLOM.Mesh.load(DATAFILE)
22+
d = pyLOM.Dataset.load(DATAFILE,ptable=m.partition_table)
23+
X = d[VARIABLE]
24+
t = d.get_variable('time')
25+
26+
27+
## Load POD spatial modes and its energy computed with the training snapshots
28+
PSI, S = pyLOM.POD.load('POD_trai_%s.h5' % VARIABLE, vars=['U','S'])
29+
30+
31+
## Load the predicted POD coefficients from SHRED
32+
V = pyLOM.POD.load('POD_predicted_%s.h5' % VARIABLE, vars='V')[0]
33+
34+
35+
## Reconstruct the flow and readd the temporal mean (SHRED was trained without)
36+
X_POD = pyLOM.POD.reconstruct(PSI,S,V)
37+
mean = pyLOM.math.temporal_mean(X)
38+
X_PODm = pyLOM.math.subtract_mean(X_POD, -1*mean)
39+
40+
41+
## Compute RMSE with original data
42+
rmse = pyLOM.math.RMSE(X_PODm,X)
43+
pyLOM.pprint(0,'RMSE = %e'%rmse)
44+
45+
46+
## Dump reconstruction to ParaView
47+
d.add_field('VELOR',1,X_PODm)
48+
pyLOM.io.pv_writer(m,d,'flow',basedir='out/flow',instants=np.arange(t.shape[0],dtype=np.int32),times=t,vars=['VELOX','VELOR'],fmt='vtkh5')
49+
50+
51+
pyLOM.cr_info()

0 commit comments

Comments
 (0)