{ "cells": [ { "cell_type": "markdown", "source": [ "# 14: Adaptive mesh refinement" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "**Note:** To improve responsiveness via caching, the notebooks are updated only once a week. They are only\n", "available for the latest stable release of Trixi.jl at the time of caching." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Adaptive mesh refinement (AMR) is a method of adapting the resolution of the numerical method\n", "to the solution features such as turbulent regions or shocks. In those critical regions\n", "of the domain, we want the simulation to use elements with smaller mesh sizes compared to other\n", "regions. This should be automatically and dynamically adapted during the run of the simulation." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "# Implementation in Trixi.jl\n", "In [Trixi.jl](https://github.com/trixi-framework/Trixi.jl), AMR is possible for the mesh types\n", "`TreeMesh` and `P4estMesh`. Both meshes are organized in a tree structure\n", "and therefore, each element can be refined independently. In Trixi.jl, AMR is restricted\n", "to a 2:1 refinement ratio between neighbor elements. This means that the maximum resolution\n", "difference of neighboring elements is a factor of two." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The implementation of AMR is divided into different steps. The basic refinement setting contains\n", "an indicator and a controller. These are added to the simulation by using an AMR callback." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "### Indicators\n", "An indicator estimates the current accuracy of the numerical approximation. It indicates which regions\n", "of the domain need finer or coarser resolutions. In Trixi.jl, you can use for instance\n", "`IndicatorLöhner` and `IndicatorHennemannGassner`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "`IndicatorLöhner` (also callable with `IndicatorLoehner`) is an interpretation and adaptation of\n", "a FEM indicator by [Löhner (1987)](https://doi.org/10.1016/0045-7825(87)90098-3) and estimates a\n", "weighted second derivative of a specified variable locally.\n", "````julia\n", "amr_indicator = IndicatorLöhner(semi, variable=variable)\n", "````\n", "All indicators have the parameter `variable` which is used to specify the variable for the\n", "indicator calculation. You can use for instance `density`, `pressure` or `density_pressure`\n", "for the compressible Euler equations. Moreover, you have the option to use simply the first\n", "conservation variable with `first` for any equations. This might be a good choice for a starting\n", "example." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "`IndicatorHennemannGassner`, also used as a shock-capturing indicator, was developed by\n", "[Hennemann et al. (2021)](https://doi.org/10.1016/j.jcp.2020.109935) and is explained in detail\n", "in the tutorial about shock-capturing. It can be constructed as follows.\n", "````julia\n", "amr_indicator = IndicatorHennemannGassner(semi,\n", " alpha_max=0.5,\n", " alpha_min=0.001,\n", " alpha_smooth=true,\n", " variable=variable)\n", "````" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Another indicator is the very basic `IndicatorMax`. It indicates the maximal value of a variable\n", "and is therefore mostly used for verification and testing. But it might be useful for the basic\n", "understanding of the implementation of indicators and AMR in Trixi.jl.\n", "````julia\n", "amr_indicator = IndicatorMax(semi, variable=variable)\n", "````" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "### Controllers\n", "The spatial discretization into elements is tree-based for both AMR supporting mesh types `TreeMesh`\n", "and `P4estMesh`. Thus, the higher the level in the tree the higher the level of refinement.\n", "For instance, a mesh element of level `3` has double resolution in each direction compared to\n", "another element with level `2`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "To map specific indicator values to a desired level of refinement, Trixi.jl uses controllers.\n", "They are build in three levels: There is a base level of refinement `base_level`, which is the\n", "minimum allowed refinement level. Then, there is a medium level `med_level`, which corresponds\n", "to the initial level of refinement, for indicator values above the threshold `med_threshold`\n", "and equally, a maximal level `max_level` for values above `max_threshold`.\n", "This variant of controller is called `ControllerThreeLevel` in Trixi.jl.\n", "````julia\n", "amr_controller = ControllerThreeLevel(semi, amr_indicator;\n", " base_level=4,\n", " med_level=5, med_threshold=0.1,\n", " max_level=6, max_threshold=0.6)\n", "````\n", "You can also set `med_level=0` to use the current level as target, see the docstring of\n", "`ControllerThreeLevel`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "An extension is `ControllerThreeLevelCombined`, which uses two different indicators.\n", "The primary indicator works the same as the single indicator for `ControllerThreeLevel`.\n", "The second indicator with its own maximum threshold adds the property, that the target level is set to\n", "`max_level` additionally if this indicator's value is greater than `max_threshold_secondary`.\n", "This is for instance used to assure that a shock has always the maximum refinement level.\n", "````julia\n", "amr_controller = ControllerThreeLevelCombined(semi, indicator_primary, indicator_secondary;\n", " base_level=2,\n", " med_level=6, med_threshold=0.0003,\n", " max_level=8, max_threshold=0.003,\n", " max_threshold_secondary=0.3)\n", "````\n", "This controller is for instance used in\n", "[`elixir_euler_astro_jet_amr.jl`](https://github.com/trixi-framework/Trixi.jl/blob/main/examples/tree_2d_dgsem/elixir_euler_astro_jet_amr.jl)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "### Callback\n", "The AMR indicator and controller are added to the simulation through the callback `AMRCallback`.\n", "It contains a semidiscretization `semi`, the controller `amr_controller` and the parameters `interval`,\n", "`adapt_initial_condition`, and `adapt_initial_condition_only_refine`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Adaptive mesh refinement will be performed every `interval` time steps. `adapt_initial_condition` indicates\n", "whether the initial condition already should be adapted before the first time step. And with\n", "`adapt_initial_condition_only_refine=true` the mesh is only refined at the beginning but not coarsened.\n", "````julia\n", "amr_callback = AMRCallback(semi, amr_controller,\n", " interval=5,\n", " adapt_initial_condition=true,\n", " adapt_initial_condition_only_refine=true)\n", "````" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "# Exemplary simulation\n", "Here, we want to implement a simple AMR simulation of the 2D linear advection equation for a Gaussian pulse." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using OrdinaryDiffEq\n", "using Trixi\n", "\n", "advection_velocity = (0.2, -0.7)\n", "equations = LinearScalarAdvectionEquation2D(advection_velocity)\n", "\n", "initial_condition = initial_condition_gauss\n", "solver = DGSEM(polydeg=3, surface_flux=flux_lax_friedrichs)\n", "\n", "coordinates_min = (-5.0, -5.0)\n", "coordinates_max = ( 5.0, 5.0)\n", "mesh = TreeMesh(coordinates_min, coordinates_max,\n", " initial_refinement_level=4,\n", " n_cells_max=30_000)\n", "\n", "semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver)\n", "\n", "\n", "tspan = (0.0, 10.0)\n", "ode = semidiscretize(semi, tspan);" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "For the best understanding about indicators and controllers, we use the simple AMR indicator\n", "`IndicatorMax`. As described before, it returns the maximal value of the specified variable\n", "(here the only conserved variable). Therefore, regions with a high maximum are refined.\n", "This is not really useful numerical application, but a nice demonstration example." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "amr_indicator = IndicatorMax(semi, variable=first)" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "These values are transferred to a refinement level with the `ControllerThreeLevel`, such that\n", "every element with maximal value greater than `0.1` is refined once and elements with maximum\n", "above `0.6` are refined twice." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "amr_controller = ControllerThreeLevel(semi, amr_indicator,\n", " base_level=4,\n", " med_level=5, med_threshold=0.1,\n", " max_level=6, max_threshold=0.6)\n", "\n", "amr_callback = AMRCallback(semi, amr_controller,\n", " interval=5,\n", " adapt_initial_condition=true,\n", " adapt_initial_condition_only_refine=true)\n", "\n", "stepsize_callback = StepsizeCallback(cfl=0.9)\n", "\n", "callbacks = CallbackSet(amr_callback, stepsize_callback);" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "Running the simulation." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "sol = solve(ode, CarpenterKennedy2N54(williamson_condition=false),\n", " dt=1.0, # solve needs some value here but it will be overwritten by the stepsize_callback\n", " save_everystep=false, callback=callbacks);" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "We plot the solution and add the refined mesh at the end of the simulation." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using Plots\n", "pd = PlotData2D(sol)\n", "plot(pd)\n", "plot!(getmesh(pd))" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "# More examples\n", "Trixi.jl provides many elixirs using AMR. We want to give some examples for different mesh types:\n", "- `elixir_euler_blast_wave_amr.jl` for [`TreeMesh`](https://github.com/trixi-framework/Trixi.jl/blob/main/examples/tree_2d_dgsem/elixir_euler_blast_wave_amr.jl)\n", " and [`P4estMesh`](https://github.com/trixi-framework/Trixi.jl/blob/main/examples/p4est_2d_dgsem/elixir_euler_blast_wave_amr.jl)\n", "- [`elixir_euler_kelvin_helmholtz_instability_amr.jl`](https://github.com/trixi-framework/Trixi.jl/blob/main/examples/tree_2d_dgsem/elixir_euler_kelvin_helmholtz_instability_amr.jl) for `TreeMesh`\n", "- [`elixir_euler_double_mach_amr.jl`](https://github.com/trixi-framework/Trixi.jl/blob/main/examples/p4est_2d_dgsem/elixir_euler_double_mach_amr.jl) for `P4estMesh`" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Animations of more interesting and complicated AMR simulations can be found below and on Trixi.jl's youtube channel\n", "[\"Trixi Framework\"](https://www.youtube.com/channel/UCpd92vU2HjjTPup-AIN0pkg)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "First, we give a [purely hyperbolic simulation of a Sedov blast wave with self-gravity](https://www.youtube.com/watch?v=dxgzgteJdOA).\n", "This simulation uses the mesh type `TreeMesh` as we did and the AMR indicator `IndicatorHennemannGassner`.\n", "\n", " \n", "
\n", "\n", "Source: Trixi.jl's YouTube channel [`Trixi Framework`](https://www.youtube.com/channel/UCpd92vU2HjjTPup-AIN0pkg)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The next example is a numerical simulation of an [ideal MHD rotor on an unstructured AMR mesh](https://www.youtube.com/watch?v=Iei7e9oQ0hs).\n", "The used mesh type is a `P4estMesh`.\n", "\n", " \n", "
\n", "\n", "Source: Trixi.jl's YouTube channel [`Trixi Framework`](https://www.youtube.com/channel/UCpd92vU2HjjTPup-AIN0pkg)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "For more information, please have a look at the respective links." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Package versions" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "These results were obtained using the following versions." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using InteractiveUtils\n", "versioninfo()\n", "\n", "using Pkg\n", "Pkg.status([\"Trixi\", \"OrdinaryDiffEq\", \"Plots\"],\n", " mode=PKGMODE_MANIFEST)" ], "metadata": {}, "execution_count": null } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.10.3" }, "kernelspec": { "name": "julia-1.10", "display_name": "Julia 1.10.3", "language": "julia" } }, "nbformat": 4 }