diff --git a/.github/workflows/build_documentation.yml b/.github/workflows/build_documentation.yml index 441e56ac..1fc28dbf 100644 --- a/.github/workflows/build_documentation.yml +++ b/.github/workflows/build_documentation.yml @@ -16,11 +16,24 @@ permissions: jobs: build-documentation: runs-on: ubuntu-latest + container: + image: ghcr.io/4c-multiphysics/4c-minimal:main + options: --user root --env OMPI_ALLOW_RUN_AS_ROOT=1 --env OMPI_ALLOW_RUN_AS_ROOT_CONFIRM=1 defaults: run: shell: bash -l {0} steps: - uses: actions/checkout@v4 + - name: Mark repo as safe for git + run: | + git config --global --add safe.directory "$GITHUB_WORKSPACE" + - name: Install rsync + run: | + sudo apt-get update + sudo apt-get install -y rsync + - name: Create links to 4C + run: | + ln -s /home/user/4C/bin/ config/4C_build - name: Create Python environment id: environment uses: ./.github/actions/create_python_environment @@ -31,14 +44,19 @@ jobs: $PYTHON_PACKAGE_MANAGER activate queens pip install -e .[tutorial] python -m ipykernel install --user --name queens --display-name "Python (queens)" + - name: Install xvfb + run: | + apt-get update + apt-get install -y xvfb - name: Sphinx build env: PYTHON_PACKAGE_MANAGER: ${{steps.environment.outputs.ppm}} run: | + set -euxo pipefail $PYTHON_PACKAGE_MANAGER activate queens sphinx-apidoc -o doc/source src/ -fMT cd doc - sphinx-build -b html -d build/doctrees source build/html -W + xvfb-run -a sphinx-build -b html -d build/doctrees source build/html -W - name: Upload html uses: actions/upload-pages-artifact@v3 with: diff --git a/doc/source/_ext/create_documentation_files.py b/doc/source/_ext/create_documentation_files.py index 4e87bd97..667e8282 100644 --- a/doc/source/_ext/create_documentation_files.py +++ b/doc/source/_ext/create_documentation_files.py @@ -286,14 +286,22 @@ def download_file_from_url(url, file_name): def copy_tutorials(): - """Copy tutorials from source to doc.""" - for tutorial in relative_path_from_root("tutorials").glob("*.ipynb"): - destination = relative_to_doc_source("tutorials/" + tutorial.name) + """Copy tutorials and other util and input files from source to doc.""" + source_root = relative_path_from_root("tutorials") + destination_root = relative_to_doc_source("tutorials") + allowed_patterns = ("*.ipynb", "*.py", "*.exo", "*.yaml") - if destination.exists(): - destination.unlink() + for pattern in allowed_patterns: + for source in source_root.rglob(pattern): + rel = source.relative_to(source_root) + destination = destination_root / rel - shutil.copyfile(tutorial, destination) + destination.parent.mkdir(parents=True, exist_ok=True) + + if destination.exists(): + destination.unlink() + + shutil.copyfile(source, destination) def main(): diff --git a/tutorial-requirements.in b/tutorial-requirements.in index cd88f906..91de3b7a 100644 --- a/tutorial-requirements.in +++ b/tutorial-requirements.in @@ -4,3 +4,4 @@ # Tutorials scikit-fem +pyvista diff --git a/tutorial-requirements.txt b/tutorial-requirements.txt index 87afb3d4..8b88b093 100644 --- a/tutorial-requirements.txt +++ b/tutorial-requirements.txt @@ -4,14 +4,104 @@ # # pip-compile --constraint=requirements.txt --output-file=tutorial-requirements.txt tutorial-requirements.in # +certifi==2024.8.30 + # via + # -c requirements.txt + # requests +charset-normalizer==3.4.0 + # via + # -c requirements.txt + # requests +contourpy==1.3.0 + # via + # -c requirements.txt + # matplotlib +cycler==0.12.1 + # via + # -c requirements.txt + # matplotlib +fonttools==4.54.1 + # via + # -c requirements.txt + # matplotlib +idna==3.10 + # via + # -c requirements.txt + # requests +kiwisolver==1.4.7 + # via + # -c requirements.txt + # matplotlib +matplotlib==3.9.2 + # via + # -c requirements.txt + # pyvista + # vtk numpy==1.26.4 # via # -c requirements.txt + # contourpy + # matplotlib + # pyvista # scikit-fem # scipy +packaging==24.1 + # via + # -c requirements.txt + # matplotlib + # pooch +pillow==11.0.0 + # via + # -c requirements.txt + # matplotlib + # pyvista +platformdirs==4.3.6 + # via + # -c requirements.txt + # pooch +pooch==1.8.2 + # via + # -c requirements.txt + # pyvista +pyparsing==3.2.0 + # via + # -c requirements.txt + # matplotlib +python-dateutil==2.9.0.post0 + # via + # -c requirements.txt + # matplotlib +pyvista==0.44.1 + # via + # -c requirements.txt + # -r tutorial-requirements.in +requests==2.32.3 + # via + # -c requirements.txt + # pooch scikit-fem==11.0.0 - # via -r tuto-requirements.in + # via -r tutorial-requirements.in scipy==1.14.1 # via # -c requirements.txt # scikit-fem +scooby==0.10.0 + # via + # -c requirements.txt + # pyvista +six==1.16.0 + # via + # -c requirements.txt + # python-dateutil +typing-extensions==4.12.2 + # via + # -c requirements.txt + # pyvista +urllib3==2.2.3 + # via + # -c requirements.txt + # requests +vtk==9.3.1 + # via + # -c requirements.txt + # pyvista diff --git a/tutorials/3_orchestrating_4c_simulations/3_orchestrating_4c_simulations.ipynb b/tutorials/3_orchestrating_4c_simulations/3_orchestrating_4c_simulations.ipynb new file mode 100644 index 00000000..3ae6dd05 --- /dev/null +++ b/tutorials/3_orchestrating_4c_simulations/3_orchestrating_4c_simulations.ipynb @@ -0,0 +1,906 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# 3. Orchestrating 4C Simulations\n", + "\n", + "Until now, all the simulation models we have used in QUEENS were Python-based. For a lot of us, this is not the case. From [fenics](https://fenicsproject.org/) to [abaqus](https://www.3ds.com/de/products/simulia/abaqus), we all have our own preferred solvers. One main advantage of QUEENS is its independence of solver type. This means, to use your QUEENS model with your favorite solver, the only thing you have to do is write an interface.\n", + "\n", + "Since QUEENS and the multiphysics solver [4C](https://github.com/4C-multiphysics/4C) share some developers, the QUEENS community already provides a ready-to-go interface. So let's start there." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## The 4C model\n", + "\n", + "In this example, we'll have a look at a nonlinear solid mechanics problem. The momentum equation is given by\n", + "\n", + "$$\\rho \\ddot{\\boldsymbol{u}} = \\nabla \\cdot \\boldsymbol{\\sigma} + \\boldsymbol{b} \\text{ in } \\Omega \\subset \\mathbb{R}^3$$\n", + "\n", + "where $\\boldsymbol{u}$ are the displacements, $\\boldsymbol{\\sigma}$ the Cauchy stresses, $\\boldsymbol{b}$ the volumetric body forces, and $\\Omega$ the domain of the solid. For the following examples, $\\boldsymbol{b}=\\boldsymbol{0}$, and we prescribe Dirichlet boundary conditions\n", + "$$ \\boldsymbol{u} = \\boldsymbol{0} \\text{ on } \\Gamma_w = \\{z=-2.5\\}$$\n", + "$$ \\boldsymbol{u} = [0, 0, 0.55]^T \\text{ on } \\Gamma_e =\\{z=2.5\\}$$\n", + "\n", + "and a Neumann boundary condition:\n", + "$$ \\boldsymbol{\\sigma} \\cdot \\boldsymbol{n} = \\boldsymbol{0} \\text{ on } \\Gamma \\setminus \\{\\Gamma_w \\cup \\Gamma_e\\}$$\n", + "\n", + "We assume a static problem, so $\\ddot{\\boldsymbol{u}} = 0$, and use a Saint Venant–Kirchhoff material model\n", + "$$ \\boldsymbol{S} = \\lambda \\text{tr}(\\boldsymbol{E})\\boldsymbol{I} + 2 \\mu \\boldsymbol{E}$$\n", + "where $\\boldsymbol{S}$ is the second Piola-Kirchhoff stress tensor and $\\boldsymbol{E}$ the Green-Lagrange strain tensor with the Lamé parameters\n", + "$$\\lambda = \\frac{E\\nu}{(1+\\nu)(1-2\\nu)}$$\n", + "$$\\mu = \\frac{E}{2(1+\\nu)}$$\n", + "where $E$ is the Young's modulus and $\\nu$ the Poisson's ratio.\n", + "\n", + "The problem is discretized in space using 40 linear finite elements.\n", + "\n", + "> Note: This tutorial is based on the 4C tutorial [Preconditioners for Iterative Solvers](https://4c-multiphysics.github.io/4C/documentation/tutorials/tutorial_preconditioning.html)\n", + "\n", + "> Note: In the following, we use PyVista for plotting; however, this might not work for some setups. In that case, set `pyvista_available` to `False` in the next cell, and please refer to the [online tutorial](https://queens-py.github.io/queens/tutorials/4_orchestrating_4C_simulations_with_QUEENS.html) in the documentation for figures." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Set this to False if pyvista is causing the kernel to crash\n", + "pyvista_available = True" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Define some paths\n", + "from pathlib import Path\n", + "from utils import find_repo_root\n", + "\n", + "home = Path.home()\n", + "NOTEBOOK_DIR = Path.cwd()\n", + "\n", + "QUEENS_BASE_DIR = find_repo_root(NOTEBOOK_DIR)\n", + "QUEENS_EXPERIMENTS_DIR = home / \"queens-experiments\"\n", + "TUTORIAL_DIR = QUEENS_BASE_DIR / \"tutorials\" / \"3_orchestrating_4c_simulations\"" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Let's look at the mesh!\n", + "import os\n", + "\n", + "if pyvista_available:\n", + " os.environ[\"PYVISTA_OFF_SCREEN\"] = \"true\"\n", + " os.environ[\"VTK_DEFAULT_RENDER_WINDOW_OFFSCREEN\"] = \"1\"\n", + " os.environ[\"PYVISTA_JUPYTER_BACKEND\"] = \"static\"\n", + "\n", + " import pyvista as pv\n", + "\n", + " fe_mesh = pv.read(TUTORIAL_DIR / \"beam_coarse.exo\")\n", + "\n", + " pl = pv.Plotter(off_screen=True, title=\"Finite Element mesh\", window_size=(1400, 500))\n", + " pl.add_mesh(fe_mesh, show_edges=True, opacity=0.8)\n", + " pl.add_mesh(pv.Plane(center=[0, 0, -2.50001], direction=[0, 0, 1]), color=\"blue\") # Gamma_w\n", + " pl.add_mesh(pv.Plane(center=[0, 0, 2.50001], direction=[0, 0, 1]), color=\"red\") # Gamma_e\n", + " pl.add_axes(line_width=5)\n", + " pl.camera_position = [\n", + " (1.1321899035097993, -6.851600196807601, 2.7649096132703574),\n", + " (0.0, 0.0, 0.275),\n", + " (-0.97637930372977, -0.08995062285804697, 0.19644933366041056),\n", + " ]\n", + "\n", + " pl.show()\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The red plane is $\\Gamma_e$, the blue one is $\\Gamma_w$." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Simulation Analytics\n", + "\n", + "Now, we want to have a look at the Cauchy stresses $\\boldsymbol{\\sigma}$ for various Young's moduli $E$ and Poisson's ratios $\\nu$. For the purpose of demonstration, let's assume a uniform distribution for each material parameter.\n", + "\n", + "Of course we use QUEENS, so let's get started!" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from queens.distributions import Uniform\n", + "from queens.parameters import Parameters\n", + "\n", + "# Let's define the distributions\n", + "young = Uniform(10, 10000)\n", + "poisson_ratio = Uniform(0, 0.5)\n", + "\n", + "# Let's define the parameters\n", + "parameters = Parameters(young=young, poisson_ratio=poisson_ratio)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now that we have the parameters ready, we need a `Driver`. In QUEENS terms, drivers call the underlying simulation models, provide the simulation input and handle logging and error management. They are the interfaces to arbitrary simulation software and essentially do everything that is needed to start a **single** simulation.\n", + "\n", + "As mentioned before, for 4C, these are already provided by the QUEENS community.\n", + "\n", + "> Note: This tutorial assumes a working local 4C installation. Please create a symbolic link to your 4C build directory and store it under `/config/4C_build` via this command:\n", + "> ```\n", + "> ln -s /config/4C_build\n", + "> ```" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from queens_interfaces.fourc.driver import Fourc\n", + "\n", + "fourc_executable = QUEENS_BASE_DIR / \"config\" / \"4C_build\" / \"4C\"\n", + "input_template = TUTORIAL_DIR / \"input_template.4C.yaml\"\n", + "input_mesh = TUTORIAL_DIR / \"beam_coarse.exo\"\n", + "\n", + "fourc_driver = Fourc(\n", + " parameters=parameters, # parameters to run\n", + " input_templates=input_template, # input template\n", + " executable=fourc_executable, # 4C executable\n", + " files_to_copy=[input_mesh], # copy mesh file to experiment directory\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "4C simulations are controlled via yaml input files. They need to be created for every simulation. Therefore, we provide a [jinja2](https://jinja.palletsprojects.com/en/stable/) input **template** to the driver. Jinja2 is a fast, expressive, and extensible templating engine. It starts from a template, for example\n", + "\n", + "```yaml\n", + "MATERIALS:\n", + " - MAT: 1\n", + " MAT_Struct_StVenantKirchhoff:\n", + " YOUNG: {{ young }}\n", + " NUE: {{ poisson_ratio }}\n", + " DENS: 1\n", + "```\n", + "\n", + "and inserts the desired values (here: `young = 1000` and `poisson_ratio = 0.3`) into the respective placeholders (here: `{{ young }}` and `{{ poisson_ratio }}`). This results in:\n", + "```yaml\n", + "MATERIALS:\n", + " - MAT: 1\n", + " MAT_Struct_StVenantKirchhoff:\n", + " YOUNG: 1000\n", + " NUE: 0.3\n", + " DENS: 1\n", + "```\n", + "\n", + "For 4C, this allows us to create inputs for simulation runs without the need for a 4C API! Let's run this example!" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "\n", + "driver_output = fourc_driver.run(\n", + " sample=np.array([1000, 0.3]),\n", + " job_id=\"fourc_young_1000_poisson_ratio_03\",\n", + " num_procs=1,\n", + " experiment_dir=Path.cwd(),\n", + " experiment_name=\"fourc_run\",\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Some information:\n", + "\n", + "- `sample`: is the inputs for which we want to evaluate the model. The order of the data is given by the order in which we define the parameters.\n", + "- `job_id`: is the identifier for the 4C run and also creates a folder with the same name where the data is stored.\n", + "- `num_procs`: 4C is called using MPI. This is the number of processes that 4C will be called with.\n", + "- `experiment_dir`: is the directory where we find all necessary data for 4C (input template or any other required files).\n", + "- `experiment_name`: is an identifier of a QUEENS run. This will come in handy when we evaluate 4C multiple times." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's look at the output folder structure of the QUEENS run:\n", + "\n", + "```bash\n", + "fourc_young_1000_poisson_ratio_03/\n", + "├── input_file.yaml\n", + "├── jobscript.sh\n", + "├── metadata.yaml\n", + "└── output\n", + " ├── output.control\n", + " ├── output.log\n", + " ├── output.mesh.structure.s0\n", + " ├── output-structure.pvd\n", + " └── output-vtk-files\n", + " ├── structure-00000-0.vtu\n", + " ├── structure-00000.pvtu\n", + " ├── structure-00001-0.vtu\n", + " └── structure-00001.pvtu\n", + "```\n", + "\n", + "What are those files?\n", + "\n", + "- `input_file.yaml`: is the 4C input file with $E=1000$ and $\\nu=0.3$\n", + "- `jobscript.sh`: contains the command with which 4C was called\n", + "- `metadata.yaml`: contains metadata of the 4C run (timings, inputs, outputs, ...)\n", + "- `output`-folder: all the files in here are output files of 4C! The sole exception is `output.log`, to which QUEENS redirected the logging data\n", + "\n", + "Go have a look at the `yaml`, `sh` and log files." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now, let us look at the Cauchy stresses $\\boldsymbol{\\sigma}$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# This is just needed for some nice plots\n", + "\n", + "from utils import plot_results\n", + "if pyvista_available:\n", + " plotter = pv.Plotter()\n", + " plot_results(\n", + " \"fourc_young_1000_poisson_ratio_03/output/output-vtk-files/structure-00001.pvtu\",\n", + " plotter,\n", + " )\n", + " plotter.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We plotted the deformed mesh. The original mesh is displayed as the blue 'wireframe'. Let us look into the driver output in Python." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "print(\"Driver output:\", driver_output)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "As we can see, it's `{}`. But why is it `{}`?\n", + "\n", + "Well, we told the driver how to run a 4C simulation, but never told it which information we want to extract. For this we need a data processor, i.e., an object that is able to extract data from the simulation results." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from queens.data_processors import PvdFile\n", + "\n", + "\n", + "# Create custom data processor from the available PVD file reader\n", + "class CustomDataProcessor(PvdFile):\n", + " def __init__(\n", + " self,\n", + " field_name,\n", + " failure_value,\n", + " file_name_identifier=None,\n", + " files_to_be_deleted_regex_lst=None,\n", + " ):\n", + " super().__init__(\n", + " field_name,\n", + " file_name_identifier,\n", + " file_options_dict={},\n", + " files_to_be_deleted_regex_lst=files_to_be_deleted_regex_lst,\n", + " point_data=False, # We use cell data\n", + " )\n", + " self.failure_value = failure_value\n", + "\n", + " def get_data_from_file(self, base_dir_file):\n", + "\n", + " # Get the data\n", + " data = super().get_data_from_file(base_dir_file)\n", + "\n", + " # Got some data\n", + " if data is not None:\n", + " data = data.flatten()\n", + " # Simulation started but only the first step was exported\n", + " if np.abs(data).sum() < 1e-8:\n", + " return self.failure_value\n", + "\n", + " # Return result\n", + " return data\n", + "\n", + " # Simulation could not even start\n", + " else:\n", + " return self.failure_value\n", + "\n", + "\n", + "custom_data_processor = CustomDataProcessor(\n", + " field_name=\"element_cauchy_stresses_xyz\", # Name of the field in the pvtu file\n", + " file_name_identifier=\"output-structure.pvd\", # Name of the file in the output folder\n", + " failure_value=np.nan * np.ones(40 * 6), # NaN for all elements\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We have defined our data processor, so let's continue." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Create a new driver\n", + "fourc_driver = Fourc(\n", + " parameters=parameters,\n", + " input_templates=input_template,\n", + " executable=fourc_executable,\n", + " data_processor=custom_data_processor, # tell the driver to extract the data\n", + " files_to_copy=[input_mesh], # copy mesh file to experiment directory\n", + ")\n", + "\n", + "# Run the driver again\n", + "driver_output = fourc_driver.run(\n", + " sample=np.array([1000, 0.3]),\n", + " job_id=\"fourc_young_1000_poisson_ratio_03\",\n", + " num_procs=1,\n", + " experiment_dir=Path.cwd(),\n", + " experiment_name=\"fourc_run\",\n", + ")\n", + "\n", + "print(\"Driver output\", driver_output)\n", + "print(\"Key(s) in driver output:\", list(driver_output.keys()))\n", + "print(\n", + " \"Number of entries in the \\'result\\' key:\",\n", + " len(driver_output[\"result\"]),\n", + " \"i.e., 40 elements x 6 Cauchy stress values\",\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The driver output can contain two keys: 'result' and 'gradient'. The first one is supposed to be the output of the model, and the second one the model gradient.\n", + "\n", + "> Note: The gradient value is only available for certain models; this particular 4C model does not provide it.\n", + "\n", + "What we did:\n", + "\n", + "- We defined a driver using QUEENS to run a 4C simulation via a Python function call\n", + "- We provided a custom way of extracting simulation results\n", + "- Although we selected a single process, setting `num_procs>1` allows us to run a simulation in parallel without having to change anything in the interface!\n", + "\n", + "**The complexity of working with a sophisticated simulation code (without an API) is hidden from the user!**\n", + "\n", + "**Instead, it is reduced to** `fourc_driver.run()`." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now in the day-to-day of computational mechanics research, we often want to study convergence of a model depending on certain parameters. In this example, we'll have a look at parameter combinations of $E \\in [10,10000]$ and $\\nu \\in [0,0.5]$ which lead to a failing simulation. Let's evaluate the 4C model on a grid of $4 \\times 21$ points." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Make it less verbose\n", + "import logging\n", + "\n", + "logging.getLogger(\"distributed\").setLevel(logging.WARN)\n", + "\n", + "from queens.schedulers import Local\n", + "from queens.iterators import Grid\n", + "from queens.models import Simulation\n", + "from queens.global_settings import GlobalSettings\n", + "from queens.main import run_iterator\n", + "from queens.utils.io import load_result\n", + "\n", + "n_grid_points_young = 4\n", + "n_grid_points_poisson_ratio = 21\n", + "\n", + "grid_output = None\n", + "if __name__ == \"__main__\":\n", + " Path(\"grid\").mkdir(exist_ok=True)\n", + " with GlobalSettings(\n", + " experiment_name=f\"grid_{n_grid_points_young}x{n_grid_points_poisson_ratio}\",\n", + " output_dir=\"grid\",\n", + " ) as gs:\n", + " scheduler = Local(\n", + " gs.experiment_name,\n", + " num_jobs=4, # 4 simulations in parallel\n", + " num_procs=1, # each simulation uses 1 process\n", + " )\n", + " model = Simulation(scheduler, fourc_driver)\n", + "\n", + " iterator = Grid(\n", + " model=model,\n", + " parameters=parameters,\n", + " global_settings=gs,\n", + " result_description={\"write_results\": True},\n", + " grid_design={\n", + " \"young\": {\n", + " \"num_grid_points\": n_grid_points_young,\n", + " \"axis_type\": \"log10\",\n", + " \"data_type\": \"FLOAT\",\n", + " },\n", + " \"poisson_ratio\": {\n", + " \"num_grid_points\": n_grid_points_poisson_ratio,\n", + " \"axis_type\": \"lin\",\n", + " \"data_type\": \"FLOAT\",\n", + " },\n", + " },\n", + " )\n", + "\n", + " run_iterator(iterator, gs)\n", + "\n", + " grid_output = load_result(gs.result_file(\"pickle\"))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Let's plot the data\n", + "grid_points = grid_output[\"input_data\"]\n", + "cauchy_stresses = grid_output[\"raw_output_data\"][\"result\"]\n", + "\n", + "import matplotlib.pyplot as plt\n", + "\n", + "converged_samples = []\n", + "failed_samples = []\n", + "failed = []\n", + "for input_sample, cauchy_output in zip(grid_points, cauchy_stresses, strict=True):\n", + " if np.isnan(cauchy_output).any():\n", + " failed_samples.append(input_sample)\n", + " failed.append(1)\n", + " else:\n", + " converged_samples.append(input_sample)\n", + " failed.append(0)\n", + "\n", + "converged_samples = np.array(converged_samples)\n", + "failed_samples = np.array(failed_samples)\n", + "failed = np.array(failed)\n", + "\n", + "fig, ax = plt.subplots(figsize=(12, 8))\n", + "ax.contourf(\n", + " grid_points[:, 0].reshape(21, 4),\n", + " grid_points[:, 1].reshape(21, 4),\n", + " failed.reshape(21, 4),\n", + " levels=[0, 0.000001],\n", + ")\n", + "\n", + "for i, gp in enumerate(grid_points):\n", + " plt.text(gp[0] * 1.04, gp[1], i)\n", + "\n", + "ax.scatter(converged_samples[:, 0], converged_samples[:, 1], c=\"b\", label=\"Converged\")\n", + "ax.scatter(failed_samples[:, 0], failed_samples[:, 1], c=\"r\", label=\"Failed\")\n", + "\n", + "ax.set_xlim(10 * 0.9, 10000 * 1.2)\n", + "ax.set_yticks(np.linspace(0, 0.5, 11))\n", + "ax.set_ylim(-0.01, 0.515)\n", + "ax.set_xlabel(\"Young's modulus $E$\")\n", + "ax.set_xscale(\"log\")\n", + "ax.set_ylabel(\"Poisson ratio $\\\\nu$\")\n", + "fig.suptitle(\"Grid study 4C\")\n", + "ax.legend()\n", + "fig.tight_layout()\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Each point represents one simulation, i.e., one 4C run with different parameters. The number next to it is the job id, so the unique identifier for this simulation for QUEENS. The color of the points:\n", + "\n", + "- blue: The simulation ran through and returned a value.\n", + "- red: The simulation failed and NaN was returned.\n", + "\n", + "A colored background indicates the region of the parameter space leading to converging results." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Previously, we called 4C using the driver directly. However, when using QUEENS iterators, this won't work anymore. Instead, we need a `Model`. The difference is that a driver handles a single simulation, whereas a model accepts a list of relevant input samples / configurations.\n", + "\n", + "When working with drivers, the `Simulation` model is used. Besides a driver, this requires a `Scheduler`. As the name suggests, schedulers schedule evaluations of the driver. The scheduler also handles the evaluation of multiple driver calls at the same time. Keep in mind that each driver call can itself be parallelized with MPI, such that QUEENS natively supports nested parallelism!\n", + "\n", + "Since we are doing our computations locally, we use a `Local` scheduler. However, there is also a `Cluster` scheduler which is able to submit jobs on high performance clusters using `SLURM`, `PBS` etc.\n", + "\n", + "The scheduler also creates the required folder structure for an experiment. The default value is set to `~/queens-experiments/`. For the grid iterator, it looks like this:\n", + "```bash\n", + "grid_4x21/\n", + "├── 0\n", + "├── 1\n", + "...\n", + "├── 80\n", + "├── 81\n", + "├── 82\n", + "├── 83\n", + "├── 9\n", + "└── input_template.4C.yaml\n", + "```\n", + "\n", + "Each number is a job with its output:\n", + "```bash\n", + "10\n", + "├── input_file.yaml\n", + "├── jobscript.sh\n", + "├── metadata.yaml\n", + "└── output\n", + " ├── output.control\n", + " ├── output.log\n", + " ├── output.mesh.structure.s0\n", + " ├── output-structure.pvd\n", + " └── output-vtk-files\n", + " ├── structure-00000-0.vtu\n", + " ├── structure-00000.pvtu\n", + " ├── structure-00001-0.vtu\n", + " └── structure-00001.pvtu\n", + "```\n", + "\n", + "Once more, have a look at the files. In particular look at the log files of the failed simulations." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Let's look at some simulation outputs!\n", + "if pyvista_available:\n", + " job_ids = [0, 41, 50, 75]\n", + "\n", + " plotter = pv.Plotter(shape=(2, 2))\n", + "\n", + " for i in range(2):\n", + " for j in range(2):\n", + " plotter.subplot(i, j)\n", + " job_id = job_ids[i * 2 + j]\n", + " file_path = QUEENS_EXPERIMENTS_DIR / f\"grid_4x21/{job_id}/output/output-vtk-files/structure-00001.pvtu\"\n", + " plotter.add_text(f\"Job: {job_id}\")\n", + " plot_results(file_path, plotter, f\"E={grid_points[job_id][0]}, nu={np.round(grid_points[job_id][1],decimals=4)}: Cauchy stress zz\\n\",)\n", + " plotter.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In summary, we employed QUEENS to do some simulation analytics of a 4C model. We identified failed simulations, looked at the outputs, and extracted data." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Calibration\n", + "\n", + "Now, let's imagine we have some experimental measurements of the Cauchy stresses at the element centers $\\boldsymbol{\\sigma}_\\text{ex}$ of a real-life cantilever. We can now use this to calibrate our parameters $E$ and $\\nu$ from experiment.\n", + "\n", + "> Note: We assume that the 4C model is able to predict the real-word accurately, i.e, we have the right boundary conditions, constitutive law etc. \n", + "\n", + "To do so, we'll use deterministic optimization to minimize the difference between simulated and measured stresses:\n", + "$$(E, \\nu) \\in \\underset{E, \\nu}{\\text{argmin}} \\ \\frac{1}{n_\\text{ele}} \\sum_{e=1}^{n_\\text{ele}} (\\boldsymbol{\\sigma}_e(E, \\nu)-\\boldsymbol{\\sigma}_{\\text{ex},e}):(\\boldsymbol{\\sigma}_e(E, \\nu)-\\boldsymbol{\\sigma}_{\\text{ex},e})$$\n", + "\n", + "where $\\boldsymbol{\\sigma}_e(E,\\nu)$ are the elementwise stresses obtained from the 4C model evaluated at $E, \\nu$, where the subscript $e$ is the element id." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Before we start the optimization, we need to change our parameter space. As we previously identified, the simulation model does not converge for $\\nu > 0.45$. So in order to make meaningful predictions, let's change the parameter space." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Reduced parameter space for the poisson ratio\n", + "poisson_ratio = Uniform(0, 0.45)\n", + "\n", + "parameters = Parameters(young=young, poisson_ratio=poisson_ratio)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Let's implement the loss function\n", + "\n", + "from queens.models._model import Model\n", + "\n", + "\n", + "class LossFunction(Model):\n", + " def __init__(self, experimental_cauchy_stresses, forward_model):\n", + " self.experimental_cauchy_stresses = experimental_cauchy_stresses\n", + " self.forward_model = forward_model\n", + " super().__init__()\n", + "\n", + " def _evaluate(self, samples):\n", + " element_cauchy_stresses = self.forward_model.evaluate(samples)[\"result\"]\n", + "\n", + " difference = element_cauchy_stresses - np.tile(\n", + " self.experimental_cauchy_stresses, (len(samples), 1)\n", + " )\n", + " squared_norm = (difference**2).sum(axis=1) + (difference[:, 3:] ** 2).sum(\n", + " axis=1\n", + " ) # the output is only a symmetrical tensor, so we have to add the offdiagonal terms once more\n", + "\n", + " squared_norm /= 40 # 40 elements\n", + "\n", + " return {\"result\": squared_norm.reshape(-1)}\n", + "\n", + " def grad(self, samples, upstream_gradient):\n", + " pass\n", + "\n", + "\n", + "from experimental_data import experimental_elementwise_cauchy_stress" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Ok, now let's set up the optimization using the [L-BFGS-B](https://docs.scipy.org/doc/scipy/reference/optimize.minimize-lbfgsb.html) algorithm." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from scipy.optimize import Bounds\n", + "from queens.iterators import Optimization\n", + "\n", + "optimum = None\n", + "\n", + "if __name__ == \"__main__\":\n", + " Path(\"optimization\").mkdir(exist_ok=True)\n", + " with GlobalSettings(\n", + " experiment_name=\"optimization_young_poisson_ratio\", output_dir=\"optimization\"\n", + " ) as gs:\n", + " scheduler = Local(\n", + " gs.experiment_name,\n", + " num_jobs=4, # 4 simulations in parallel\n", + " num_procs=1, # each simulation uses 1 process\n", + " )\n", + " custom_data_processor = CustomDataProcessor(\n", + " field_name=\"element_cauchy_stresses_xyz\", # Name of the field in the pvtu file\n", + " file_name_identifier=\"output-structure.pvd\", # Name of the file in the output folder\n", + " failure_value=np.nan * np.ones(40 * 6), # NaN for all elements\n", + " files_to_be_deleted_regex_lst=[\n", + " \"output.control\",\n", + " \"*.mesh.s0\",\n", + " \"output-vtk-files/*\",\n", + " ], # Delete these files after we've extracted the data to save storage space\n", + " )\n", + " fourc_driver = Fourc(\n", + " parameters=parameters,\n", + " input_templates=input_template,\n", + " executable=fourc_executable,\n", + " data_processor=custom_data_processor, # tell the driver to extract the data\n", + " files_to_copy=[input_mesh], # copy mesh file to experiment directory\n", + " )\n", + " model = Simulation(scheduler, fourc_driver)\n", + " loss_function = LossFunction(experimental_elementwise_cauchy_stress, model)\n", + " optimization = Optimization(\n", + " loss_function,\n", + " parameters,\n", + " gs,\n", + " initial_guess=[200, 0.4],\n", + " result_description={\"write_results\": True},\n", + " objective_and_jacobian=True,\n", + " bounds=Bounds([10, 0], [10000, 0.45]),\n", + " algorithm=\"L-BFGS-B\",\n", + " max_feval=1000,\n", + " )\n", + "\n", + " run_iterator(optimization, gs)\n", + "\n", + " optimum = optimization.solution.x" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Let's do a forward simulation of the optimal parameters:\n", + "\n", + "fourc_driver.data_processor = None\n", + "\n", + "fourc_driver.run(\n", + " optimum,\n", + " experiment_dir=Path.cwd(),\n", + " job_id=\"calibrated_value\",\n", + " num_procs=1,\n", + " experiment_name=None,\n", + ")\n", + "\n", + "\n", + "def tensor_norm(tensor_data):\n", + " data = np.hstack(\n", + " (tensor_data.tolist(), tensor_data[:, 3:].tolist())\n", + " ) # duplicate offdiagonal terms\n", + " return np.sqrt((data**2).sum(axis=1)).reshape(-1, 1)\n", + "\n", + "\n", + "def plot_tensor_norm(output_mesh, field_name, data, plotter, color_bar_title):\n", + " output_mesh[field_name] = data\n", + "\n", + " plotter.add_mesh(\n", + " output_mesh,\n", + " scalars=field_name,\n", + " scalar_bar_args={\n", + " \"title\": color_bar_title,\n", + " \"title_font_size\": 15,\n", + " \"label_font_size\": 15,\n", + " },\n", + " )\n", + " plotter.add_axes(line_width=5)\n", + " plotter.camera_position = [\n", + " (1.1321899035097993, -6.851600196807601, 2.7649096132703574),\n", + " (0.0, 0.0, 0.2749999999999999),\n", + " (-0.97637930372977, -0.08995062285804697, 0.19644933366041056),\n", + " ]\n", + "\n", + "if pyvista_available:\n", + " output_mesh = pv.read(\"calibrated_value/output/output-vtk-files/structure-00001.pvtu\")\n", + " output_calibrated_value = output_mesh[\"element_cauchy_stresses_xyz\"]\n", + " difference = output_calibrated_value - experimental_elementwise_cauchy_stress.reshape(\n", + " -1, 6\n", + " )\n", + " relative_error = tensor_norm(difference) / tensor_norm(\n", + " experimental_elementwise_cauchy_stress.reshape(-1, 6)\n", + " )\n", + " plotter = pv.Plotter()\n", + " plot_tensor_norm(\n", + " output_mesh,\n", + " \"relative\",\n", + " relative_error,\n", + " plotter,\n", + " \"Relative error of the Cauchy stress norm difference between calibrated value and experimental data\\n\",\n", + " )\n", + " plotter.add_text(\n", + " f\"E={np.round(optimum[0], decimals=2)}, nu={np.round(optimum[1], decimals=2)}\",\n", + " position=\"upper_right\",\n", + " )\n", + " plotter.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We have successfully conducted a calibration of the Young's modulus $E$ and Poisson's ratio $\\nu$! The largest relative error\n", + "$$\\frac{\\sqrt{(\\boldsymbol{\\sigma}_e(E, \\nu)-\\boldsymbol{\\sigma}_{\\text{ex},e}):(\\boldsymbol{\\sigma}_e(E, \\nu)-\\boldsymbol{\\sigma}_{\\text{ex},e})}}{\\sqrt{\\boldsymbol{\\sigma}_{\\text{ex},e}:\\boldsymbol{\\sigma}_{\\text{ex},e}}}$$\n", + "between calibrated model output and experimental data is $< 1.7\\cdot10^{-8}$.\n", + "\n", + "This means we have a good state estimation, i.e., the model output matches the provided data! In fact, the data was not experimental at all, but was generated with $E=123$ and $\\nu=0.4$, which the optimizer identified quickly." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Why QUEENS?\n", + "\n", + "- Hide complexity in handling complex simulation codes:\n", + " - `Driver`: make 4C callable via a Python function call\n", + " - `Scheduler`: schedule and handle simulations in parallel\n", + " - `Model`: hide the scheduler and driver logic for the QUEENS user\n", + "- Nested parallelism\n", + " - Multiple processes per 4C run\n", + " - Multiple parallel 4C runs at the same time\n", + "- Modularity:\n", + " - Switched from grid evaluation to optimization simply by changing the iterator\n", + " - 4C simulation model did not change!" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Let's play around\n", + "\n", + "Let your creativity flow and try out stuff.\n", + "\n", + "### Inspirations\n", + "- Change the parameter ranges\n", + "- Look at the 4C input \n", + "- Look at the driver outputs\n", + "- Why did the simulations fail? (hint: two different reasons)\n", + "- Plot the displacement field" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python (queens)", + "language": "python", + "name": "queens" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.0" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/tutorials/3_orchestrating_4c_simulations/beam_coarse.exo b/tutorials/3_orchestrating_4c_simulations/beam_coarse.exo new file mode 100644 index 00000000..184ecb46 Binary files /dev/null and b/tutorials/3_orchestrating_4c_simulations/beam_coarse.exo differ diff --git a/tutorials/3_orchestrating_4c_simulations/experimental_data.py b/tutorials/3_orchestrating_4c_simulations/experimental_data.py new file mode 100644 index 00000000..7bd67c43 --- /dev/null +++ b/tutorials/3_orchestrating_4c_simulations/experimental_data.py @@ -0,0 +1,261 @@ +# +# SPDX-License-Identifier: LGPL-3.0-or-later +# Copyright (c) 2024-2025, QUEENS contributors. +# +# This file is part of QUEENS. +# +# QUEENS is free software: you can redistribute it and/or modify it under the terms of the GNU +# Lesser General Public License as published by the Free Software Foundation, either version 3 of +# the License, or (at your option) any later version. QUEENS is distributed in the hope that it will +# be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details. You +# should have received a copy of the GNU Lesser General Public License along with QUEENS. If not, +# see . +# +"""Data for tutorial 4.""" +import numpy as np + +experimental_elementwise_cauchy_stress = np.array( + [ + 4.965463975769041, + 4.965463975769033, + 17.51894212636035, + 0.03002653745373526, + -1.3302958186504616, + -1.3302958186504539, + 0.028146592901132292, + 0.028146592901133374, + 18.636951800218696, + 0.012885974118424563, + 0.17919017027298328, + 0.17919017027297202, + -0.042180090097681235, + -0.04218009009768057, + 18.650502468387778, + 0.03464416037371552, + 0.052885158232057194, + 0.05288515823204195, + -0.009775365334722787, + -0.00977536533473154, + 18.643362157159533, + 0.0016409123275933793, + 0.006633901693475661, + 0.006633901693473053, + -0.0010330150043624002, + -0.0010330150043614053, + 18.641487576354407, + 0.00021339106204969465, + 0.0005422584723292606, + 0.0005422584723294157, + -0.0010330150043516542, + -0.0010330150043471765, + 18.641487576354415, + 0.00021339106205215902, + -0.0005422584723329064, + -0.0005422584723276938, + -0.009775365334704078, + -0.00977536533470786, + 18.643362157159565, + 0.0016409123275951348, + -0.0066339016934708, + -0.006633901693459221, + -0.04218009009768332, + -0.04218009009768949, + 18.65050246838775, + 0.034644160373714546, + -0.05288515823203891, + -0.05288515823201958, + 0.028146592901169887, + 0.028146592901162837, + 18.6369518002187, + 0.012885974118419446, + -0.17919017027296463, + -0.17919017027297335, + 4.965463975769042, + 4.965463975769036, + 17.518942126360358, + 0.03002653745372676, + 1.330295818650467, + 1.3302958186504656, + 4.965463975768985, + 4.965463975768996, + 17.518942126360265, + -0.030026537453727027, + 1.3302958186504579, + -1.3302958186504907, + 0.028146592901138225, + 0.02814659290114923, + 18.636951800218714, + -0.012885974118421925, + -0.1791901702729771, + 0.17919017027296902, + -0.04218009009770295, + -0.04218009009770612, + 18.650502468387746, + -0.03464416037371499, + -0.05288515823204254, + 0.05288515823201454, + -0.009775365334677315, + -0.009775365334677804, + 18.64336215715961, + -0.0016409123275940257, + -0.006633901693475254, + 0.00663390169346114, + -0.0010330150043375255, + -0.0010330150043386197, + 18.641487576354447, + -0.00021339106205258012, + -0.0005422584723254254, + 0.0005422584723254766, + -0.0010330150043826982, + -0.0010330150043885683, + 18.6414875763544, + -0.00021339106205113003, + 0.0005422584723280942, + -0.0005422584723277691, + -0.009775365334733528, + -0.009775365334734625, + 18.643362157159515, + -0.001640912327592868, + 0.006633901693463425, + -0.006633901693471567, + -0.04218009009769092, + -0.04218009009769621, + 18.65050246838773, + -0.03464416037371196, + 0.05288515823201943, + -0.0528851582320438, + 0.028146592901187772, + 0.028146592901194423, + 18.636951800218775, + -0.01288597411842288, + 0.17919017027295925, + -0.1791901702729603, + 4.965463975768974, + 4.965463975768976, + 17.518942126360244, + -0.030026537453731753, + -1.3302958186504656, + 1.3302958186504699, + 4.965463975769046, + 4.965463975769051, + 17.51894212636036, + -0.030026537453726874, + -1.3302958186504585, + 1.3302958186504716, + 0.028146592901088938, + 0.028146592901090135, + 18.636951800218597, + -0.012885974118420528, + 0.1791901702729778, + -0.17919017027297107, + -0.04218009009768932, + -0.04218009009768941, + 18.650502468387778, + -0.03464416037371437, + 0.05288515823203733, + -0.05288515823202167, + -0.009775365334708455, + -0.009775365334703481, + 18.643362157159565, + -0.001640912327593154, + 0.006633901693467521, + -0.00663390169346678, + -0.001033015004370161, + -0.0010330150043544404, + 18.64148757635443, + -0.00021339106205401455, + 0.00054225847232167, + -0.0005422584723292693, + -0.001033015004361405, + -0.0010330150043472762, + 18.64148757635443, + -0.0002133910620527629, + -0.0005422584723280483, + 0.0005422584723261902, + -0.009775365334742384, + -0.009775365334735817, + 18.6433621571595, + -0.0016409123275928927, + -0.006633901693463455, + 0.006633901693462081, + -0.04218009009767108, + -0.042180090097671305, + 18.650502468387813, + -0.03464416037371468, + -0.05288515823202985, + 0.05288515823203297, + 0.028146592901229596, + 0.028146592901220215, + 18.63695180021883, + -0.012885974118419248, + -0.17919017027295764, + 0.1791901702729505, + 4.965463975769022, + 4.965463975769019, + 17.518942126360287, + -0.030026537453728092, + 1.3302958186504648, + -1.3302958186504628, + 4.965463975769042, + 4.9654639757690475, + 17.51894212636036, + 0.03002653745372974, + 1.330295818650454, + 1.3302958186504577, + 0.02814659290113046, + 0.02814659290113048, + 18.636951800218633, + 0.012885974118422452, + -0.1791901702729813, + -0.17919017027297185, + -0.04218009009767682, + -0.04218009009768283, + 18.65050246838779, + 0.034644160373716835, + -0.05288515823204267, + -0.05288515823202154, + -0.009775365334698504, + -0.009775365334699599, + 18.64336215715957, + 0.001640912327594638, + -0.006633901693466627, + -0.006633901693474796, + -0.0010330150043444905, + -0.0010330150043428988, + 18.64148757635444, + 0.00021339106205077853, + -0.0005422584723261475, + -0.0005422584723291357, + -0.0010330150043227, + -0.001033015004327078, + 18.641487576354503, + 0.00021339106205158948, + 0.0005422584723278418, + 0.0005422584723235945, + -0.009775365334732734, + -0.009775365334741087, + 18.6433621571595, + 0.0016409123275929834, + 0.006633901693458997, + 0.00663390169346926, + -0.042180090097676856, + -0.04218009009767581, + 18.650502468387796, + 0.03464416037371584, + 0.05288515823202486, + 0.052885158232036794, + 0.028146592901174397, + 0.028146592901190655, + 18.63695180021874, + 0.012885974118425469, + 0.17919017027296008, + 0.17919017027297757, + 4.965463975769032, + 4.965463975769038, + 17.518942126360354, + 0.030026537453729938, + -1.3302958186504696, + -1.3302958186504688, + ] +) diff --git a/tutorials/3_orchestrating_4c_simulations/input_template.4C.yaml b/tutorials/3_orchestrating_4c_simulations/input_template.4C.yaml new file mode 100644 index 00000000..db6fd21c --- /dev/null +++ b/tutorials/3_orchestrating_4c_simulations/input_template.4C.yaml @@ -0,0 +1,56 @@ +TITLE: + QUEENS injected: + poisson_ratio: {{ poisson_ratio }} + young: {{ young }} +PROBLEM TYPE: + PROBLEMTYPE: "Structure" +STRUCTURAL DYNAMIC: + INT_STRATEGY: "Standard" + DYNAMICTYPE: "Statics" + NUMSTEP: 1 + MAXTIME: 1 + MAXITER: 7 + RESULTSEVERY: 0 + RESTARTEVERY: 0 + DIVERCONT: stop + LINEAR_SOLVER: 1 + TIMESTEP: 1 +IO: + VERBOSITY: "standard" + STRUCT_STRESS: Cauchy +IO/RUNTIME VTK OUTPUT: + INTERVAL_STEPS: 1 +IO/RUNTIME VTK OUTPUT/STRUCTURE: + OUTPUT_STRUCTURE: true + DISPLACEMENT: true + STRESS_STRAIN: true +SOLVER 1: + SOLVER: "UMFPACK" + NAME: "Direct Solver" +MATERIALS: + - MAT: 1 + MAT_Struct_StVenantKirchhoff: + YOUNG: {{ young }} + NUE: {{ poisson_ratio }} + DENS: 1 +DESIGN SURF DIRICH CONDITIONS: + - E: 2 + ENTITY_TYPE: node_set_id + NUMDOF: 3 + ONOFF: [1, 1, 1] + VAL: [0, 0, 0] + FUNCT: [null, null, null] + - E: 1 + ENTITY_TYPE: node_set_id + NUMDOF: 3 + ONOFF: [1, 1, 1] + VAL: [0, 0, 0.55] + FUNCT: [null, null, null] +STRUCTURE GEOMETRY: + FILE: ../beam_coarse.exo + ELEMENT_BLOCKS: + - ID: 1 + SOLID: + HEX8: + MAT: 1 + KINEM: nonlinear diff --git a/tutorials/3_orchestrating_4c_simulations/utils.py b/tutorials/3_orchestrating_4c_simulations/utils.py new file mode 100644 index 00000000..6ddeb74a --- /dev/null +++ b/tutorials/3_orchestrating_4c_simulations/utils.py @@ -0,0 +1,72 @@ +# +# SPDX-License-Identifier: LGPL-3.0-or-later +# Copyright (c) 2024-2025, QUEENS contributors. +# +# This file is part of QUEENS. +# +# QUEENS is free software: you can redistribute it and/or modify it under the terms of the GNU +# Lesser General Public License as published by the Free Software Foundation, either version 3 of +# the License, or (at your option) any later version. QUEENS is distributed in the hope that it will +# be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details. You +# should have received a copy of the GNU Lesser General Public License along with QUEENS. If not, +# see . +# +"""Utility functions for tutorial 4.""" +import os +from pathlib import Path + +os.environ["PYVISTA_OFF_SCREEN"] = "true" +os.environ["VTK_DEFAULT_RENDER_WINDOW_OFFSCREEN"] = "1" +os.environ["PYVISTA_JUPYTER_BACKEND"] = "static" + +import pyvista as pv + +fe_mesh = pv.read("beam_coarse.exo") + + +def plot_results( + result_file: str, + plotter: pv.Plotter, + color_bar_title: str = "zz-component of Cauchy stress tensor\n", +) -> None: + """Plots results. + + Args: + result_file (str): Path to pvtu file. + plotter (pyvista Plotter): pyvista plotter. + color_bar_title (string): title of the colorbar. + """ + outputs = pv.read(result_file).warp_by_vector("displacement") + plotter.add_mesh(fe_mesh.copy(), style="wireframe", color="blue") + outputs["cauchy_zz"] = outputs["element_cauchy_stresses_xyz"][:, 2] + plotter.add_mesh( + outputs, + scalars="cauchy_zz", + scalar_bar_args={ + "title": color_bar_title, + "title_font_size": 15, + "label_font_size": 15, + }, + ) + plotter.add_axes(line_width=5) # type: ignore[call-arg] + plotter.camera_position = [ + (1.1321899035097993, -6.851600196807601, 2.7649096132703574), + (0.0, 0.0, 0.2749999999999999), + (-0.97637930372977, -0.08995062285804697, 0.19644933366041056), + ] + + +def find_repo_root(start: Path) -> Path: + """Find the repository root by going upwards from a starting directory. + + Args: + start (Path): Starting directory to begin the upward search. + + Returns: + Path: Path to the detected repository root (containing ``pyproject.toml``). + """ + for p in [start, *start.parents]: + if (p / "pyproject.toml").exists(): + return p + raise RuntimeError(f"Could not find repo root from {start}")