{ "cells": [ { "cell_type": "markdown", "source": [ "# 16: Unstructured meshes with HOHQMesh.jl" ], "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": [ "Trixi.jl supports numerical approximations on unstructured quadrilateral meshes\n", "with the `UnstructuredMesh2D` mesh type." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The purpose of this tutorial is to demonstrate how to use the `UnstructuredMesh2D`\n", "functionality of Trixi.jl. This begins by running and visualizing an available unstructured\n", "quadrilateral mesh example. Then, the tutorial will demonstrate how to\n", "conceptualize a problem with curved boundaries, generate\n", "a curvilinear mesh using the available software in the Trixi.jl ecosystem,\n", "and then run a simulation using Trixi.jl on said mesh." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Unstructured quadrilateral meshes can be made\n", "with the [High-Order Hex-Quad Mesh (HOHQMesh) generator](https://github.com/trixi-framework/HOHQMesh)\n", "created and developed by David Kopriva.\n", "HOHQMesh is a mesh generator specifically designed for spectral element methods.\n", "It provides high-order boundary curve information (needed to accurately set boundary conditions)\n", "and elements can be larger (due to the high accuracy of the spatial approximation)\n", "compared to traditional finite element mesh generators.\n", "For more information about the design and features of HOHQMesh one can refer to its\n", "[official documentation](https://trixi-framework.github.io/HOHQMesh/)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "HOHQMesh is incorporated into the Trixi.jl framework via the registered Julia package\n", "[HOHQMesh.jl](https://github.com/trixi-framework/HOHQMesh.jl).\n", "This package provides a Julia wrapper for the HOHQMesh generator that allows users to easily create mesh\n", "files without the need to build HOHQMesh from source. To install the HOHQMesh package execute\n", "```julia\n", "import Pkg; Pkg.add(\"HOHQMesh\")\n", "```\n", "Now we are ready to generate an unstructured quadrilateral mesh that can be used by Trixi.jl." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Running and visualizing an unstructured simulation" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Trixi.jl supports solving hyperbolic problems on several mesh types.\n", "There is a default example for this mesh type that can be executed by" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using Trixi\n", "trixi_include(default_example_unstructured())" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "This will compute a smooth, manufactured solution test case for the 2D compressible Euler equations\n", "on the curved quadrilateral mesh described in the\n", "[Trixi.jl documentation](https://trixi-framework.github.io/Trixi.jl/stable/meshes/unstructured_quad_mesh/)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Apart from the usual error and timing output provided by the Trixi.jl run, it is useful to visualize and inspect\n", "the solution. One option available in the Trixi.jl framework to visualize the solution on\n", "unstructured quadrilateral meshes is post-processing the\n", "Trixi.jl output file(s) with the [`Trixi2Vtk`](https://github.com/trixi-framework/Trixi2Vtk.jl) tool\n", "and plotting them with [ParaView](https://www.paraview.org/download/)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "To convert the HDF5-formatted `.h5` output file(s) from Trixi.jl into VTK format execute the following" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using Trixi2Vtk\n", "trixi2vtk(\"out/solution_000180.h5\", output_directory=\"out\")" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "Note this step takes about 15-30 seconds as the package `Trixi2Vtk` must be precompiled and executed for the first time\n", "in your REPL session. The `trixi2vtk` command above will convert the solution file at the final time into a `.vtu` file\n", "which can be read in and visualized with ParaView. Optional arguments for `trixi2vtk` are: (1) Pointing to the `output_directory`\n", "where the new files will be saved; it defaults to the current directory. (2) Specifying a higher number of\n", "visualization nodes. For instance, if we want to use 12 uniformly spaced nodes for visualization we can execute" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "trixi2vtk(\"out/solution_000180.h5\", output_directory=\"out\", nvisnodes=12)" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "By default `trixi2vtk` sets `nvisnodes` to be the same as the number of nodes specified in\n", "the `elixir` file used to run the simulation." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Finally, if you want to convert all the solution files to VTK execute" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "trixi2vtk(\"out/solution_000*.h5\", output_directory=\"out\", nvisnodes=12)" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "then it is possible to open the `.pvd` file with ParaView and create a video of the simulation." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Creating a mesh using HOHQMesh" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The creation of an unstructured quadrilateral mesh using HOHQMesh.jl is driven by a **control file**. In this file the user dictates\n", "the domain to be meshed, prescribes any desired boundary curvature, the polynomial order of said boundaries, etc.\n", "In this tutorial we cover several basic features of the possible control inputs. For a complete discussion\n", "on this topic see the [HOHQMesh control file documentation](https://trixi-framework.github.io/HOHQMesh/the-control-file/)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "To begin, we provide a complete control file in this tutorial. After this we give a breakdown\n", "of the control file components to explain the chosen parameters." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Suppose we want to create a mesh of a domain with straight sided\n", "outer boundaries and a curvilinear \"ice cream cone\" shaped object at its center." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "![mesh_boundary_cartoon](https://user-images.githubusercontent.com/25242486/129603954-9788500d-bba8-49be-8e6f-7555099dbf7c.png)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The associated `ice_cream_straight_sides.control` file is created below." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "open(\"out/ice_cream_straight_sides.control\", \"w\") do io\n", " println(io, raw\"\"\"\n", "\\begin{CONTROL_INPUT}\n", " \\begin{RUN_PARAMETERS}\n", " mesh file name = ice_cream_straight_sides.mesh\n", " plot file name = ice_cream_straight_sides.tec\n", " stats file name = none\n", " mesh file format = ISM-v2\n", " polynomial order = 4\n", " plot file format = skeleton\n", " \\end{RUN_PARAMETERS}\n", "\n", " \\begin{BACKGROUND_GRID}\n", " x0 = [-8.0, -8.0, 0.0]\n", " dx = [1.0, 1.0, 0.0]\n", " N = [16,16,1]\n", " \\end{BACKGROUND_GRID}\n", "\n", " \\begin{SPRING_SMOOTHER}\n", " smoothing = ON\n", " smoothing type = LinearAndCrossBarSpring\n", " number of iterations = 25\n", " \\end{SPRING_SMOOTHER}\n", "\n", "\\end{CONTROL_INPUT}\n", "\n", "\\begin{MODEL}\n", "\n", " \\begin{INNER_BOUNDARIES}\n", "\n", " \\begin{CHAIN}\n", " name = IceCreamCone\n", " \\begin{END_POINTS_LINE}\n", " name = LeftSlant\n", " xStart = [-2.0, 1.0, 0.0]\n", " xEnd = [ 0.0, -3.0, 0.0]\n", " \\end{END_POINTS_LINE}\n", "\n", " \\begin{END_POINTS_LINE}\n", " name = RightSlant\n", " xStart = [ 0.0, -3.0, 0.0]\n", " xEnd = [ 2.0, 1.0, 0.0]\n", " \\end{END_POINTS_LINE}\n", "\n", " \\begin{CIRCULAR_ARC}\n", " name = IceCream\n", " units = degrees\n", " center = [ 0.0, 1.0, 0.0]\n", " radius = 2.0\n", " start angle = 0.0\n", " end angle = 180.0\n", " \\end{CIRCULAR_ARC}\n", " \\end{CHAIN}\n", "\n", " \\end{INNER_BOUNDARIES}\n", "\n", "\\end{MODEL}\n", "\\end{FILE}\n", "\"\"\")\n", "end" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "The first three blocks of information are wrapped within a `CONTROL_INPUT` environment block as they define the\n", "core components of the quadrilateral mesh that will be generated." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The first block of information in `RUN_PARAMETERS` is\n", "```\n", "\\begin{RUN_PARAMETERS}\n", " mesh file name = ice_cream_straight_sides.mesh\n", " plot file name = ice_cream_straight_sides.tec\n", " stats file name = none\n", " mesh file format = ISM-v2\n", " polynomial order = 4\n", " plot file format = skeleton\n", "\\end{RUN_PARAMETERS}\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The mesh and plot file names will be the files created by HOHQMesh once successfully executed. The stats file name is\n", "available if you wish to also save a collection of mesh statistics. For this example it is deactivated.\n", "These file names given within `RUN_PARAMETERS` **should match** that of the control file, and although this is not required by\n", "HOHQMesh, it is a useful style convention.\n", "The mesh file format `ISM-v2` in the format currently required by Trixi.jl. The `polynomial order` prescribes the order\n", "of an interpolant constructed on the Chebyshev-Gauss-Lobatto nodes that is used to represent any curved boundaries on a particular element.\n", "The plot file format of `skeleton` means that visualizing the plot file will only draw the element boundaries (and no internal nodes).\n", "Alternatively, the format can be set to `sem` to visualize the interior nodes of the approximation as well." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The second block of information in `BACKGROUND_GRID` is\n", "```\n", "\\begin{BACKGROUND_GRID}\n", " x0 = [-8.0, -8.0, 0.0]\n", " dx = [1.0, 1.0, 0.0]\n", " N = [16,16,1]\n", "\\end{BACKGROUND_GRID}\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "This lays a grid of Cartesian elements for the domain beginning at the point `x0` as its bottom-left corner.\n", "The value of `dx`, which could differ in each direction if desired, controls the step size taken in each Cartesian direction.\n", "The values in `N` set how many Cartesian box elements are set in each coordinate direction.\n", "The above parameters define a $16\\times 16$ element square mesh on $[-8,8]^2$.\n", "Further, this sets up four outer boundaries of the domain that are given the default names: `Top, Left, Bottom, Right`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The third block of information in `SPRING_SMOOTHER` is\n", "```\n", "\\begin{SPRING_SMOOTHER}\n", " smoothing = ON\n", " smoothing type = LinearAndCrossBarSpring\n", " number of iterations = 25\n", "\\end{SPRING_SMOOTHER}\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Once HOHQMesh generates the mesh, a spring-mass-dashpot model is created to smooth the mesh and create \"nicer\" quadrilateral elements.\n", "The [default parameters of Hooke's law](https://trixi-framework.github.io/HOHQMesh/the-control-input/#the-smoother)\n", "for the spring-mass-dashpot model have been selected after a fair amount of experimentation across many meshes.\n", "If you wish to deactivate this feature you can set `smoothing = OFF` (or remove this block from the control file)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "After the `CONTROL_INPUT` environment block comes the `MODEL` environment block. It is here where the user prescribes curved boundary information with either:\n", "* An `OUTER_BOUNDARY` (covered in the next section of this tutorial).\n", "* One or more `INNER_BOUNDARIES`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "There are several options to describe the boundary curve data to HOHQMesh like splines or parametric curves." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "For the example `ice_cream_straight_sides.control` we define three internal boundaries; two straight-sided and\n", "one as a circular arc.\n", "Within the HOHQMesh control input each curve must be assigned to a `CHAIN` as shown below in the complete\n", "`INNER_BOUNDARIES` block.\n", "```\n", "\\begin{INNER_BOUNDARIES}\n", "\n", " \\begin{CHAIN}\n", " name = IceCreamCone\n", " \\begin{END_POINTS_LINE}\n", " name = LeftSlant\n", " xStart = [-2.0, 1.0, 0.0]\n", " xEnd = [ 0.0, -3.0, 0.0]\n", " \\end{END_POINTS_LINE}\n", "\n", " \\begin{END_POINTS_LINE}\n", " name = RightSlant\n", " xStart = [ 0.0, -3.0, 0.0]\n", " xEnd = [ 2.0, 1.0, 0.0]\n", " \\end{END_POINTS_LINE}\n", "\n", " \\begin{CIRCULAR_ARC}\n", " name = IceCream\n", " units = degrees\n", " center = [ 0.0, 1.0, 0.0]\n", " radius = 2.0\n", " start angle = 0.0\n", " end angle = 180.0\n", " \\end{CIRCULAR_ARC}\n", " \\end{CHAIN}\n", "\n", "\\end{INNER_BOUNDARIES}\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "It is important to note there are two `name` quantities one for the `CHAIN` and one for the `PARAMETRIC_EQUATION_CURVE`.\n", "The name for the `CHAIN` is used internally by HOHQMesh, so if you have multiple `CHAIN`s they **must be given a unique name**.\n", "The name for the `PARAMETRIC_EQUATION_CURVE` will be printed to the appropriate boundaries within the `.mesh` file produced by\n", "HOHQMesh." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We create the mesh file `ice_cream_straight_sides.mesh` and its associated file for plotting\n", "`ice_cream_straight_sides.tec` by using HOHQMesh.jl's function `generate_mesh`." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using HOHQMesh\n", "control_file = joinpath(\"out\", \"ice_cream_straight_sides.control\")\n", "output = generate_mesh(control_file);" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "The mesh file `ice_cream_straight_sides.mesh` and its associated file for plotting\n", "`ice_cream_straight_sides.tec` are placed in the `out` folder.\n", "The resulting mesh generated by HOHQMesh.jl is given in the following figure." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "![mesh_straight_sides](https://user-images.githubusercontent.com/25242486/129603958-08e4b874-53d5-4511-9a54-6daf4c21edca.png)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We note that Trixi.jl uses the boundary name information from the control file\n", "to assign boundary conditions in an elixir file.\n", "Therefore, the name should start with a letter and consist only of alphanumeric characters and underscores. Please note that the name will be treated as case sensitive." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Example simulation on `ice_cream_straight_sides.mesh`" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "With this newly generated mesh we are ready to run a Trixi.jl simulation on an unstructured quadrilateral mesh.\n", "For this we must create a new elixir file." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The elixir file given below creates an initial condition for a\n", "uniform background flow state with a free stream Mach number of 0.3.\n", "A focus for this part of the tutorial is to specify the boundary conditions and to construct the new mesh from the\n", "file that was generated in the previous exercise." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "It is straightforward to set the different boundary\n", "condition types in an elixir by assigning a particular function to a boundary name inside a\n", "Julia dictionary, `Dict`, variable. Observe that the names of these boundaries match those provided by HOHQMesh\n", "either by default, e.g. `Bottom`, or user assigned, e.g. `IceCream`. For this problem setup use\n", "* Freestream boundary conditions on the four box edges.\n", "* Free slip wall boundary condition on the interior curved boundaries." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "To construct the unstructured quadrilateral mesh from the HOHQMesh file we point to the appropriate location\n", "with the variable `mesh_file` and then feed this into the constructor for the `UnstructuredMesh2D` type in Trixi.jl" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "```julia\n", "# create the unstructured mesh from your mesh file\n", "using Trixi\n", "mesh_file = joinpath(\"out\", \"ice_cream_straight_sides.mesh\")\n", "mesh = UnstructuredMesh2D(mesh_file);\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The complete elixir file for this simulation example is given below." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using OrdinaryDiffEq, Trixi\n", "\n", "equations = CompressibleEulerEquations2D(1.4) # set gas gamma = 1.4\n", "\n", "# freestream flow state with Ma_inf = 0.3\n", "@inline function uniform_flow_state(x, t, equations::CompressibleEulerEquations2D)\n", "\n", " # set the freestream flow parameters\n", " rho_freestream = 1.0\n", " u_freestream = 0.3\n", " p_freestream = inv(equations.gamma)\n", "\n", " theta = 0.0 # zero angle of attack\n", " si, co = sincos(theta)\n", " v1 = u_freestream * co\n", " v2 = u_freestream * si\n", "\n", " prim = SVector(rho_freestream, v1, v2, p_freestream)\n", " return prim2cons(prim, equations)\n", "end\n", "\n", "# initial condition\n", "initial_condition = uniform_flow_state\n", "\n", "# boundary condition types\n", "boundary_condition_uniform_flow = BoundaryConditionDirichlet(uniform_flow_state)\n", "\n", "# boundary condition dictionary\n", "boundary_conditions = Dict( :Bottom => boundary_condition_uniform_flow,\n", " :Top => boundary_condition_uniform_flow,\n", " :Right => boundary_condition_uniform_flow,\n", " :Left => boundary_condition_uniform_flow,\n", " :LeftSlant => boundary_condition_slip_wall,\n", " :RightSlant => boundary_condition_slip_wall,\n", " :IceCream => boundary_condition_slip_wall );\n", "\n", "# DGSEM solver.\n", "# 1) polydeg must be >= the polynomial order set in the HOHQMesh control file to guarantee\n", "# freestream preservation. As a extra task try setting polydeg=3\n", "# 2) VolumeIntegralFluxDifferencing with central volume flux is activated\n", "# for dealiasing\n", "volume_flux = flux_ranocha\n", "solver = DGSEM(polydeg=4, surface_flux=flux_hll,\n", " volume_integral=VolumeIntegralFluxDifferencing(volume_flux))\n", "\n", "# create the unstructured mesh from your mesh file\n", "mesh_file = joinpath(\"out\", \"ice_cream_straight_sides.mesh\")\n", "mesh = UnstructuredMesh2D(mesh_file)\n", "\n", "# Create semidiscretization with all spatial discretization-related components\n", "semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver,\n", " boundary_conditions=boundary_conditions)\n", "\n", "# Create ODE problem from semidiscretization with time span from 0.0 to 2.0\n", "tspan = (0.0, 2.0)\n", "ode = semidiscretize(semi, tspan)\n", "\n", "\n", "# Create the callbacks to output solution files and adapt the time step\n", "summary_callback = SummaryCallback()\n", "save_solution = SaveSolutionCallback(interval=10,\n", " save_initial_solution=true,\n", " save_final_solution=true)\n", "stepsize_callback = StepsizeCallback(cfl=1.0)\n", "\n", "callbacks = CallbackSet(summary_callback, save_solution, stepsize_callback)\n", "\n", "# Evolve ODE problem in time using `solve` from OrdinaryDiffEq\n", "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);\n", "# print the timer summary\n", "summary_callback()" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "Visualization of the solution is carried out in a similar way as above. That is, one converts the `.h5`\n", "output files with `trixi2vtk` and then plot the solution in ParaView. An example plot of the pressure\n", "at the final time is shown below." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "![simulation_straight_sides](https://user-images.githubusercontent.com/25242486/129733926-6ef80676-779b-4f1e-9826-3ebf750cf382.png)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Making a mesh with a curved outer boundary" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Let us modify the mesh from the previous task and place a circular outer boundary instead\n", "of straight-sided outer boundaries.\n", "Note, the \"ice cream cone\" shape is still placed at the center of the domain." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We create the new control file `ice_cream_curved_sides.control` file below and will then highlight the\n", "major differences compared to `ice_cream_straight_sides.control`." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "open(\"out/ice_cream_curved_sides.control\", \"w\") do io\n", " println(io, raw\"\"\"\n", "\\begin{CONTROL_INPUT}\n", " \\begin{RUN_PARAMETERS}\n", " mesh file name = ice_cream_curved_sides.mesh\n", " plot file name = ice_cream_curved_sides.tec\n", " stats file name = none\n", " mesh file format = ISM-v2\n", " polynomial order = 4\n", " plot file format = skeleton\n", " \\end{RUN_PARAMETERS}\n", "\n", " \\begin{BACKGROUND_GRID}\n", " background grid size = [1.0, 1.0, 0.0]\n", " \\end{BACKGROUND_GRID}\n", "\n", " \\begin{SPRING_SMOOTHER}\n", " smoothing = ON\n", " smoothing type = LinearAndCrossBarSpring\n", " number of iterations = 25\n", " \\end{SPRING_SMOOTHER}\n", "\n", "\\end{CONTROL_INPUT}\n", "\n", "\\begin{MODEL}\n", "\n", " \\begin{OUTER_BOUNDARY}\n", " \\begin{PARAMETRIC_EQUATION_CURVE}\n", " name = OuterCircle\n", " xEqn = x(t) = 8.0*sin(2.0*pi*t)\n", " yEqn = y(t) = 8.0*cos(2.0*pi*t)\n", " zEqn = z(t) = 0.0\n", " \\end{PARAMETRIC_EQUATION_CURVE}\n", "\n", " \\end{OUTER_BOUNDARY}\n", "\n", " \\begin{INNER_BOUNDARIES}\n", "\n", " \\begin{CHAIN}\n", " name = IceCreamCone\n", " \\begin{END_POINTS_LINE}\n", " name = LeftSlant\n", " xStart = [-2.0, 1.0, 0.0]\n", " xEnd = [ 0.0, -3.0, 0.0]\n", " \\end{END_POINTS_LINE}\n", "\n", " \\begin{END_POINTS_LINE}\n", " name = RightSlant\n", " xStart = [ 0.0, -3.0, 0.0]\n", " xEnd = [ 2.0, 1.0, 0.0]\n", " \\end{END_POINTS_LINE}\n", "\n", " \\begin{CIRCULAR_ARC}\n", " name = IceCream\n", " units = degrees\n", " center = [ 0.0, 1.0, 0.0]\n", " radius = 2.0\n", " start angle = 0.0\n", " end angle = 180.0\n", " \\end{CIRCULAR_ARC}\n", " \\end{CHAIN}\n", "\n", " \\end{INNER_BOUNDARIES}\n", "\n", "\\end{MODEL}\n", "\\end{FILE}\n", "\"\"\")\n", "end" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "The first alteration is that we have altered the second block of information\n", "`BACKGROUND_GRID` within the `CONTROL_INPUT` to be\n", "```\n", "\\begin{BACKGROUND_GRID}\n", " background grid size = [1.0, 1.0, 0.0]\n", "\\end{BACKGROUND_GRID}\n", "```\n", "This mesh control file has an outer boundary that determines the extent of the domain to be meshed.\n", "Therefore, we only need to supply the `background grid size` to the `BACKGROUND_GRID` control input." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The second alteration is that the `MODEL` now contains information for an `OUTER_BOUNDARY`.\n", "In this case it is a circle of radius `8` centered at `[0.0, 0.0, 0.0]` written as a set of\n", "`PARAMETRIC_EQUATION_CURVE`s.\n", "```\n", " \\begin{OUTER_BOUNDARY}\n", "\n", " \\begin{PARAMETRIC_EQUATION_CURVE}\n", " name = OuterCircle\n", " xEqn = x(t) = 8.0*sin(2.0*pi*t)\n", " yEqn = y(t) = 8.0*cos(2.0*pi*t)\n", " zEqn = z(t) = 0.0\n", " \\end{PARAMETRIC_EQUATION_CURVE}\n", "\n", " \\end{OUTER_BOUNDARY}\n", "```\n", "Just as with the inner boundary curves, we must assign a name to the `OUTER_BOUNDARY`. It will be included\n", "in the generated `.mesh` file and is used within the Trixi.jl elixir file to set boundary conditions." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Again, we create the `.mesh` and `.tec` files with HOHQMesh.jl's function `generate_mesh`" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "control_file = joinpath(\"out\", \"ice_cream_curved_sides.control\")\n", "output = generate_mesh(control_file);" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "The files are placed in the `out` folder." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "The resulting mesh generated by HOHQMesh.jl is given in the following figure." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "![mesh_curved_sides](https://user-images.githubusercontent.com/25242486/129603957-6a92618f-9ed8-4072-b6ab-05533bea746a.png)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Running Trixi.jl on `ice_cream_curved_sides.mesh`" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We can reuse much of the elixir file to setup the uniform flow over an ice cream cone from the\n", "previous part of this tutorial. The only component of the elixir file that must be changed is the boundary condition\n", "dictionary because we now have a boundary named `OuterCircle` instead of four edges of a bounding box." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "# boundary condition dictionary\n", "boundary_conditions = Dict( :OuterCircle => boundary_condition_uniform_flow,\n", " :LeftSlant => boundary_condition_slip_wall,\n", " :RightSlant => boundary_condition_slip_wall,\n", " :IceCream => boundary_condition_slip_wall );" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "Also, we must update the construction of the mesh from our new mesh file `ice_cream_curved_sides.mesh` that\n", "is located in the `out` folder." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "# create the unstructured mesh from your mesh file\n", "mesh_file = joinpath(\"out\", \"ice_cream_curved_sides.mesh\")\n", "mesh = UnstructuredMesh2D(mesh_file);" ], "metadata": {}, "execution_count": null }, { "cell_type": "markdown", "source": [ "We can then post-process the solution file at the final time on the new mesh with `Trixi2Vtk` and visualize with ParaView." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "![simulation_curved_sides](https://user-images.githubusercontent.com/25242486/129733924-778795c1-9119-419a-8b89-bcbe13e33cd7.png)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "## Setting up a simulation with AMR via `P4estMesh`\n", "The above explained mesh file format of `ISM-V2` only works with `UnstructuredMesh2D` and so does\n", "not support AMR. On the other hand, the mesh type `P4estMesh` allows AMR. The mesh\n", "constructor for the `P4estMesh` imports an unstructured, conforming mesh from an Abaqus mesh file\n", "(`.inp`)." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "As described above, the first block of the HOHQMesh control file contains the parameter\n", "`mesh file format`. If you set `mesh file format = ABAQUS` instead of `ISM-V2`,\n", "HOHQMesh.jl's function `generate_mesh` creates an Abaqus mesh file `.inp`.\n", "```julia\n", "using HOHQMesh\n", "control_file = joinpath(\"out\", \"ice_cream_straight_sides.control\")\n", "output = generate_mesh(control_file);\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Now, you can create a `P4estMesh` from your mesh file. It is described in detail in the\n", "[P4est-based mesh](https://trixi-framework.github.io/Trixi.jl/stable/meshes/p4est_mesh/#P4est-based-mesh)\n", "part of the Trixi.jl docs.\n", "```julia\n", "using Trixi\n", "mesh_file = joinpath(\"out\", \"ice_cream_straight_sides.inp\")\n", "mesh = P4estMesh{2}(mesh_file)\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "Since `P4estMesh` supports AMR, we just have to extend the setup from the first example by the\n", "standard AMR procedure. For more information about AMR in Trixi.jl, see the matching tutorial." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "```julia\n", "amr_indicator = IndicatorLöhner(semi, variable=density)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "amr_controller = ControllerThreeLevel(semi, amr_indicator,\n", " base_level=0,\n", " med_level =1, med_threshold=0.05,\n", " max_level =3, max_threshold=0.1)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "amr_callback = AMRCallback(semi, amr_controller,\n", " interval=5,\n", " adapt_initial_condition=true,\n", " adapt_initial_condition_only_refine=true)" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "callbacks = CallbackSet(..., amr_callback)\n", "```" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We can then post-process the solution file at the final time on the new mesh with `Trixi2Vtk` and visualize\n", "with ParaView, see the appropriate [visualization section](https://trixi-framework.github.io/Trixi.jl/stable/visualization/#Trixi2Vtk)\n", "for details." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "![simulation_straight_sides_p4est_amr](https://user-images.githubusercontent.com/74359358/168049930-8abce6ac-cd47-4d04-b40b-0fa459bbd98d.png)" ], "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\", \"Trixi2Vtk\", \"HOHQMesh\"],\n", " mode=PKGMODE_MANIFEST)" ], "metadata": {}, "execution_count": null } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.9.4" }, "kernelspec": { "name": "julia-1.9", "display_name": "Julia 1.9.4", "language": "julia" } }, "nbformat": 4 }