{ "cells": [ { "cell_type": "markdown", "id": "36dccb97", "metadata": {}, "source": [ "# Introduction to ROM methods: SVD based techniques\n", "\n", "This notebook presents the basis concepts of Reduced Order Modelling (ROM), starting from the most basic techniques based on the Singular Value Decomposition (SVD) for matrices. This notebook will also present the import class from OpenFOAM simulation, implemented within the *pyforce* package.\n", "\n", "In this tutorial, you will learn:\n", "\n", "- How to import snapshots from OpenFOAM simulations through the `ReadFromOF` class;\n", "- How to compress fluid dynamics fields (pressure and velocity) using randomised SVD and the Proper Orthogonal Decomposition (POD) technique.\n", "- How to create a simple surrogate model for the reduced dynamics using interpolation and regression techniques." ] }, { "cell_type": "markdown", "id": "a6d0a5a7", "metadata": {}, "source": [ "The data for this tutorial comes from an OpenFOAM tutorial (com version 2012) for *pimpleFoam* in laminar condition, called [`cylinder2D`](https://gitlab.com/openfoam/openfoam/-/tree/maintenance-v2012/tutorials/incompressible/pimpleFoam/laminar/cylinder2D?ref_type=heads)." ] }, { "cell_type": "markdown", "id": "31cf26d0", "metadata": {}, "source": [ "## Read OpenFOAM cases\n", "The *pyforce* package comes with the `ReadFromOF` class, which allows users to easily import simulation data from OpenFOAM cases. This class must be initialized with the path to the OpenFOAM case folder." ] }, { "cell_type": "code", "execution_count": 1, "id": "b4fa9220", "metadata": {}, "outputs": [], "source": [ "from pyforce.tools.write_read import ReadFromOF\n", "import warnings\n", "warnings.filterwarnings(\"ignore\", category=RuntimeWarning)\n", "\n", "io_data = ReadFromOF('Datasets/LaminarFlowCyl_OF2012/', \n", " skip_zero_time=True # optional\n", " )" ] }, { "cell_type": "markdown", "id": "948d5476", "metadata": {}, "source": [ "The variables to import are the pressure (scalar field) and the velocity (vector field). The snapshots can be loaded using `pyvista` (default) or `fluidfoam` or `foamlib` libraries: the latter is more efficient for large datasets, but it only supports cell data. The method `import_field` can be used to load the data from the OpenFOAM case and returns the snapshots as a `FunctionsList` object, including a vector with the time instances (from the OpenFOAM directories)." ] }, { "cell_type": "code", "execution_count": 2, "id": "c5dc54a4", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Importing p using pyvista\n", "Importing U using pyvista\n" ] } ], "source": [ "import numpy as np\n", "\n", "var_names = ['p', 'U']\n", "\n", "dataset = dict()\n", "\n", "for field in var_names:\n", "\n", " print('Importing '+field+f' using pyvista')\n", " dataset[field], times = io_data.import_field(field,\n", " import_mode='pyvista', verbose=False # optional\n", " )" ] }, { "cell_type": "markdown", "id": "0c980f23", "metadata": {}, "source": [ "Let us also extract the grid needed by *pyforce* using the `mesh` method, to make a plot of the snapshots" ] }, { "cell_type": "code", "execution_count": 3, "id": "b3d7e852", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "
| UnstructuredGrid | Information |
|---|---|
| N Cells | 34764 |
| N Points | 70396 |
| X Bounds | -6.075e-01, 1.120e+00 |
| Y Bounds | -6.000e-01, 6.000e-01 |
| Z Bounds | -7.500e-03, 7.500e-03 |
| N Arrays | 0 |