Commit a1e7656e authored by Gurecky, William's avatar Gurecky, William
Browse files

adds a new segmented pod rom helper class to handle joining per-zone pod models

parent 6a9a18f9
Loading
Loading
Loading
Loading
+132 −9
Original line number Diff line number Diff line
@@ -16,6 +16,7 @@ import matplotlib.pyplot as plt
import argparse
import time
import os
from typing import List
from scipy.interpolate import interp1d, RBFInterpolator
import h5py
try:
@@ -28,9 +29,9 @@ v_i_name = "Velocity[i] (m/s)"
v_j_name = "Velocity[j] (m/s)"
v_k_name = "Velocity[k] (m/s)"
p_name = "Static Pressure (Pa)"
v_i_name = "Velocity in Rotating[i] (m/s)"
v_j_name = "Velocity in Rotating[j] (m/s)"
v_k_name = "Velocity in Rotating[k] (m/s)"
#v_i_name = "Velocity in Rotating[i] (m/s)"
#v_j_name = "Velocity in Rotating[j] (m/s)"
#v_k_name = "Velocity in Rotating[k] (m/s)"


def construct_pod_model(snapshot_h5file, data_name="Static Pressure (Pa)", n_pod_modes=5):
@@ -91,6 +92,44 @@ class PodCfdModel(object):
        return self._xyz_points


class SegmentedPodCfdModel(object):
    """
    Contains several sub-POD ROM models that can be combined to form one output.
    Requires that each sub-POD ROM has the exact same inputs. Particularly useful when
    sharp boundaries in the flow or geometry exist in the full model: A different POD
    model can be fitted to the flow field in each region and later rejoined
    using this class.
    """
    def __init__(self, snapshot_h5files: List[str], data_name: str="Static Pressure (Pa)", n_pod_modes: int=7):
        self._pod_cfd_models = []
        for snapshot_h5file in snapshot_h5files:
            self._pod_cfd_models.append(PodCfdModel(snapshot_h5file, data_name, n_pod_modes))

    def predict(self, vane_angle: float, min_vane_angle: float=4., max_vane_angle: float=14.):
        """
        Evaluate the POD ROM.

        Args:
            vane_angle:  turbine guide vane angle in degrees
            min_vane_angle:  min vane angle (model only works
                             between min and max CFD vane angle data)
            min_vane_angle:  max vane angle (model only works
                             between min and max CFD vane angle data)
        """
        va = np.clip(vane_angle, min_vane_angle, max_vane_angle)
        seg_results = []
        for pod_model in self._pod_cfd_models:
            seg_results.append(pod_model.predict(va, min_vane_angle, max_vane_angle))
        return np.concatenate(seg_results, axis=0)

    @property
    def xyz_points(self):
        seg_xyz = []
        for pod_model in self._pod_cfd_models:
            seg_xyz.append(pod_model._xyz_points)
        return np.concatenate(seg_xyz, axis=0)


def test_pod_cfd_model(use_plotly=True,
                       snapshot_h5file="cfd_pod_snapshots.h5",
                       n_pod_modes=7):
