Select Git revision
ita_openGL2Matlab.m
Code owners
Assign users and groups as approvers for specific file changes. Learn more.
09_VisualizeReconstructed.py 5.48 KiB
import splinepy as sp
import numpy as np
import pandas as pd
import xml.etree.ElementTree as ET
from os import path
import shutil
import subprocess
import os
from splinepy.helpme.integrate import _get_integral_measure, _get_quadrature_information
from miscellaneous.DataPreparation import save_numpy_to_xml
from miscellaneous.error_evaluation import integrate, integrate_multipatch, compute_integral_error, get_solution_vectors, load_geometry, show_multipatch_field
from models.LinearRegression import LinearRegressionModel
from models.RBFInterpolator import RadialBasisRegressionModel
from models.GaussianProcessRegressor import GaussianProcessRegressionModel
VELOCITY_SIZE = 3066
PRESSURE_SIZE = 1399
DEGREE_ELEVATIONS = 1
common_show_options = {
"cmap": "jet",
"scalarbar": True,
"lighting": "off",
"control_points": False,
"knots": False,
}
path_to_executable = 'gismo/stokes_recreation'
models = {
'linregMod': LinearRegressionModel(),
# 'gaussianProcessMod_DotWhite': GaussianProcessRegressionModel(),
# 'gaussianProcessMod_RBF': GaussianProcessRegressionModel(),
# 'gaussianProcessMod_Matern': GaussianProcessRegressionModel(),
'gaussianProcessMod_RationalQuadratic': GaussianProcessRegressionModel(),
# 'gaussianProcessMod_ExpSineSquared': GaussianProcessRegressionModel(),
# 'radialBasisMod_Linear': RadialBasisRegressionModel(),
'radialBasisMod_thinplatespline': RadialBasisRegressionModel(),
'radialBasisMod_cubic': RadialBasisRegressionModel(),
'radialBasisMod_quintic': RadialBasisRegressionModel(),
# 'radialBasisMod_multiquadric': RadialBasisRegressionModel(),
# 'radialBasisMod_inversemultiquadric': RadialBasisRegressionModel(),
# 'radialBasisMod_inversequadratic': RadialBasisRegressionModel(),
# 'radialBasisMod_gaussian': RadialBasisRegressionModel(),
}
df = pd.read_excel(f'Data/test/parameter_input.xlsx').transpose()
SAMPLES_TEST = df.shape[1]
for model_name in models:
print(f'Model: \t\t\t{model_name}')
# Create temp folder
if not os.path.exists(f'temp_{model_name}'):
os.mkdir(f'temp_{model_name}')
# Read POD modes
data_training = np.load(f'models/trained_models/{model_name}/MetaData.npz')
# R = data_training['R']
R = [1, 11, 21]
for R_ in R:
print(f'Singular value: \t{R_}')
model = models[model_name]
model_velocity, model_pressure = model.load_model(f'{model_name}/velocity_{int(R_)}.pkl', f'{model_name}/pressure_{int(R_)}.pkl')
# Predict values
# Predicting test values
test_velocity_predict, test_pressure_predict = model.predict_test(model_velocity, model_pressure, rescale=True)
# for sample in range(SAMPLES_TEST):
# for sample in range(10):
for sample in range(5):
save_numpy_to_xml(test_velocity_predict[:,sample], f'temp_{model_name}/predict_velocity_field_{sample}.xml', rows=VELOCITY_SIZE, cols='1')
save_numpy_to_xml(test_pressure_predict[:,sample], f'temp_{model_name}/predict_pressure_field_{sample}.xml', rows=PRESSURE_SIZE, cols='1')
# Get all the input files for the evaluation
geometry_file = f'Data/test/microstructure/index_{sample}.xml'
velocity_solution = f'Data/test/velocities/velocity_field_{sample}.xml'
pressure_solution = f'Data/test/pressure/pressure_field_{sample}.xml'
velocity_predicted = f'temp_{model_name}/predict_velocity_field_{sample}.xml'
pressure_predicted = f'temp_{model_name}/predict_pressure_field_{sample}.xml'
print(f'Pressure ranges from {np.min(test_pressure_predict[:,sample])} to {np.max(test_pressure_predict[:,sample])}')
print(f'Velocity ranges from {np.min(test_velocity_predict[:,sample])} to {np.max(test_velocity_predict[:,sample])}')
# Create patches for the geometry with gismo script
program_output = subprocess.Popen(
[f'./{path_to_executable}', f'-f {geometry_file}', f'-v {velocity_solution}', f'-p {pressure_solution}', f'-w {velocity_predicted}', f'-q {pressure_predicted}', '-e 1', '--no-plot'],
stdout=subprocess.PIPE,
text=True
).communicate()[0]
# Recreate the geometry
PRESSURE_FILE = 'pressure_field_patches.xml'
VELOCITY_FILE = 'velocity_field_patches.xml'
PRESSURE_REC_FILE = 'pressure_field_rec_patches.xml'
VELOCITY_REC_FILE = 'velocity_field_rec_patches.xml'
GEOMETRY_FILE = geometry_file
pressure_data = get_solution_vectors(file_path=PRESSURE_FILE)
velocity_data = get_solution_vectors(file_path=VELOCITY_FILE, two_dimensional=True)
pressure_rec_data = get_solution_vectors(file_path=PRESSURE_REC_FILE)
velocity_rec_data = get_solution_vectors(file_path=VELOCITY_REC_FILE, two_dimensional=True)
# Load the geometry
microstructure, ms_vel = load_geometry(GEOMETRY_FILE, degree_elevations=DEGREE_ELEVATIONS)
# Show pressure and velocity field
mp_pressure = show_multipatch_field(microstructure, pressure_data, data_name="pressure")
mp_pressure.show()
# sp.io.default.export('test.svg', mp_pressure)
mp_velocity = show_multipatch_field(ms_vel, velocity_data, data_name="velocity")
mp_velocity.show()
shutil.rmtree(f'temp_{model_name}')
os.remove(PRESSURE_FILE)
os.remove(VELOCITY_FILE)
os.remove(PRESSURE_REC_FILE)
os.remove(VELOCITY_REC_FILE)