@@ -155,6 +194,7 @@ def test_pod_cfd_model(use_plotly=True,
        print("Time to make plot: %0.4e (sec)" % (te-ts))

    else:
        # use matplotlib
        # create 3d plot figure of predicted pressure distribution
        ts = time.time()
        fig = plt.figure()
@@ -169,6 +209,83 @@ def test_pod_cfd_model(use_plotly=True,
        print("Time to make plot: %0.4e (sec)" % (te-ts))


def test_pod_segmented_cfd_model(snapshot_h5files, n_pod_modes=7):
    assert len(snapshot_h5files) > 1
    # create the reduced order model(s)
    pod_model = SegmentedPodCfdModel(
            snapshot_h5files, data_name=p_name,
            n_pod_modes=n_pod_modes)
    xyz_points = pod_model.xyz_points

    pod_model_vx = SegmentedPodCfdModel(
            snapshot_h5files, data_name=v_i_name, n_pod_modes=n_pod_modes)
    pod_model_vy = SegmentedPodCfdModel(
            snapshot_h5files, data_name=v_j_name, n_pod_modes=n_pod_modes)
    pod_model_vz = SegmentedPodCfdModel(
            snapshot_h5files, data_name=v_k_name, n_pod_modes=n_pod_modes)

    # evaluate the pod model at an vane angle not in training set
    ts = time.time()
    current_vane_angle = 12.5
    pred_p = pod_model.predict(current_vane_angle)
    pred_vx = pod_model_vx.predict(current_vane_angle)
    pred_vy = pod_model_vy.predict(current_vane_angle)
    pred_vz = pod_model_vz.predict(current_vane_angle)
    pred_vmag = np.sqrt(pred_vx**2 + pred_vy**2 + pred_vz**2)
    te = time.time()
    print("Time to eval POD ROM: %0.4e (sec)" % (te-ts))

    X = np.asarray(xyz_points[:, 0], dtype=np.float32)
    Y = np.asarray(xyz_points[:, 1], dtype=np.float32)
    Z = np.asarray(xyz_points[:, 2], dtype=np.float32)
    values = np.asarray(pred_p, dtype=np.float32)

    # make plots
    import plotly.graph_objects as go
    import plotly
    ts = time.time()
    fig_p = go.Figure(data=go.Scatter3d(
        x=X.flatten(),
        y=Y.flatten(),
        z=Z.flatten(),
        mode='markers',
        marker=dict(size=4, color=values, opacity=0.95,
                    colorbar=dict(thickness=5, title="Static Pressure (Pa)")),
        customdata=values,
        hovertemplate=('%{customdata}%'),
        ))
    fig_p.show()

    fig_v = go.Figure(data=go.Cone(
        x=X.flatten(),
        y=Y.flatten(),
        z=Z.flatten(),
        u=pred_vx,
        v=pred_vy,
        w=pred_vz,
        sizemode="absolute",
        sizeref=150,
        colorbar=dict(thickness=8, title="Velocity (m/s)"),
        ))
    fig_v.show()
#     fig_vs = go.Figure(data=go.Streamtube(
#         x=X.flatten(),
#         y=Y.flatten(),
#         z=Z.flatten(),
#         u=pred_vx,
#         v=pred_vy,
#         w=pred_vz,
#         starts= dict(
#             x=[3.32,],
#             y=[1.0,],
#             z=[-1.75,]
#             ),
#         ))
#     fig_vs.show()
    te = time.time()
    print("Time to make plot: %0.4e (sec)" % (te-ts))


def test_dynamic_pod_model(snapshot_h5file="cfd_pod_snapshots.h5", n_pod_modes=7):
    # create the reduced order model(s)
    pod_model = PodCfdModel(
@@ -258,14 +375,20 @@ def test_dynamic_pod_model(snapshot_h5file="cfd_pod_snapshots.h5", n_pod_modes=7

if __name__ == "__main__":
    parser = argparse.ArgumentParser()
    parser.add_argument("-i", help="Input POD snapshot file.", type=str, default="cfd_pod_snapshots.h5")
    parser.add_argument("-i", help="Input POD snapshot file(s).", type=str, nargs='+', default=["cfd_pod_snapshots.h5"])
    parser.add_argument("--dynamic", help="Run in dynamic test mode with moving vane angle", type=int, default=0)
    parser.add_argument("--n_pod_modes", help="Run in dynamic test mode with moving vane angle", type=int, default=7)
    args = parser.parse_args()
    if len(args.i) == 1:
        if bool(args.dynamic):
            test_dynamic_pod_model(
                args.i, n_pod_modes=args.n_pod_modes)
                    snapshot_h5file=args.i[0],
                    n_pod_modes=args.n_pod_modes)
        else:
            test_pod_cfd_model(
                use_plotly=True, snapshot_h5file=args.i,
                    use_plotly=True, snapshot_h5file=args.i[0],
                    n_pod_modes=args.n_pod_modes)
    else:
        test_pod_segmented_cfd_model(
                snapshot_h5files=args.i,
                n_pod_modes=args.n_pod_modes)
+66 −37
Original line number Diff line number Diff line
@@ -37,46 +37,65 @@ run_parameters = [
        ]

# NEW CFD DATA FILES
run_dir = "./alder_u12/"
# run_dir = "./alder_u12_fixed/"
run_dir = "./alder_u12_fixed/"
run_csv_files = [
        "alpha-01_runner.csv",
        "alpha-02_runner.csv",
        "alpha-03_runner.csv",
        "alpha-04_runner.csv",
        "alpha-05_runner.csv",
        "alpha-06_runner.csv",
        "alpha-07_runner.csv",
        "alpha-08_runner.csv",
        "alpha-09_runner.csv",
#         "alpha-01_runner.csv",
#         "alpha-02_runner.csv",
#         "alpha-03_runner.csv",
#         "alpha-04_runner.csv",
#         "alpha-05_runner.csv",
#         "alpha-06_runner.csv",
#         "alpha-07_runner.csv",
#         "alpha-08_runner.csv",
#         "alpha-09_runner.csv",
#         "alpha-10_runner.csv",
        "alpha-9_runner.csv",
        "alpha-10_runner.csv",
        ]
run_csv_files = [
        "alpha-01_spiral-casing-and-draft-tube.csv",
        "alpha-02_spiral-casing-and-draft-tube.csv",
        "alpha-03_spiral-casing-and-draft-tube.csv",
        "alpha-04_spiral-casing-and-draft-tube.csv",
        "alpha-05_spiral-casing-and-draft-tube.csv",
        "alpha-06_spiral-casing-and-draft-tube.csv",
        "alpha-07_spiral-casing-and-draft-tube.csv",
        "alpha-08_spiral-casing-and-draft-tube.csv",
        "alpha-09_spiral-casing-and-draft-tube.csv",
        "alpha-10_spiral-casing-and-draft-tube.csv",
        ]
# run_csv_files = [
# #         "alpha-01_spiral-casing-and-draft-tube.csv",
# #         "alpha-02_spiral-casing-and-draft-tube.csv",
# #         "alpha-03_spiral-casing-and-draft-tube.csv",
# #         "alpha-04_spiral-casing-and-draft-tube.csv",
# #         "alpha-05_spiral-casing-and-draft-tube.csv",
# #         "alpha-06_spiral-casing-and-draft-tube.csv",
# #         "alpha-07_spiral-casing-and-draft-tube.csv",
# #         "alpha-08_spiral-casing-and-draft-tube.csv",
# #         "alpha-09_spiral-casing-and-draft-tube.csv",
# #         "alpha-10_spiral-casing-and-draft-tube.csv",
#           "alpha-9_casing-draft_tube.csv",
#           "alpha-10_casing-draft_tube.csv",
#          ]
# run_csv_files = [
#         # "alpha-01_distributor.csv",
#         # "alpha-02_distributor.csv",
#         # "alpha-03_distributor.csv",
#         # "alpha-04_distributor.csv",
#         # "alpha-05_distributor.csv",
#         # "alpha-06_distributor.csv",
#         # "alpha-07_distributor.csv",
#         # "alpha-08_distributor.csv",
#         # "alpha-09_distributor.csv",
#         # "alpha-10_distributor.csv",
#         "alpha-9_distributor.csv",
#         "alpha-10_distributor.csv",
#         ]
run_parameters = [
        [1.3875],
        [2.775],
        [4.1625],
        [5.55],
        [6.9375],
        [8.325],
        [9.7125],
        [11.1],
#         [1.3875],
#         [2.775],
#         [4.1625],
#         [5.55],
#         [6.9375],
#         [8.325],
#         [9.7125],
#         [11.1],
        [12.4875],
        [13.875],
        ]


def plot_cfd_vel_field(mesh_x, mesh_y, mesh_z, v_i, v_j, v_k):
def plot_cfd_vel_field(mesh_x, mesh_y, mesh_z, v_i, v_j, v_k, cbar_label=""):
    import plotly.graph_objects as go
    import plotly

@@ -96,8 +115,8 @@ def plot_cfd_vel_field(mesh_x, mesh_y, mesh_z, v_i, v_j, v_k):
        v=v_j[subsample_idxs],
        w=v_k[subsample_idxs],
        sizemode="absolute",
        sizeref=150,
        colorbar=dict(thickness=8, title="CFD Velocity (m/s)"),
        sizeref=25,
        colorbar=dict(thickness=8, title=cbar_label+"CFD Velocity (m/s)"),
        ))
    fig_v.show()
    fig_p = go.Figure(data=go.Scatter3d(
@@ -106,7 +125,7 @@ def plot_cfd_vel_field(mesh_x, mesh_y, mesh_z, v_i, v_j, v_k):
        z=mesh_z[subsample_idxs].flatten(),
        mode='markers',
        marker=dict(size=4, color=v_mag[subsample_idxs], opacity=0.95,
                    colorbar=dict(thickness=5, title="CFD Vel Magnitude (m/s)")),
                    colorbar=dict(thickness=5, title=cbar_label+"CFD Vel Magnitude (m/s)")),
        customdata=v_mag[subsample_idxs],
        hovertemplate=('%{customdata}%'),
        ))
@@ -168,8 +187,12 @@ def conv2hdf(csv_fname, outfile_name="fixed_cfd.h5", case_params=[], dset_names=
        v_i = hdf["/Velocity[i] (m/s)"][:]
        v_j = hdf["/Velocity[j] (m/s)"][:]
        v_k = hdf["/Velocity[k] (m/s)"][:]
        vr_i = hdf["Velocity in Rotating[i] (m/s)"][:]
        vr_j = hdf["Velocity in Rotating[j] (m/s)"][:]
        vr_k = hdf["Velocity in Rotating[k] (m/s)"][:]
        p_s = hdf["/Static Pressure (Pa)"][:]
        plot_cfd_vel_field(xx, yy, zz, v_i, v_j, v_k)
        plot_cfd_vel_field(xx, yy, zz, vr_i, vr_j, vr_k, cbar_label="Rotating Ref ")
        plot_cfd_pres_field(xx, yy, zz, p_s)

    # clean up
@@ -179,7 +202,8 @@ def conv2hdf(csv_fname, outfile_name="fixed_cfd.h5", case_params=[], dset_names=
def reinterpolate_cfd(
        h5_fname_base, h5_fnames=[], n_subsample=100000,
        field_names=["Pressure (Pa)"],
        out_fname="cfd_pod_snapshots.h5"):
        out_fname="cfd_pod_snapshots.h5",
        vel_mult=1.0):
    """
    Re-interpolates all data to be on a consistent mesh for POD
    modeling.
@@ -248,6 +272,7 @@ def reinterpolate_cfd(

            # Summary data
            if "Velocity" in field_name:
                case_sub_field_data = vel_mult * case_sub_field_data
                summary_field_data = np.abs(case_sub_field_data)
            else:
                summary_field_data = case_sub_field_data
@@ -287,6 +312,7 @@ def reinterpolate_cfd(
if __name__ == "__main__":
    parser = argparse.ArgumentParser()
    parser.add_argument("-o", help="output h5 file name", type=str, default='cfd_pod_snapshots.h5')
    parser.add_argument("--vel_mult", help="Velocity vector multiplier", type=float, default=1.0)
    args = parser.parse_args()
    h5_fnames = []
#     field_names = [
@@ -311,10 +337,13 @@ if __name__ == "__main__":
        print("Preprocessing file: %s" % csv_file)
        out_file = csv_file + ".h5"
        make_plot = False
        if i == 7:
        if i == 1:
            make_plot = True
        conv2hdf(run_dir + csv_file, outfile_name=out_file, case_params=params, make_plot=make_plot,
                 csv_col_names=csv_col_names)
        h5_fnames.append(out_file)
        i += 1
    reinterpolate_cfd(h5_fnames[0], h5_fnames=h5_fnames, field_names=field_names, out_fname=args.o)
    base_id = max(0, int(len(h5_fnames) - len(h5_fnames) / 2))
    reinterpolate_cfd(h5_fnames[base_id], h5_fnames=h5_fnames,
                      field_names=field_names, out_fname=args.o,
                      vel_mult=args.vel_mult)