From d53fdb59fe71c1ba956fc758b5874c3914304ddf Mon Sep 17 00:00:00 2001 From: Johannes Markert <10619309+jmark@users.noreply.github.com> Date: Wed, 31 Jan 2024 16:49:35 +0100 Subject: [PATCH] Feature: T8codeMesh backend with MPI support (#1803) * Initial commit for the new feature using t8code as meshing backend. * Delete t8code_2d_dgsem * Added new examples and tests. Testing updates for T8code.jl. * Worked in the comments. * Fixed spelling. * Update src/auxiliary/auxiliary.jl Co-authored-by: Hendrik Ranocha * Added whitespace in Unions. * Adapted commented out code block reporting the no. of elements per level. * Added dummy save mesh support for . * Added test . * Added to method signature. * Deleted unnecessary comments. * Removed commented out tests. * Fixed Morton ordering bug in 2D at mortar interfaces. * Disabled `save_solution` callbacks and added more tests. * Added more tests. * Updated code according to the review. * Update src/auxiliary/t8code.jl Co-authored-by: Hendrik Ranocha * Update src/auxiliary/t8code.jl Co-authored-by: Hendrik Ranocha * Update src/auxiliary/t8code.jl Co-authored-by: Hendrik Ranocha * Update src/auxiliary/t8code.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/solvers/dgsem_t8code/containers_2d.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Code cleanup. * Updated to T8code@0.3.0 * Fixing minor issues. * Fixed typo. * Code cleanup. * Enabled `set_ghost` in examples. * Generalized type info in function signature. * Added namespace qualifier. * Updated comments. * Refactored code and deleted lots of it. * Removed a copy operation. * Initial commit. * Fxinig minor bugs. * Fixed minor typo. * Added first 3d example and fixed segfault. * Added many 3D examples and tests. * Backup. * Fixed merging issues. * Adding more tests. * Fixed some merging issues and formatting. * Fixed spelling. * Fixed spelling and changed assert macro. * Applied automatic formatting. * Applied automatic formatting. * Backup. * Removed superfluous outer constructor for T8codeMesh. * Added return statement for consistency. * Fixed wrong indentation by autoformatter. * Added comments. * Made sure an exception is thrown. * Changed flags for sc_init for t8code initialization. * Updated formatting. * Workaround for error about calling MPI routines after MPI has been finalized. * Upped to T8code v0.4.1. * Added mpi_finailize_hook for proper memory cleanup. * Added t8code to test_threaded.jl * Added a `save_mesh_file` call in order to satisfy code coverage. * Improved finalizer logic for T8coeMesh. * Refined code. * Restructured to do blocks. * Moved save_mesh_file call to test file. * Fixed spelling error. * Made sc_finalize optional. * Fixed spelling. * Cleaned up examples. * Updated and cleaned t8code solver codes. * Updated tests for t8code 3D code. * Fixed spelling. * Update elixir_euler_source_terms_nonconforming_unstructured_curved.jl * Update elixir_euler_source_terms_nonconforming_unstructured_curved.jl * Fixed indentation. * Update src/solvers/dgsem_structured/dg_3d.jl Co-authored-by: Hendrik Ranocha * Update src/solvers/dgsem_t8code/containers_3d.jl Co-authored-by: Andrew Winters * Update src/callbacks_step/amr_dg3d.jl Co-authored-by: Andrew Winters * Update examples/t8code_3d_dgsem/elixir_euler_ec.jl Co-authored-by: Andrew Winters * Update examples/t8code_3d_dgsem/elixir_advection_unstructured_curved.jl Co-authored-by: Andrew Winters * Update examples/t8code_3d_dgsem/elixir_advection_amr_unstructured_curved.jl Co-authored-by: Andrew Winters * Update src/solvers/dgsem_structured/dg_3d.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/callbacks_step/analysis_dg3d.jl Co-authored-by: Hendrik Ranocha * Update examples/t8code_3d_dgsem/elixir_euler_free_stream.jl Co-authored-by: Andrew Winters * Removed NDIMS from T8codeMesh construction in case of p4est/p8est connectivity input. * Aligned T8codeMesh constructur with other mesh constructors. * Update examples/t8code_3d_dgsem/elixir_euler_sedov.jl Co-authored-by: Andrew Winters * Update examples/t8code_3d_dgsem/elixir_euler_sedov.jl Co-authored-by: Andrew Winters * Cleanup up. * Added @allocated test. * Fixed formatting. * Applied formatter. * Added ParallelT8codeMesh to constructor. * Code cleanup. * Extended T8codeMesh backend to support MPI interfaces datastructures. * Update NEWS.md (#1780) * set version to v0.6.5 * set development version to v0.6.6-pre * hotfix: restrict DiffEqBase.jl to let CI pass (#1788) * hotfix: restrict DiffEqBase.jl to let CI pass * restrict DiffEqBase.jl in main Project.toml * Update Project.toml * Parabolic Mortar for AMR `P4estMesh{3}` (#1765) * Clean branch * Un-Comment * un-comment * test coarsen * remove redundancy * Remove support for passive terms * expand resize * comments * format * Avoid code duplication * Update src/callbacks_step/amr_dg1d.jl Co-authored-by: Michael Schlottke-Lakemper * comment * comment & format * Try to increase coverage * Slightly more expressive names * Apply suggestions from code review * add specifier for 1d * Structs for resizing parabolic helpers * check if mortars are present * reuse `reinitialize_containers!` * resize calls for parabolic helpers * update analysis callbacks * Velocities for compr euler * Init container * correct copy-paste error * resize each dim * add dispatch * Add AMR for shear layer * USe only amr shear layer * first steps towards p4est parabolic amr * Add tests * remove plots * Format * remove redundant line * platform independent tests * No need for different flux_viscous comps after adding container_viscous to p4est * Laplace 3d * Longer times to allow converage to hit coarsen! * Increase testing of Laplace 3D * Add tests for velocities * remove comment * add elixir for amr testing * adding commented out mortar routines in 2D * Adding Mortar to 2d dg parabolic term * remove testing snippet * fix comments * add more arguments for dispatch * add some temporary todo notes * some updates for AP and KS * specialize mortar_fluxes_to_elements * BUGFIX: apply_jacobian_parabolic! was incorrect for P4estMesh * fixed rhs_parabolic! for mortars * more changes to elixir * indexing bug * comments * Adding the example for nonperiodic BCs with amr * hopefully this fixes AMR boundaries for parabolic terms * add elixir * Example with non periodic bopundary conditions * remove cruft * 3D parabolic amr * TGV elixir * Creating test for AMR 3D parabolic * Formatting * test formatting * Update src/Trixi.jl * Update src/equations/compressible_euler_1d.jl * Update src/equations/compressible_euler_2d.jl * Update src/equations/compressible_euler_3d.jl * Update src/solvers/dgsem_tree/container_viscous_2d.jl * Update src/solvers/dgsem_tree/container_viscous_2d.jl * Update src/solvers/dgsem_tree/container_viscous_2d.jl * Update src/solvers/dgsem_tree/container_viscous_3d.jl * Update src/solvers/dgsem_tree/container_viscous_3d.jl * Update src/solvers/dgsem_tree/container_viscous_3d.jl * Update src/solvers/dgsem_p4est/dg_2d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update src/solvers/dgsem_tree/dg_2d_parabolic.jl * Update test/test_parabolic_3d.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_tree/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_2d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_2d.jl * Update src/solvers/dgsem_p4est/dg_2d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_2d_parabolic.jl * Update test/test_parabolic_3d.jl * Update test/test_parabolic_3d.jl * Update test/test_parabolic_3d.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_taylor_green_vortex_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_taylor_green_vortex_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_taylor_green_vortex_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_taylor_green_vortex_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_taylor_green_vortex_amr.jl * Update dg_3d_parabolic.jl * Update test_parabolic_3d.jl * Update elixir_navierstokes_3d_blast_wave_amr.jl * Update elixir_navierstokes_taylor_green_vortex_amr.jl * Update dg_3d_parabolic.jl * Update test_parabolic_3d.jl * Update test_parabolic_3d.jl * Update examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Update elixir_navierstokes_3d_blast_wave_amr.jl * Update elixir_navierstokes_taylor_green_vortex_amr.jl * Update dg_3d_parabolic.jl * Update test_parabolic_3d.jl * Delete examples/p4est_3d_dgsem/elixir_navierstokes_3d_blast_wave_amr.jl * Create elixir_navierstokes_blast_wave_amr.jl * Update test_parabolic_3d.jl * Update NEWS.md * Update NEWS.md * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl * Update src/solvers/dgsem_p4est/dg_3d_parabolic.jl --------- Co-authored-by: Daniel_Doehring Co-authored-by: Daniel Doehring Co-authored-by: Michael Schlottke-Lakemper Co-authored-by: Jesse Chan <1156048+jlchan@users.noreply.github.com> Co-authored-by: Jesse Chan * reset the timer also on non-root MPI processes (#1787) Co-authored-by: Michael Schlottke-Lakemper * Add HLLC flux for non-cartesian meshes to CompressibleEulerEquations{2,3}D (#1790) * add HLLC flux for non-cartesian meshes * add tests for HLLC flux * Add 2D test with HLLC * Update test_p4est_3d.jl * Update test_p4est_2d.jl * Update test_p4est_3d.jl * Update src/equations/compressible_euler_3d.jl Co-authored-by: Hendrik Ranocha * Update src/equations/compressible_euler_2d.jl Co-authored-by: Hendrik Ranocha * Update compressible_euler_2d.jl * Update compressible_euler_3d.jl * Update test_p4est_2d.jl * Update test_p4est_3d.jl * Update compressible_euler_2d.jl * Update compressible_euler_2d.jl --------- Co-authored-by: Daniel Doehring Co-authored-by: Hendrik Ranocha * Bump crate-ci/typos from 1.16.23 to 1.16.26 (#1793) Bumps [crate-ci/typos](https://github.com/crate-ci/typos) from 1.16.23 to 1.16.26. - [Release notes](https://github.com/crate-ci/typos/releases) - [Changelog](https://github.com/crate-ci/typos/blob/master/CHANGELOG.md) - [Commits](https://github.com/crate-ci/typos/compare/v1.16.23...v1.16.26) --- updated-dependencies: - dependency-name: crate-ci/typos dependency-type: direct:production update-type: version-update:semver-patch ... Signed-off-by: dependabot[bot] Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> * Extend `CompressibleEulerQuasi1D` and `ShallowWaterQuasi1D` to `DGMulti` (#1797) * adding DGMulti versions of fluxes * remove incorrect factor of 2 * add example and test * formatting * add comment * revert removing factor of 2 * formatting * add SWE quasi-1D test d * enable quasi1D SWE for DGMulti * add docstrings * formatting * Update src/equations/compressible_euler_quasi_1d.jl Co-authored-by: Hendrik Ranocha * adding comments explaining why `normal_direction` is included in 1D * Apply suggestions from code review Co-authored-by: Daniel Doehring --------- Co-authored-by: Hendrik Ranocha Co-authored-by: Daniel Doehring * Fix boundary_condition_slip_wall for SWE (#1798) * fix_wall_bc * add test * apply formatter * adjust some comments * update elixir and test; add topography * comments explaining usage of `ForwardDiff.jacobian` (#1800) * Bump actions/upload-artifact from 3 to 4 (#1795) * Bump actions/upload-artifact from 3 to 4 Bumps [actions/upload-artifact](https://github.com/actions/upload-artifact) from 3 to 4. - [Release notes](https://github.com/actions/upload-artifact/releases) - [Commits](https://github.com/actions/upload-artifact/compare/v3...v4) --- updated-dependencies: - dependency-name: actions/upload-artifact dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] * bump download-artifact to v4 --------- Signed-off-by: dependabot[bot] Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Co-authored-by: Hendrik Ranocha * Code cleanup. * Applied formatter. * Registered mpi t8code tests. * Added load-balancing feature to T8codeMesh. * Applied formatter. * Updated some test values. * Fixed typo. * Removed unused member variable. * Apply suggestions from code review Co-authored-by: Daniel Doehring * suggestions from review * fix format (strange?) * Added comments to help interpreting the source code. * Update src/callbacks_step/amr_dg3d.jl Co-authored-by: Benedict <135045760+bgeihe@users.noreply.github.com> * Adhered to unified mesh constructor calling scheme. * Applied formatter. * Switched to Float64 instead of Cdouble. * Update src/meshes/t8code_mesh.jl Co-authored-by: Daniel Doehring * Refactored negative volume check. * Applied formatter. * Fixed typo resp. bug. * Apply suggestions from code review Co-authored-by: Hendrik Ranocha * add missing allocation checks * Some refactoring. * Deleted msh file. * Fixed a bug. * Code cleanup. * Ignore gmsh files. * Removed adapt! from global namespace. * Added documentation. * Added @test_warn to test. * Applied formatter. * Apply suggestions from code review Co-authored-by: Hendrik Ranocha Co-authored-by: Joshua Lampert <51029046+JoshuaLampert@users.noreply.github.com> * Turned @warn to @info. * Code cleanup and added @deprecated routines in order to avoid breaking release. * Applied formatter. * Added formatter pragmas to avoid ugly formatting. * Applied formatter. * Code cleanup. Documenting. * Applied formatter. * Removed remnants of development. * Removed dubious commented out line. * minors in comments skip ci * Added `retrieve` function which downloads missing mesh file without race conditions in case of MPI environment. * Added Downloads as dependency. * Added missing namespace qualifier. * add compat bound * Apply suggestions from code review Co-authored-by: Daniel Doehring Co-authored-by: Joshua Lampert <51029046+JoshuaLampert@users.noreply.github.com> * use retrieve everywhere for downloading * try to improve formatting * and what julia formatter makes out of it * try to consolidate max_dt computation * format * add default allocation tests * remove left over vtk_write * Revert "try to consolidate max_dt computation" This reverts commit d0e9a2cb2882f426a3966aac03e6ce497e1b720f. * fix variable name * less consolidation * remove unused variables * Update src/auxiliary/auxiliary.jl Co-authored-by: Hendrik Ranocha * Update src/Trixi.jl Co-authored-by: Hendrik Ranocha * revert b068d1183 (manually) * use Trixi.download * Revert "use retrieve everywhere for downloading" This reverts commit a438547e23bc489f12beea625dcc0905a50d18ce. * change to new Trixi.download signature * format * remove merge leftover * format * Code refactor. Applied review comments. * Fixed bug. * Fixed typo. * Update src/auxiliary/auxiliary.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Update src/meshes/t8code_mesh.jl Co-authored-by: Hendrik Ranocha * Applying review comments. * Added return early if there is nothing to do in adapt routines. --------- Signed-off-by: dependabot[bot] Co-authored-by: Johannes Markert Co-authored-by: Hendrik Ranocha Co-authored-by: Andrew Winters Co-authored-by: Jesse Chan <1156048+jlchan@users.noreply.github.com> Co-authored-by: Hendrik Ranocha Co-authored-by: Ahmad Peyvan <115842305+apey236@users.noreply.github.com> Co-authored-by: Daniel_Doehring Co-authored-by: Daniel Doehring Co-authored-by: Michael Schlottke-Lakemper Co-authored-by: Jesse Chan Co-authored-by: Benedict <135045760+bgeihe@users.noreply.github.com> Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Co-authored-by: Patrick Ersing <114223904+patrickersing@users.noreply.github.com> Co-authored-by: Benedict Geihe Co-authored-by: Joshua Lampert <51029046+JoshuaLampert@users.noreply.github.com> --- Project.toml | 2 + .../elixir_2d_euler_vortex_unstructured.jl | 9 +- examples/dgmulti_2d/elixir_euler_hohqmesh.jl | 9 +- .../elixir_advection_amr_unstructured_flag.jl | 7 +- .../elixir_advection_unstructured_flag.jl | 7 +- .../elixir_euler_blast_wave_amr.jl | 7 +- .../elixir_euler_double_mach_amr.jl | 8 +- .../elixir_euler_forward_step_amr.jl | 8 +- .../elixir_euler_free_stream.jl | 7 +- ...e_terms_nonconforming_unstructured_flag.jl | 7 +- .../elixir_euler_supersonic_cylinder.jl | 8 +- .../elixir_euler_wall_bc_amr.jl | 8 +- examples/p4est_2d_dgsem/elixir_mhd_rotor.jl | 7 +- ...lixir_advection_amr_unstructured_curved.jl | 7 +- .../elixir_advection_unstructured_curved.jl | 7 +- examples/p4est_3d_dgsem/elixir_euler_ec.jl | 7 +- .../elixir_euler_free_stream.jl | 7 +- .../elixir_euler_free_stream_extruded.jl | 7 +- ...terms_nonconforming_unstructured_curved.jl | 7 +- ...euler_source_terms_nonperiodic_hohqmesh.jl | 8 +- .../elixir_shallowwater_well_balanced.jl | 1 - .../elixir_advection_amr_unstructured_flag.jl | 7 +- .../elixir_advection_unstructured_flag.jl | 7 +- .../elixir_euler_free_stream.jl | 7 +- ...e_terms_nonconforming_unstructured_flag.jl | 7 +- examples/t8code_2d_dgsem/elixir_mhd_rotor.jl | 7 +- ...lixir_advection_amr_unstructured_curved.jl | 9 +- .../elixir_advection_unstructured_curved.jl | 7 +- examples/t8code_3d_dgsem/elixir_euler_ec.jl | 7 +- .../elixir_euler_free_stream.jl | 7 +- .../elixir_euler_free_stream_extruded.jl | 7 +- ...terms_nonconforming_unstructured_curved.jl | 7 +- .../elixir_acoustics_gauss_wall.jl | 8 +- .../elixir_advection_basic.jl | 9 +- .../elixir_euler_basic.jl | 9 +- .../unstructured_2d_dgsem/elixir_euler_ec.jl | 9 +- .../elixir_euler_free_stream.jl | 9 +- .../elixir_euler_periodic.jl | 9 +- .../elixir_euler_sedov.jl | 7 +- .../elixir_euler_wall_bc.jl | 9 +- .../elixir_mhd_alfven_wave.jl | 8 +- .../unstructured_2d_dgsem/elixir_mhd_ec.jl | 8 +- .../elixir_shallowwater_dirichlet.jl | 8 +- .../elixir_shallowwater_ec.jl | 8 +- .../elixir_shallowwater_ec_shockcapturing.jl | 8 +- .../elixir_shallowwater_source_terms.jl | 8 +- ...ixir_shallowwater_three_mound_dam_break.jl | 13 +- ...lixir_shallowwater_twolayer_convergence.jl | 8 +- .../elixir_shallowwater_twolayer_dam_break.jl | 8 +- ...xir_shallowwater_twolayer_well_balanced.jl | 8 +- ...xir_shallowwater_wall_bc_shockcapturing.jl | 9 +- .../elixir_shallowwater_well_balanced.jl | 9 +- .../elixir_advection_basic.jl | 9 +- .../elixir_euler_free_stream.jl | 9 +- .../elixir_euler_source_terms.jl | 9 +- src/Trixi.jl | 1 + src/auxiliary/auxiliary.jl | 23 + src/auxiliary/t8code.jl | 235 +------ src/callbacks_step/amr.jl | 40 +- src/callbacks_step/amr_dg.jl | 11 +- src/callbacks_step/amr_dg2d.jl | 7 +- src/callbacks_step/amr_dg3d.jl | 7 +- src/callbacks_step/analysis_dg2d_parallel.jl | 6 +- src/callbacks_step/analysis_dg3d_parallel.jl | 6 +- src/callbacks_step/stepsize_dg2d.jl | 32 + src/callbacks_step/stepsize_dg3d.jl | 32 + src/meshes/p4est_mesh.jl | 4 + src/meshes/t8code_mesh.jl | 577 +++++++++++++++++- .../dgsem_p4est/containers_parallel.jl | 6 +- src/solvers/dgsem_p4est/dg_2d_parallel.jl | 21 +- src/solvers/dgsem_p4est/dg_3d_parallel.jl | 23 +- src/solvers/dgsem_p4est/dg_parallel.jl | 9 +- src/solvers/dgsem_t8code/containers.jl | 13 +- src/solvers/dgsem_t8code/containers_2d.jl | 3 +- src/solvers/dgsem_t8code/containers_3d.jl | 3 +- .../dgsem_t8code/containers_parallel.jl | 65 ++ src/solvers/dgsem_t8code/dg.jl | 7 +- src/solvers/dgsem_t8code/dg_parallel.jl | 135 ++++ src/solvers/dgsem_tree/dg_2d_parallel.jl | 3 +- test/test_mpi.jl | 4 +- test/test_mpi_p4est_2d.jl | 63 ++ test/test_mpi_p4est_3d.jl | 81 +++ test/test_mpi_t8code_2d.jl | 142 +++++ test/test_mpi_t8code_3d.jl | 180 ++++++ test/test_t8code_2d.jl | 7 +- 85 files changed, 1529 insertions(+), 640 deletions(-) create mode 100644 src/solvers/dgsem_t8code/containers_parallel.jl create mode 100644 src/solvers/dgsem_t8code/dg_parallel.jl create mode 100644 test/test_mpi_t8code_2d.jl create mode 100644 test/test_mpi_t8code_3d.jl diff --git a/Project.toml b/Project.toml index 0bbdec206d8..e99b08e0e81 100644 --- a/Project.toml +++ b/Project.toml @@ -9,6 +9,7 @@ ConstructionBase = "187b0558-2788-49d3-abe0-74a17ed4e7c9" DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e" DiffEqCallbacks = "459566f4-90b8-5000-8ac3-15dfb0a30def" +Downloads = "f43a241f-c20a-4ad4-852c-f6b1247861c6" EllipsisNotation = "da5c29d0-fa7d-589e-88eb-ea29b0a81949" FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" @@ -57,6 +58,7 @@ ConstructionBase = "1.3" DataStructures = "0.18.15" DiffEqBase = "6 - 6.143" DiffEqCallbacks = "2.25" +Downloads = "1.6" EllipsisNotation = "1.0" FillArrays = "0.13.2, 1" ForwardDiff = "0.10.18" diff --git a/benchmark/elixir_2d_euler_vortex_unstructured.jl b/benchmark/elixir_2d_euler_vortex_unstructured.jl index 082b6648abf..43e4b6559de 100644 --- a/benchmark/elixir_2d_euler_vortex_unstructured.jl +++ b/benchmark/elixir_2d_euler_vortex_unstructured.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -49,11 +48,9 @@ end initial_condition = initial_condition_isentropic_vortex solver = DGSEM(polydeg = 3, surface_flux = flux_lax_friedrichs) -default_mesh_file = joinpath(@__DIR__, "mesh_uniform_cartesian.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/ranocha/f4ea19ba3b62348968c971db43d7798b/raw/a506abb9479c020920cf6068c142670fc1a9aadc/mesh_uniform_cartesian.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/ranocha/f4ea19ba3b62348968c971db43d7798b/raw/a506abb9479c020920cf6068c142670fc1a9aadc/mesh_uniform_cartesian.mesh", + joinpath(@__DIR__, "mesh_uniform_cartesian.mesh")) + mesh = UnstructuredMesh2D(mesh_file, periodicity = true) semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver) diff --git a/examples/dgmulti_2d/elixir_euler_hohqmesh.jl b/examples/dgmulti_2d/elixir_euler_hohqmesh.jl index f534b5bc8ad..9b14a5c6827 100644 --- a/examples/dgmulti_2d/elixir_euler_hohqmesh.jl +++ b/examples/dgmulti_2d/elixir_euler_hohqmesh.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -30,12 +29,8 @@ dg = DGMulti(polydeg = 8, element_type = Quad(), approximation_type = SBP(), ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_trixi_unstructured_mesh_docs.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/52056f1487853fab63b7f4ed7f171c80/raw/9d573387dfdbb8bce2a55db7246f4207663ac07f/mesh_trixi_unstructured_mesh_docs.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/52056f1487853fab63b7f4ed7f171c80/raw/9d573387dfdbb8bce2a55db7246f4207663ac07f/mesh_trixi_unstructured_mesh_docs.mesh", + joinpath(@__DIR__, "mesh_trixi_unstructured_mesh_docs.mesh")) mesh = DGMultiMesh(dg, mesh_file) diff --git a/examples/p4est_2d_dgsem/elixir_advection_amr_unstructured_flag.jl b/examples/p4est_2d_dgsem/elixir_advection_amr_unstructured_flag.jl index 0a50b3644f0..4bfb2d3e375 100644 --- a/examples/p4est_2d_dgsem/elixir_advection_amr_unstructured_flag.jl +++ b/examples/p4est_2d_dgsem/elixir_advection_amr_unstructured_flag.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -31,10 +30,8 @@ Trixi.validate_faces(faces) mapping_flag = Trixi.transfinite_mapping(faces) # Unstructured mesh with 24 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) mesh = P4estMesh{2}(mesh_file, polydeg = 3, mapping = mapping_flag, diff --git a/examples/p4est_2d_dgsem/elixir_advection_unstructured_flag.jl b/examples/p4est_2d_dgsem/elixir_advection_unstructured_flag.jl index 37fcc547f60..1ab96925fe6 100644 --- a/examples/p4est_2d_dgsem/elixir_advection_unstructured_flag.jl +++ b/examples/p4est_2d_dgsem/elixir_advection_unstructured_flag.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -29,10 +28,8 @@ Trixi.validate_faces(faces) mapping_flag = Trixi.transfinite_mapping(faces) # Unstructured mesh with 24 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) mesh = P4estMesh{2}(mesh_file, polydeg = 3, mapping = mapping_flag, diff --git a/examples/p4est_2d_dgsem/elixir_euler_blast_wave_amr.jl b/examples/p4est_2d_dgsem/elixir_euler_blast_wave_amr.jl index 0ca4fdc2eb7..5db5f74a686 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_blast_wave_amr.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_blast_wave_amr.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -50,10 +49,8 @@ volume_integral = VolumeIntegralShockCapturingHG(indicator_sc; solver = DGSEM(basis, surface_flux, volume_integral) # Unstructured mesh with 48 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/a075f8ec39a67fa9fad8f6f84342cbca/raw/a7206a02ed3a5d3cadacd8d9694ac154f9151db7/square_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/a075f8ec39a67fa9fad8f6f84342cbca/raw/a7206a02ed3a5d3cadacd8d9694ac154f9151db7/square_unstructured_1.inp", + joinpath(@__DIR__, "square_unstructured_1.inp")) mesh = P4estMesh{2}(mesh_file, polydeg = 3, initial_refinement_level = 1) diff --git a/examples/p4est_2d_dgsem/elixir_euler_double_mach_amr.jl b/examples/p4est_2d_dgsem/elixir_euler_double_mach_amr.jl index 92928146d7b..fbc11e89185 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_double_mach_amr.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_double_mach_amr.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -99,11 +98,8 @@ solver = DGSEM(polydeg = polydeg, surface_flux = surface_flux, volume_integral = volume_integral) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "abaqus_double_mach.inp") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/a0806ef0d03cf5ea221af523167b6e32/raw/61ed0eb017eb432d996ed119a52fb041fe363e8c/abaqus_double_mach.inp", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/a0806ef0d03cf5ea221af523167b6e32/raw/61ed0eb017eb432d996ed119a52fb041fe363e8c/abaqus_double_mach.inp", + joinpath(@__DIR__, "abaqus_double_mach.inp")) mesh = P4estMesh{2}(mesh_file) diff --git a/examples/p4est_2d_dgsem/elixir_euler_forward_step_amr.jl b/examples/p4est_2d_dgsem/elixir_euler_forward_step_amr.jl index 0ec9fc222f2..654efd5e209 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_forward_step_amr.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_forward_step_amr.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -104,11 +103,8 @@ solver = DGSEM(polydeg = polydeg, surface_flux = surface_flux, volume_integral = volume_integral) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "abaqus_forward_step.inp") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/b346ee6aa5446687f128eab8b37d52a7/raw/cd1e1d43bebd8d2631a07caec45585ec8456ca4c/abaqus_forward_step.inp", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/b346ee6aa5446687f128eab8b37d52a7/raw/cd1e1d43bebd8d2631a07caec45585ec8456ca4c/abaqus_forward_step.inp", + joinpath(@__DIR__, "abaqus_forward_step.inp")) mesh = P4estMesh{2}(mesh_file) diff --git a/examples/p4est_2d_dgsem/elixir_euler_free_stream.jl b/examples/p4est_2d_dgsem/elixir_euler_free_stream.jl index 38307a7d781..ab11dc11567 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_free_stream.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_free_stream.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -31,10 +30,8 @@ end # Get the uncurved mesh from a file (downloads the file if not available locally) # Unstructured mesh with 48 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/a075f8ec39a67fa9fad8f6f84342cbca/raw/a7206a02ed3a5d3cadacd8d9694ac154f9151db7/square_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/a075f8ec39a67fa9fad8f6f84342cbca/raw/a7206a02ed3a5d3cadacd8d9694ac154f9151db7/square_unstructured_1.inp", + joinpath(@__DIR__, "square_unstructured_1.inp")) # Map the unstructured mesh with the mapping above mesh = P4estMesh{2}(mesh_file, polydeg = 3, mapping = mapping, initial_refinement_level = 1) diff --git a/examples/p4est_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl b/examples/p4est_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl index 09d018309a6..084fd699b8e 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -31,10 +30,8 @@ mapping_flag = Trixi.transfinite_mapping(faces) # Get the uncurved mesh from a file (downloads the file if not available locally) # Unstructured mesh with 24 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) mesh = P4estMesh{2}(mesh_file, polydeg = 3, mapping = mapping_flag, diff --git a/examples/p4est_2d_dgsem/elixir_euler_supersonic_cylinder.jl b/examples/p4est_2d_dgsem/elixir_euler_supersonic_cylinder.jl index 36c5624ba97..76ee96d4766 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_supersonic_cylinder.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_supersonic_cylinder.jl @@ -13,7 +13,6 @@ # # Keywords: supersonic flow, shock capturing, AMR, unstructured curved mesh, positivity preservation, compressible Euler, 2D -using Downloads: download using OrdinaryDiffEq using Trixi @@ -82,11 +81,8 @@ solver = DGSEM(polydeg = polydeg, surface_flux = surface_flux, volume_integral = volume_integral) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "abaqus_cylinder_in_channel.inp") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/a08f78f6b185b63c3baeff911a63f628/raw/addac716ea0541f588b9d2bd3f92f643eb27b88f/abaqus_cylinder_in_channel.inp", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/a08f78f6b185b63c3baeff911a63f628/raw/addac716ea0541f588b9d2bd3f92f643eb27b88f/abaqus_cylinder_in_channel.inp", + joinpath(@__DIR__, "abaqus_cylinder_in_channel.inp")) mesh = P4estMesh{2}(mesh_file) diff --git a/examples/p4est_2d_dgsem/elixir_euler_wall_bc_amr.jl b/examples/p4est_2d_dgsem/elixir_euler_wall_bc_amr.jl index 8b8d05bade8..75e60d0c78b 100644 --- a/examples/p4est_2d_dgsem/elixir_euler_wall_bc_amr.jl +++ b/examples/p4est_2d_dgsem/elixir_euler_wall_bc_amr.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -39,11 +38,8 @@ solver = DGSEM(polydeg = 5, surface_flux = flux_lax_friedrichs, volume_integral = VolumeIntegralFluxDifferencing(volume_flux)) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "abaqus_gingerbread_man.inp") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/0e9e990a04b5105d1d2e3096a6e41272/raw/0d924b1d7e7d3cc1070a6cc22fe1d501687aa6dd/abaqus_gingerbread_man.inp", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/0e9e990a04b5105d1d2e3096a6e41272/raw/0d924b1d7e7d3cc1070a6cc22fe1d501687aa6dd/abaqus_gingerbread_man.inp", + joinpath(@__DIR__, "abaqus_gingerbread_man.inp")) mesh = P4estMesh{2}(mesh_file) diff --git a/examples/p4est_2d_dgsem/elixir_mhd_rotor.jl b/examples/p4est_2d_dgsem/elixir_mhd_rotor.jl index 380db487356..089e82580c9 100644 --- a/examples/p4est_2d_dgsem/elixir_mhd_rotor.jl +++ b/examples/p4est_2d_dgsem/elixir_mhd_rotor.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -69,10 +68,8 @@ function mapping_twist(xi, eta) return SVector(x, y) end -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) mesh = P4estMesh{2}(mesh_file, polydeg = 4, diff --git a/examples/p4est_3d_dgsem/elixir_advection_amr_unstructured_curved.jl b/examples/p4est_3d_dgsem/elixir_advection_amr_unstructured_curved.jl index cd280cf5bf6..33afd2e030e 100644 --- a/examples/p4est_3d_dgsem/elixir_advection_amr_unstructured_curved.jl +++ b/examples/p4est_3d_dgsem/elixir_advection_amr_unstructured_curved.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -49,10 +48,8 @@ function mapping(xi, eta, zeta) end # Unstructured mesh with 48 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", + joinpath(@__DIR__, "cube_unstructured_2.inp")) # Mesh polydeg of 2 (half the solver polydeg) to ensure FSP (see above). mesh = P4estMesh{3}(mesh_file, polydeg = 2, diff --git a/examples/p4est_3d_dgsem/elixir_advection_unstructured_curved.jl b/examples/p4est_3d_dgsem/elixir_advection_unstructured_curved.jl index 6df9ac0b16a..83adcbf6a63 100644 --- a/examples/p4est_3d_dgsem/elixir_advection_unstructured_curved.jl +++ b/examples/p4est_3d_dgsem/elixir_advection_unstructured_curved.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -46,10 +45,8 @@ function mapping(xi, eta, zeta) end # Unstructured mesh with 68 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", + joinpath(@__DIR__, "cube_unstructured_1.inp")) mesh = P4estMesh{3}(mesh_file, polydeg = 3, mapping = mapping, diff --git a/examples/p4est_3d_dgsem/elixir_euler_ec.jl b/examples/p4est_3d_dgsem/elixir_euler_ec.jl index d9d774a7ffc..91698545052 100644 --- a/examples/p4est_3d_dgsem/elixir_euler_ec.jl +++ b/examples/p4est_3d_dgsem/elixir_euler_ec.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -46,10 +45,8 @@ function mapping(xi_, eta_, zeta_) end # Unstructured mesh with 48 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", + joinpath(@__DIR__, "cube_unstructured_2.inp")) mesh = P4estMesh{3}(mesh_file, polydeg = 5, mapping = mapping, diff --git a/examples/p4est_3d_dgsem/elixir_euler_free_stream.jl b/examples/p4est_3d_dgsem/elixir_euler_free_stream.jl index 24a781ca59e..6406a38186b 100644 --- a/examples/p4est_3d_dgsem/elixir_euler_free_stream.jl +++ b/examples/p4est_3d_dgsem/elixir_euler_free_stream.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -47,10 +46,8 @@ function mapping(xi_, eta_, zeta_) end # Unstructured mesh with 68 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", + joinpath(@__DIR__, "cube_unstructured_1.inp")) # Mesh polydeg of 2 (half the solver polydeg) to ensure FSP (see above). mesh = P4estMesh{3}(mesh_file, polydeg = 2, diff --git a/examples/p4est_3d_dgsem/elixir_euler_free_stream_extruded.jl b/examples/p4est_3d_dgsem/elixir_euler_free_stream_extruded.jl index f56fe3a429d..08307a449a7 100644 --- a/examples/p4est_3d_dgsem/elixir_euler_free_stream_extruded.jl +++ b/examples/p4est_3d_dgsem/elixir_euler_free_stream_extruded.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -36,10 +35,8 @@ function mapping(xi, eta_, zeta_) end # Unstructured mesh with 48 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", + joinpath(@__DIR__, "cube_unstructured_2.inp")) mesh = P4estMesh{3}(mesh_file, polydeg = 3, mapping = mapping, diff --git a/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl b/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl index 0de22eaea40..e7ca0cad4ba 100644 --- a/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl +++ b/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -49,10 +48,8 @@ function mapping(xi, eta, zeta) end # Unstructured mesh with 68 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", + joinpath(@__DIR__, "cube_unstructured_1.inp")) # Mesh polydeg of 2 (half the solver polydeg) to ensure FSP (see above). mesh = P4estMesh{3}(mesh_file, polydeg = 2, diff --git a/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonperiodic_hohqmesh.jl b/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonperiodic_hohqmesh.jl index 0fa3a28fe8b..7d81d6739bf 100644 --- a/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonperiodic_hohqmesh.jl +++ b/examples/p4est_3d_dgsem/elixir_euler_source_terms_nonperiodic_hohqmesh.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -19,11 +18,8 @@ boundary_conditions = Dict(:Bottom => boundary_condition, solver = DGSEM(polydeg = 4, surface_flux = flux_lax_friedrichs) # Unstructured 3D half circle mesh from HOHQMesh -default_mesh_file = joinpath(@__DIR__, "abaqus_half_circle_3d.inp") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/11461efbfb02c42e06aca338b3d0b645/raw/81deeb1ebc4945952c30af5bb75fe222a18d975c/abaqus_half_circle_3d.inp", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/11461efbfb02c42e06aca338b3d0b645/raw/81deeb1ebc4945952c30af5bb75fe222a18d975c/abaqus_half_circle_3d.inp", + joinpath(@__DIR__, "abaqus_half_circle_3d.inp")) mesh = P4estMesh{3}(mesh_file, initial_refinement_level = 0) diff --git a/examples/structured_2d_dgsem/elixir_shallowwater_well_balanced.jl b/examples/structured_2d_dgsem/elixir_shallowwater_well_balanced.jl index 61dd252fd83..a6a56aa807c 100644 --- a/examples/structured_2d_dgsem/elixir_shallowwater_well_balanced.jl +++ b/examples/structured_2d_dgsem/elixir_shallowwater_well_balanced.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi diff --git a/examples/t8code_2d_dgsem/elixir_advection_amr_unstructured_flag.jl b/examples/t8code_2d_dgsem/elixir_advection_amr_unstructured_flag.jl index f285d24fc6c..0923e328487 100644 --- a/examples/t8code_2d_dgsem/elixir_advection_amr_unstructured_flag.jl +++ b/examples/t8code_2d_dgsem/elixir_advection_amr_unstructured_flag.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -31,10 +30,8 @@ Trixi.validate_faces(faces) mapping_flag = Trixi.transfinite_mapping(faces) # Unstructured mesh with 24 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/t8code_2d_dgsem/elixir_advection_unstructured_flag.jl b/examples/t8code_2d_dgsem/elixir_advection_unstructured_flag.jl index 5ba1ab15489..ba8f1b59b80 100644 --- a/examples/t8code_2d_dgsem/elixir_advection_unstructured_flag.jl +++ b/examples/t8code_2d_dgsem/elixir_advection_unstructured_flag.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -28,10 +27,8 @@ Trixi.validate_faces(faces) mapping_flag = Trixi.transfinite_mapping(faces) # Unstructured mesh with 24 cells of the square domain [-1, 1]^n. -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/t8code_2d_dgsem/elixir_euler_free_stream.jl b/examples/t8code_2d_dgsem/elixir_euler_free_stream.jl index 37d15f38566..5e6c4193c50 100644 --- a/examples/t8code_2d_dgsem/elixir_euler_free_stream.jl +++ b/examples/t8code_2d_dgsem/elixir_euler_free_stream.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -30,10 +29,8 @@ end # Get the uncurved mesh from a file (downloads the file if not available locally) # Unstructured mesh with 48 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/a075f8ec39a67fa9fad8f6f84342cbca/raw/a7206a02ed3a5d3cadacd8d9694ac154f9151db7/square_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/a075f8ec39a67fa9fad8f6f84342cbca/raw/a7206a02ed3a5d3cadacd8d9694ac154f9151db7/square_unstructured_1.inp", + joinpath(@__DIR__, "square_unstructured_1.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/t8code_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl b/examples/t8code_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl index bcc1abc560e..e496eb76729 100644 --- a/examples/t8code_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl +++ b/examples/t8code_2d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_flag.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -30,10 +29,8 @@ mapping_flag = Trixi.transfinite_mapping(faces) # Get the uncurved mesh from a file (downloads the file if not available locally) # Unstructured mesh with 24 cells of the square domain [-1, 1]^n -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/t8code_2d_dgsem/elixir_mhd_rotor.jl b/examples/t8code_2d_dgsem/elixir_mhd_rotor.jl index adb154948fb..ff2e40ae607 100644 --- a/examples/t8code_2d_dgsem/elixir_mhd_rotor.jl +++ b/examples/t8code_2d_dgsem/elixir_mhd_rotor.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -68,10 +67,8 @@ function mapping_twist(xi, eta) return SVector(x, y) end -mesh_file = joinpath(@__DIR__, "square_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/63ff2ea224409e55ee8423b3a33e316a/raw/7db58af7446d1479753ae718930741c47a3b79b7/square_unstructured_2.inp", + joinpath(@__DIR__, "square_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/t8code_3d_dgsem/elixir_advection_amr_unstructured_curved.jl b/examples/t8code_3d_dgsem/elixir_advection_amr_unstructured_curved.jl index 617736afbdd..e7c0f4b7318 100644 --- a/examples/t8code_3d_dgsem/elixir_advection_amr_unstructured_curved.jl +++ b/examples/t8code_3d_dgsem/elixir_advection_amr_unstructured_curved.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -47,11 +46,9 @@ function mapping(xi, eta, zeta) return SVector(5 * x, 5 * y, 5 * z) end -# Unstructured mesh with 48 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", - mesh_file) +# Unstructured mesh with 48 cells of the cube domain [-1, 1]^3. +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", + joinpath(@__DIR__, "cube_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connectivity object first from which diff --git a/examples/t8code_3d_dgsem/elixir_advection_unstructured_curved.jl b/examples/t8code_3d_dgsem/elixir_advection_unstructured_curved.jl index df358435c9a..ee27ee117fe 100644 --- a/examples/t8code_3d_dgsem/elixir_advection_unstructured_curved.jl +++ b/examples/t8code_3d_dgsem/elixir_advection_unstructured_curved.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -45,10 +44,8 @@ function mapping(xi, eta, zeta) end # Unstructured mesh with 68 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", + joinpath(@__DIR__, "cube_unstructured_1.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connectivity object first from which diff --git a/examples/t8code_3d_dgsem/elixir_euler_ec.jl b/examples/t8code_3d_dgsem/elixir_euler_ec.jl index 07745c3ac56..b720bfcd375 100644 --- a/examples/t8code_3d_dgsem/elixir_euler_ec.jl +++ b/examples/t8code_3d_dgsem/elixir_euler_ec.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -45,10 +44,8 @@ function mapping(xi_, eta_, zeta_) end # Unstructured mesh with 48 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", + joinpath(@__DIR__, "cube_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connectivity object first from which diff --git a/examples/t8code_3d_dgsem/elixir_euler_free_stream.jl b/examples/t8code_3d_dgsem/elixir_euler_free_stream.jl index e135d464810..b70a6091adf 100644 --- a/examples/t8code_3d_dgsem/elixir_euler_free_stream.jl +++ b/examples/t8code_3d_dgsem/elixir_euler_free_stream.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -46,10 +45,8 @@ function mapping(xi_, eta_, zeta_) end # Unstructured mesh with 68 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", + joinpath(@__DIR__, "cube_unstructured_1.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connectivity object first from which diff --git a/examples/t8code_3d_dgsem/elixir_euler_free_stream_extruded.jl b/examples/t8code_3d_dgsem/elixir_euler_free_stream_extruded.jl index d129b59826e..6ae38d20b5a 100644 --- a/examples/t8code_3d_dgsem/elixir_euler_free_stream_extruded.jl +++ b/examples/t8code_3d_dgsem/elixir_euler_free_stream_extruded.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -35,10 +34,8 @@ function mapping(xi, eta_, zeta_) end # Unstructured mesh with 48 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_2.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/b8df0033798e4926dec515fc045e8c2c/raw/b9254cde1d1fb64b6acc8416bc5ccdd77a240227/cube_unstructured_2.inp", + joinpath(@__DIR__, "cube_unstructured_2.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/t8code_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl b/examples/t8code_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl index d4664522bea..6856be36ea1 100644 --- a/examples/t8code_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl +++ b/examples/t8code_3d_dgsem/elixir_euler_source_terms_nonconforming_unstructured_curved.jl @@ -1,4 +1,3 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -48,10 +47,8 @@ function mapping(xi, eta, zeta) end # Unstructured mesh with 68 cells of the cube domain [-1, 1]^3 -mesh_file = joinpath(@__DIR__, "cube_unstructured_1.inp") -isfile(mesh_file) || - download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", - mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/efaulhaber/d45c8ac1e248618885fa7cc31a50ab40/raw/37fba24890ab37cfa49c39eae98b44faf4502882/cube_unstructured_1.inp", + joinpath(@__DIR__, "cube_unstructured_1.inp")) # INP mesh files are only support by p4est. Hence, we # create a p4est connecvity object first from which diff --git a/examples/unstructured_2d_dgsem/elixir_acoustics_gauss_wall.jl b/examples/unstructured_2d_dgsem/elixir_acoustics_gauss_wall.jl index 8f8e867dca8..9741430d11c 100644 --- a/examples/unstructured_2d_dgsem/elixir_acoustics_gauss_wall.jl +++ b/examples/unstructured_2d_dgsem/elixir_acoustics_gauss_wall.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -14,11 +13,8 @@ equations = AcousticPerturbationEquations2D(v_mean_global = (0.0, -0.5), solver = DGSEM(polydeg = 4, surface_flux = flux_lax_friedrichs) # Create unstructured quadrilateral mesh from a file -default_mesh_file = joinpath(@__DIR__, "mesh_five_circles_in_circle.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/3c79baad6b4d73bb26ec6420b5d16f45/raw/22aefc4ec2107cf0bffc40e81dfbc52240c625b1/mesh_five_circles_in_circle.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/3c79baad6b4d73bb26ec6420b5d16f45/raw/22aefc4ec2107cf0bffc40e81dfbc52240c625b1/mesh_five_circles_in_circle.mesh", + joinpath(@__DIR__, "mesh_five_circles_in_circle.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_dgsem/elixir_advection_basic.jl b/examples/unstructured_2d_dgsem/elixir_advection_basic.jl index afef6c2c38f..c0ee453344d 100644 --- a/examples/unstructured_2d_dgsem/elixir_advection_basic.jl +++ b/examples/unstructured_2d_dgsem/elixir_advection_basic.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -16,12 +15,8 @@ solver = DGSEM(polydeg = 6, surface_flux = flux_lax_friedrichs) ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", + joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_euler_basic.jl b/examples/unstructured_2d_dgsem/elixir_euler_basic.jl index cd6a1995757..f8976120d53 100644 --- a/examples/unstructured_2d_dgsem/elixir_euler_basic.jl +++ b/examples/unstructured_2d_dgsem/elixir_euler_basic.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -25,12 +24,8 @@ solver = DGSEM(polydeg = 8, surface_flux = flux_lax_friedrichs) ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_trixi_unstructured_mesh_docs.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/52056f1487853fab63b7f4ed7f171c80/raw/9d573387dfdbb8bce2a55db7246f4207663ac07f/mesh_trixi_unstructured_mesh_docs.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/52056f1487853fab63b7f4ed7f171c80/raw/9d573387dfdbb8bce2a55db7246f4207663ac07f/mesh_trixi_unstructured_mesh_docs.mesh", + joinpath(@__DIR__, "mesh_trixi_unstructured_mesh_docs.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_dgsem/elixir_euler_ec.jl b/examples/unstructured_2d_dgsem/elixir_euler_ec.jl index 0f53aa62a18..58b4d9a1dd2 100644 --- a/examples/unstructured_2d_dgsem/elixir_euler_ec.jl +++ b/examples/unstructured_2d_dgsem/elixir_euler_ec.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -19,12 +18,8 @@ solver = DGSEM(polydeg = 6, surface_flux = flux_ranocha, ############################################################################### # Get the curved quad mesh from a file - -default_mesh_file = joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", + joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_euler_free_stream.jl b/examples/unstructured_2d_dgsem/elixir_euler_free_stream.jl index a2fec1a320a..f266a3de0b2 100644 --- a/examples/unstructured_2d_dgsem/elixir_euler_free_stream.jl +++ b/examples/unstructured_2d_dgsem/elixir_euler_free_stream.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -26,12 +25,8 @@ solver = DGSEM(polydeg = 6, surface_flux = flux_hll) ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_gingerbread_man.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/2c6440b5f8a57db131061ad7aa78ee2b/raw/1f89fdf2c874ff678c78afb6fe8dc784bdfd421f/mesh_gingerbread_man.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/2c6440b5f8a57db131061ad7aa78ee2b/raw/1f89fdf2c874ff678c78afb6fe8dc784bdfd421f/mesh_gingerbread_man.mesh", + joinpath(@__DIR__, "mesh_gingerbread_man.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_dgsem/elixir_euler_periodic.jl b/examples/unstructured_2d_dgsem/elixir_euler_periodic.jl index afd177f0740..e640001ad7f 100644 --- a/examples/unstructured_2d_dgsem/elixir_euler_periodic.jl +++ b/examples/unstructured_2d_dgsem/elixir_euler_periodic.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -19,12 +18,8 @@ solver = DGSEM(polydeg = 6, surface_flux = FluxRotated(flux_hll)) ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", + joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_euler_sedov.jl b/examples/unstructured_2d_dgsem/elixir_euler_sedov.jl index e1cc0932969..06053273b74 100644 --- a/examples/unstructured_2d_dgsem/elixir_euler_sedov.jl +++ b/examples/unstructured_2d_dgsem/elixir_euler_sedov.jl @@ -55,11 +55,8 @@ solver = DGSEM(polydeg = polydeg, surface_flux = surface_flux, volume_integral = volume_integral) # Get the curved quad mesh from a file -default_mesh_file = joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", + joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_euler_wall_bc.jl b/examples/unstructured_2d_dgsem/elixir_euler_wall_bc.jl index b2abefe7eeb..65e5eb51ce6 100644 --- a/examples/unstructured_2d_dgsem/elixir_euler_wall_bc.jl +++ b/examples/unstructured_2d_dgsem/elixir_euler_wall_bc.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -40,12 +39,8 @@ solver = DGSEM(polydeg = 4, surface_flux = flux_hll) ############################################################################### # Get the curved quad mesh from a file - -default_mesh_file = joinpath(@__DIR__, "mesh_box_around_circle.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8b9b11a1eedfa54b215c122c3d17b271/raw/0d2b5d98c87e67a6f384693a8b8e54b4c9fcbf3d/mesh_box_around_circle.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8b9b11a1eedfa54b215c122c3d17b271/raw/0d2b5d98c87e67a6f384693a8b8e54b4c9fcbf3d/mesh_box_around_circle.mesh", + joinpath(@__DIR__, "mesh_box_around_circle.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_dgsem/elixir_mhd_alfven_wave.jl b/examples/unstructured_2d_dgsem/elixir_mhd_alfven_wave.jl index 3ed3e828ca8..0c7152a6ea0 100644 --- a/examples/unstructured_2d_dgsem/elixir_mhd_alfven_wave.jl +++ b/examples/unstructured_2d_dgsem/elixir_mhd_alfven_wave.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -17,12 +16,9 @@ solver = DGSEM(polydeg = 7, volume_integral = VolumeIntegralFluxDifferencing(volume_flux)) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) -mesh_file = default_mesh_file mesh = UnstructuredMesh2D(mesh_file, periodicity = true) semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver) diff --git a/examples/unstructured_2d_dgsem/elixir_mhd_ec.jl b/examples/unstructured_2d_dgsem/elixir_mhd_ec.jl index a40f92cac02..805934e305d 100644 --- a/examples/unstructured_2d_dgsem/elixir_mhd_ec.jl +++ b/examples/unstructured_2d_dgsem/elixir_mhd_ec.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -35,11 +34,8 @@ solver = DGSEM(polydeg = 6, volume_integral = VolumeIntegralFluxDifferencing(volume_flux)) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_dirichlet.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_dirichlet.jl index 1148f25fae3..df1a69192ce 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_dirichlet.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_dirichlet.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -38,11 +37,8 @@ solver = DGSEM(polydeg = 4, surface_flux = (flux_hll, flux_nonconservative_fjord # This setup is for the curved, split form well-balancedness testing # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_outer_circle.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/9beddd9cd00e2a0a15865129eeb24928/raw/be71e67fa48bc4e1e97f5f6cd77c3ed34c6ba9be/mesh_outer_circle.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/9beddd9cd00e2a0a15865129eeb24928/raw/be71e67fa48bc4e1e97f5f6cd77c3ed34c6ba9be/mesh_outer_circle.mesh", + joinpath(@__DIR__, "mesh_outer_circle.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_ec.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_ec.jl index 8e9d396d826..9122fb8287d 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_ec.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_ec.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -25,11 +24,8 @@ solver = DGSEM(polydeg = 6, # This setup is for the curved, split form entropy conservation testing (needs periodic BCs) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_ec_shockcapturing.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_ec_shockcapturing.jl index 94202b81df0..98408db5a78 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_ec_shockcapturing.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_ec_shockcapturing.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -34,11 +33,8 @@ solver = DGSEM(basis, surface_flux, volume_integral) # This setup is for the curved, split form entropy conservation testing (needs periodic BCs) # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_source_terms.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_source_terms.jl index 07668688406..a7aa5808955 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_source_terms.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_source_terms.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -23,11 +22,8 @@ solver = DGSEM(polydeg = 6, surface_flux = surface_flux, # This setup is for the curved, split form convergence test on a periodic domain # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_three_mound_dam_break.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_three_mound_dam_break.jl index 6164f9d4a55..df321aad267 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_three_mound_dam_break.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_three_mound_dam_break.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -93,16 +92,10 @@ solver = DGSEM(basis, surface_flux, volume_integral) ############################################################################### # Get the unstructured quad mesh from a file (downloads the file if not available locally) +mesh_file = Trixi.download("https://gist.githubusercontent.com/svengoldberg/c3c87fecb3fc6e46be7f0d1c7cb35f83/raw/e817ecd9e6c4686581d63c46128f9b6468d396d3/mesh_three_mound.mesh", + joinpath(@__DIR__, "mesh_three_mound.mesh")) -default_meshfile = joinpath(@__DIR__, "mesh_three_mound.mesh") - -isfile(default_meshfile) || - download("https://gist.githubusercontent.com/svengoldberg/c3c87fecb3fc6e46be7f0d1c7cb35f83/raw/e817ecd9e6c4686581d63c46128f9b6468d396d3/mesh_three_mound.mesh", - default_meshfile) - -meshfile = default_meshfile - -mesh = UnstructuredMesh2D(meshfile) +mesh = UnstructuredMesh2D(mesh_file) # Create the semi discretization object semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, solver; diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_convergence.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_convergence.jl index 0b86095663a..fcc08b6f991 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_convergence.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_convergence.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -24,11 +23,8 @@ solver = DGSEM(polydeg = 6, surface_flux = surface_flux, # This setup is for the curved, split form convergence test on a periodic domain # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_dam_break.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_dam_break.jl index 4ad5f7e3201..821f31c52ac 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_dam_break.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_dam_break.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -48,11 +47,8 @@ solver = DGSEM(polydeg = 6, surface_flux = surface_flux, ############################################################################### # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = false) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_well_balanced.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_well_balanced.jl index 6a727df2502..ca1f54595bb 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_well_balanced.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_twolayer_well_balanced.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -42,11 +41,8 @@ solver = DGSEM(polydeg = 6, surface_flux = surface_flux, # This setup is for the curved, split form well-balancedness testing # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_wall_bc_shockcapturing.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_wall_bc_shockcapturing.jl index 76b9642d595..f115113ed27 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_wall_bc_shockcapturing.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_wall_bc_shockcapturing.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -55,12 +54,8 @@ solver = DGSEM(polydeg = polydeg, surface_flux = surface_flux, ############################################################################### # Get the unstructured quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_outer_circle.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/9beddd9cd00e2a0a15865129eeb24928/raw/be71e67fa48bc4e1e97f5f6cd77c3ed34c6ba9be/mesh_outer_circle.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/9beddd9cd00e2a0a15865129eeb24928/raw/be71e67fa48bc4e1e97f5f6cd77c3ed34c6ba9be/mesh_outer_circle.mesh", + joinpath(@__DIR__, "mesh_outer_circle.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_dgsem/elixir_shallowwater_well_balanced.jl b/examples/unstructured_2d_dgsem/elixir_shallowwater_well_balanced.jl index bf4d0be682a..6bad3a77f03 100644 --- a/examples/unstructured_2d_dgsem/elixir_shallowwater_well_balanced.jl +++ b/examples/unstructured_2d_dgsem/elixir_shallowwater_well_balanced.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -36,13 +35,9 @@ solver = DGSEM(polydeg = 6, surface_flux = surface_flux, ############################################################################### # This setup is for the curved, split form well-balancedness testing - # Get the unstructured quad mesh from a file (downloads the file if not available locally) -default_mesh_file = joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/8f8cd23df27fcd494553f2a89f3c1ba4/raw/85e3c8d976bbe57ca3d559d653087b0889535295/mesh_alfven_wave_with_twist_and_flip.mesh", + joinpath(@__DIR__, "mesh_alfven_wave_with_twist_and_flip.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_fdsbp/elixir_advection_basic.jl b/examples/unstructured_2d_fdsbp/elixir_advection_basic.jl index c181203e7a4..fe7e708f3b3 100644 --- a/examples/unstructured_2d_fdsbp/elixir_advection_basic.jl +++ b/examples/unstructured_2d_fdsbp/elixir_advection_basic.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -21,12 +20,8 @@ solver = FDSBP(D_SBP, ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", + joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/examples/unstructured_2d_fdsbp/elixir_euler_free_stream.jl b/examples/unstructured_2d_fdsbp/elixir_euler_free_stream.jl index 7ada50c0c65..25a81c16bf9 100644 --- a/examples/unstructured_2d_fdsbp/elixir_euler_free_stream.jl +++ b/examples/unstructured_2d_fdsbp/elixir_euler_free_stream.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -33,12 +32,8 @@ solver = FDSBP(D_SBP, ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_gingerbread_man.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/2c6440b5f8a57db131061ad7aa78ee2b/raw/1f89fdf2c874ff678c78afb6fe8dc784bdfd421f/mesh_gingerbread_man.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/2c6440b5f8a57db131061ad7aa78ee2b/raw/1f89fdf2c874ff678c78afb6fe8dc784bdfd421f/mesh_gingerbread_man.mesh", + joinpath(@__DIR__, "mesh_gingerbread_man.mesh")) mesh = UnstructuredMesh2D(mesh_file) diff --git a/examples/unstructured_2d_fdsbp/elixir_euler_source_terms.jl b/examples/unstructured_2d_fdsbp/elixir_euler_source_terms.jl index edcd221bf59..5f11d41ad5c 100644 --- a/examples/unstructured_2d_fdsbp/elixir_euler_source_terms.jl +++ b/examples/unstructured_2d_fdsbp/elixir_euler_source_terms.jl @@ -1,5 +1,4 @@ -using Downloads: download using OrdinaryDiffEq using Trixi @@ -22,12 +21,8 @@ solver = FDSBP(D_SBP, ############################################################################### # Get the curved quad mesh from a file (downloads the file if not available locally) - -default_mesh_file = joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh") -isfile(default_mesh_file) || - download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", - default_mesh_file) -mesh_file = default_mesh_file +mesh_file = Trixi.download("https://gist.githubusercontent.com/andrewwinters5000/12ce661d7c354c3d94c74b964b0f1c96/raw/8275b9a60c6e7ebbdea5fc4b4f091c47af3d5273/mesh_periodic_square_with_twist.mesh", + joinpath(@__DIR__, "mesh_periodic_square_with_twist.mesh")) mesh = UnstructuredMesh2D(mesh_file, periodicity = true) diff --git a/src/Trixi.jl b/src/Trixi.jl index e18b2f6415c..8d74fbc9736 100644 --- a/src/Trixi.jl +++ b/src/Trixi.jl @@ -39,6 +39,7 @@ import SciMLBase: get_du, get_tmp_cache, u_modified!, get_proposed_dt, set_proposed_dt!, terminate!, remake, add_tstop!, has_tstop, first_tstop +using Downloads: Downloads using CodeTracking: CodeTracking using ConstructionBase: ConstructionBase using DiffEqCallbacks: PeriodicCallback, PeriodicCallbackAffect diff --git a/src/auxiliary/auxiliary.jl b/src/auxiliary/auxiliary.jl index 1f7d30d6aa8..92da9a5ba8b 100644 --- a/src/auxiliary/auxiliary.jl +++ b/src/auxiliary/auxiliary.jl @@ -345,4 +345,27 @@ function register_error_hints() return nothing end + +""" + Trixi.download(src_url, file_path) + +Download a file from given `src_url` to given `file_path` if +`file_path` is not already a file. This function just returns +`file_path`. +This is a small wrapper of `Downloads.download(src_url, file_path)` +that avoids race conditions when multiple MPI ranks are used. +""" +function download(src_url, file_path) + # Note that `mpi_isroot()` is also `true` if running + # in serial (without MPI). + if mpi_isroot() + isfile(file_path) || Downloads.download(src_url, file_path) + end + + if mpi_isparallel() + MPI.Barrier(mpi_comm()) + end + + return file_path +end end # @muladd diff --git a/src/auxiliary/t8code.jl b/src/auxiliary/t8code.jl index db01476bb86..7c1399fc803 100644 --- a/src/auxiliary/t8code.jl +++ b/src/auxiliary/t8code.jl @@ -46,230 +46,6 @@ function init_t8code() return nothing end -function trixi_t8_unref_forest(forest) - t8_forest_unref(Ref(forest)) -end - -function t8_free(ptr) - T8code.Libt8.sc_free(t8_get_package_id(), ptr) -end - -function trixi_t8_count_interfaces(forest) - # Check that forest is a committed, that is valid and usable, forest. - @assert t8_forest_is_committed(forest) != 0 - - # Get the number of local elements of forest. - num_local_elements = t8_forest_get_local_num_elements(forest) - # Get the number of ghost elements of forest. - num_ghost_elements = t8_forest_get_num_ghosts(forest) - # Get the number of trees that have elements of this process. - num_local_trees = t8_forest_get_num_local_trees(forest) - - current_index = t8_locidx_t(0) - - local_num_conform = 0 - local_num_mortars = 0 - local_num_boundary = 0 - - for itree in 0:(num_local_trees - 1) - tree_class = t8_forest_get_tree_class(forest, itree) - eclass_scheme = t8_forest_get_eclass_scheme(forest, tree_class) - - # Get the number of elements of this tree. - num_elements_in_tree = t8_forest_get_tree_num_elements(forest, itree) - - for ielement in 0:(num_elements_in_tree - 1) - element = t8_forest_get_element_in_tree(forest, itree, ielement) - - level = t8_element_level(eclass_scheme, element) - - num_faces = t8_element_num_faces(eclass_scheme, element) - - for iface in 0:(num_faces - 1) - pelement_indices_ref = Ref{Ptr{t8_locidx_t}}() - pneighbor_leafs_ref = Ref{Ptr{Ptr{t8_element}}}() - pneigh_scheme_ref = Ref{Ptr{t8_eclass_scheme}}() - - dual_faces_ref = Ref{Ptr{Cint}}() - num_neighbors_ref = Ref{Cint}() - - forest_is_balanced = Cint(1) - - t8_forest_leaf_face_neighbors(forest, itree, element, - pneighbor_leafs_ref, iface, dual_faces_ref, - num_neighbors_ref, - pelement_indices_ref, pneigh_scheme_ref, - forest_is_balanced) - - num_neighbors = num_neighbors_ref[] - neighbor_ielements = unsafe_wrap(Array, pelement_indices_ref[], - num_neighbors) - neighbor_leafs = unsafe_wrap(Array, pneighbor_leafs_ref[], num_neighbors) - neighbor_scheme = pneigh_scheme_ref[] - - if num_neighbors > 0 - neighbor_level = t8_element_level(neighbor_scheme, neighbor_leafs[1]) - - # Conforming interface: The second condition ensures we only visit the interface once. - if level == neighbor_level && current_index <= neighbor_ielements[1] - local_num_conform += 1 - elseif level < neighbor_level - local_num_mortars += 1 - end - else - local_num_boundary += 1 - end - - t8_free(dual_faces_ref[]) - t8_free(pneighbor_leafs_ref[]) - t8_free(pelement_indices_ref[]) - end # for - - current_index += 1 - end # for - end # for - - return (interfaces = local_num_conform, - mortars = local_num_mortars, - boundaries = local_num_boundary) -end - -function trixi_t8_fill_mesh_info(forest, elements, interfaces, mortars, boundaries, - boundary_names) - # Check that forest is a committed, that is valid and usable, forest. - @assert t8_forest_is_committed(forest) != 0 - - # Get the number of local elements of forest. - num_local_elements = t8_forest_get_local_num_elements(forest) - # Get the number of ghost elements of forest. - num_ghost_elements = t8_forest_get_num_ghosts(forest) - # Get the number of trees that have elements of this process. - num_local_trees = t8_forest_get_num_local_trees(forest) - - current_index = t8_locidx_t(0) - - local_num_conform = 0 - local_num_mortars = 0 - local_num_boundary = 0 - - for itree in 0:(num_local_trees - 1) - tree_class = t8_forest_get_tree_class(forest, itree) - eclass_scheme = t8_forest_get_eclass_scheme(forest, tree_class) - - # Get the number of elements of this tree. - num_elements_in_tree = t8_forest_get_tree_num_elements(forest, itree) - - for ielement in 0:(num_elements_in_tree - 1) - element = t8_forest_get_element_in_tree(forest, itree, ielement) - - level = t8_element_level(eclass_scheme, element) - - num_faces = t8_element_num_faces(eclass_scheme, element) - - for iface in 0:(num_faces - 1) - - # Compute the `orientation` of the touching faces. - if t8_element_is_root_boundary(eclass_scheme, element, iface) == 1 - cmesh = t8_forest_get_cmesh(forest) - itree_in_cmesh = t8_forest_ltreeid_to_cmesh_ltreeid(forest, itree) - iface_in_tree = t8_element_tree_face(eclass_scheme, element, iface) - orientation_ref = Ref{Cint}() - - t8_cmesh_get_face_neighbor(cmesh, itree_in_cmesh, iface_in_tree, C_NULL, - orientation_ref) - orientation = orientation_ref[] - else - orientation = zero(Cint) - end - - pelement_indices_ref = Ref{Ptr{t8_locidx_t}}() - pneighbor_leafs_ref = Ref{Ptr{Ptr{t8_element}}}() - pneigh_scheme_ref = Ref{Ptr{t8_eclass_scheme}}() - - dual_faces_ref = Ref{Ptr{Cint}}() - num_neighbors_ref = Ref{Cint}() - - forest_is_balanced = Cint(1) - - t8_forest_leaf_face_neighbors(forest, itree, element, - pneighbor_leafs_ref, iface, dual_faces_ref, - num_neighbors_ref, - pelement_indices_ref, pneigh_scheme_ref, - forest_is_balanced) - - num_neighbors = num_neighbors_ref[] - dual_faces = unsafe_wrap(Array, dual_faces_ref[], num_neighbors) - neighbor_ielements = unsafe_wrap(Array, pelement_indices_ref[], - num_neighbors) - neighbor_leafs = unsafe_wrap(Array, pneighbor_leafs_ref[], num_neighbors) - neighbor_scheme = pneigh_scheme_ref[] - - if num_neighbors > 0 - neighbor_level = t8_element_level(neighbor_scheme, neighbor_leafs[1]) - - # Conforming interface: The second condition ensures we only visit the interface once. - if level == neighbor_level && current_index <= neighbor_ielements[1] - local_num_conform += 1 - - faces = (iface, dual_faces[1]) - interface_id = local_num_conform - - # Write data to interfaces container. - interfaces.neighbor_ids[1, interface_id] = current_index + 1 - interfaces.neighbor_ids[2, interface_id] = neighbor_ielements[1] + 1 - - # Save interfaces.node_indices dimension specific in containers_3d.jl. - init_interface_node_indices!(interfaces, faces, orientation, - interface_id) - # Non-conforming interface. - elseif level < neighbor_level - local_num_mortars += 1 - - faces = (dual_faces[1], iface) - - mortar_id = local_num_mortars - - # Last entry is the large element. - mortars.neighbor_ids[end, mortar_id] = current_index + 1 - - # Fill in the `mortars.neighbor_ids` array and reorder if necessary. - init_mortar_neighbor_ids!(mortars, faces[2], faces[1], - orientation, neighbor_ielements, - mortar_id) - - # Fill in the `mortars.node_indices` array. - init_mortar_node_indices!(mortars, faces, orientation, mortar_id) - - # else: "level > neighbor_level" is skipped since we visit the mortar interface only once. - end - - # Domain boundary. - else - local_num_boundary += 1 - boundary_id = local_num_boundary - - boundaries.neighbor_ids[boundary_id] = current_index + 1 - - init_boundary_node_indices!(boundaries, iface, boundary_id) - - # One-based indexing. - boundaries.name[boundary_id] = boundary_names[iface + 1, itree + 1] - end - - t8_free(dual_faces_ref[]) - t8_free(pneighbor_leafs_ref[]) - t8_free(pelement_indices_ref[]) - end # for iface = ... - - current_index += 1 - end # for - end # for - - return (interfaces = local_num_conform, - mortars = local_num_mortars, - boundaries = local_num_boundary) -end - function trixi_t8_get_local_element_levels(forest) # Check that forest is a committed, that is valid and usable, forest. @assert t8_forest_is_committed(forest) != 0 @@ -341,23 +117,16 @@ function adapt_callback(forest, end function trixi_t8_adapt_new(old_forest, indicators) - # Check that forest is a committed, that is valid and usable, forest. - @assert t8_forest_is_committed(old_forest) != 0 - - # Init new forest. new_forest_ref = Ref{t8_forest_t}() t8_forest_init(new_forest_ref) new_forest = new_forest_ref[] - let set_from = C_NULL, recursive = 0, set_for_coarsening = 0, no_repartition = 0, - do_ghost = 1 - + let set_from = C_NULL, recursive = 0, no_repartition = 1, do_ghost = 1 t8_forest_set_user_data(new_forest, pointer(indicators)) t8_forest_set_adapt(new_forest, old_forest, @t8_adapt_callback(adapt_callback), recursive) t8_forest_set_balance(new_forest, set_from, no_repartition) - t8_forest_set_partition(new_forest, set_from, set_for_coarsening) - t8_forest_set_ghost(new_forest, do_ghost, T8_GHOST_FACES) # Note: MPI support not available yet so it is a dummy call. + t8_forest_set_ghost(new_forest, do_ghost, T8_GHOST_FACES) t8_forest_commit(new_forest) end diff --git a/src/callbacks_step/amr.jl b/src/callbacks_step/amr.jl index 5854c8617c3..6f57d6647fc 100644 --- a/src/callbacks_step/amr.jl +++ b/src/callbacks_step/amr.jl @@ -726,7 +726,7 @@ function (amr_callback::AMRCallback)(u_ode::AbstractVector, mesh::P4estMesh, return has_changed end -function (amr_callback::AMRCallback)(u_ode::AbstractVector, mesh::SerialT8codeMesh, +function (amr_callback::AMRCallback)(u_ode::AbstractVector, mesh::T8codeMesh, equations, dg::DG, cache, semi, t, iter; only_refine = false, only_coarsen = false, @@ -754,29 +754,29 @@ function (amr_callback::AMRCallback)(u_ode::AbstractVector, mesh::SerialT8codeMe @trixi_timeit timer() "adapt" begin difference = @trixi_timeit timer() "mesh" trixi_t8_adapt!(mesh, indicators) - @trixi_timeit timer() "solver" adapt!(u_ode, adaptor, mesh, equations, dg, - cache, difference) - end + # Store whether there were any cells coarsened or refined and perform load balancing. + has_changed = any(difference .!= 0) - # Store whether there were any cells coarsened or refined and perform load balancing. - has_changed = any(difference .!= 0) + # Check if mesh changed on other processes + if mpi_isparallel() + has_changed = MPI.Allreduce!(Ref(has_changed), |, mpi_comm())[] + end - # TODO: T8codeMesh for MPI not implemented yet. - # Check if mesh changed on other processes - # if mpi_isparallel() - # has_changed = MPI.Allreduce!(Ref(has_changed), |, mpi_comm())[] - # end + if has_changed + @trixi_timeit timer() "solver" adapt!(u_ode, adaptor, mesh, equations, dg, + cache, difference) + end + end if has_changed - # TODO: T8codeMesh for MPI not implemented yet. - # if mpi_isparallel() && amr_callback.dynamic_load_balancing - # @trixi_timeit timer() "dynamic load balancing" begin - # global_first_quadrant = unsafe_wrap(Array, mesh.p4est.global_first_quadrant, mpi_nranks() + 1) - # old_global_first_quadrant = copy(global_first_quadrant) - # partition!(mesh) - # rebalance_solver!(u_ode, mesh, equations, dg, cache, old_global_first_quadrant) - # end - # end + if mpi_isparallel() && amr_callback.dynamic_load_balancing + @trixi_timeit timer() "dynamic load balancing" begin + old_global_first_element_ids = get_global_first_element_ids(mesh) + partition!(mesh) + rebalance_solver!(u_ode, mesh, equations, dg, cache, + old_global_first_element_ids) + end + end reinitialize_boundaries!(semi.boundary_conditions, cache) end diff --git a/src/callbacks_step/amr_dg.jl b/src/callbacks_step/amr_dg.jl index 1dcfdccdea8..0a7055af409 100644 --- a/src/callbacks_step/amr_dg.jl +++ b/src/callbacks_step/amr_dg.jl @@ -6,11 +6,14 @@ #! format: noindent # Redistribute data for load balancing after partitioning the mesh -function rebalance_solver!(u_ode::AbstractVector, mesh::ParallelP4estMesh, equations, +function rebalance_solver!(u_ode::AbstractVector, + mesh::Union{ParallelP4estMesh, ParallelT8codeMesh}, + equations, dg::DGSEM, cache, old_global_first_quadrant) - # mpi ranks are 0-based, this array uses 1-based indices - global_first_quadrant = unsafe_wrap(Array, mesh.p4est.global_first_quadrant, - mpi_nranks() + 1) + + # MPI ranks are 0-based. This array uses 1-based indices. + global_first_quadrant = get_global_first_element_ids(mesh) + if global_first_quadrant[mpi_rank() + 1] == old_global_first_quadrant[mpi_rank() + 1] && global_first_quadrant[mpi_rank() + 2] == diff --git a/src/callbacks_step/amr_dg2d.jl b/src/callbacks_step/amr_dg2d.jl index b816bc06e65..94524b23a3a 100644 --- a/src/callbacks_step/amr_dg2d.jl +++ b/src/callbacks_step/amr_dg2d.jl @@ -385,7 +385,12 @@ function adapt!(u_ode::AbstractVector, adaptor, mesh::T8codeMesh{2}, equations, # Return early if there is nothing to do. if !any(difference .!= 0) - return nothing + if mpi_isparallel() + # MPICache init uses all-to-all communication -> reinitialize even if there is nothing to do + # locally (there still might be other MPI ranks that have refined elements) + reinitialize_containers!(mesh, equations, dg, cache) + end + return end # Number of (local) cells/elements. diff --git a/src/callbacks_step/amr_dg3d.jl b/src/callbacks_step/amr_dg3d.jl index 392cbba9e28..3f67951bafe 100644 --- a/src/callbacks_step/amr_dg3d.jl +++ b/src/callbacks_step/amr_dg3d.jl @@ -316,7 +316,12 @@ function adapt!(u_ode::AbstractVector, adaptor, mesh::T8codeMesh{3}, equations, # Return early if there is nothing to do. if !any(difference .!= 0) - return nothing + if mpi_isparallel() + # MPICache init uses all-to-all communication -> reinitialize even if there is nothing to do + # locally (there still might be other MPI ranks that have refined elements) + reinitialize_containers!(mesh, equations, dg, cache) + end + return end # Number of (local) cells/elements. diff --git a/src/callbacks_step/analysis_dg2d_parallel.jl b/src/callbacks_step/analysis_dg2d_parallel.jl index a04bf732604..000daa015dc 100644 --- a/src/callbacks_step/analysis_dg2d_parallel.jl +++ b/src/callbacks_step/analysis_dg2d_parallel.jl @@ -91,7 +91,8 @@ function calc_error_norms_per_element(func, u, t, analyzer, end function calc_error_norms(func, u, t, analyzer, - mesh::ParallelP4estMesh{2}, equations, + mesh::Union{ParallelP4estMesh{2}, ParallelT8codeMesh{2}}, + equations, initial_condition, dg::DGSEM, cache, cache_analysis) @unpack vandermonde, weights = analyzer @unpack node_coordinates, inverse_jacobian = cache.elements @@ -171,7 +172,8 @@ function integrate_via_indices(func::Func, u, end function integrate_via_indices(func::Func, u, - mesh::ParallelP4estMesh{2}, equations, + mesh::Union{ParallelP4estMesh{2}, ParallelT8codeMesh{2}}, + equations, dg::DGSEM, cache, args...; normalize = true) where {Func} @unpack weights = dg.basis diff --git a/src/callbacks_step/analysis_dg3d_parallel.jl b/src/callbacks_step/analysis_dg3d_parallel.jl index d8756d91c9d..de777be406d 100644 --- a/src/callbacks_step/analysis_dg3d_parallel.jl +++ b/src/callbacks_step/analysis_dg3d_parallel.jl @@ -6,7 +6,8 @@ #! format: noindent function calc_error_norms(func, u, t, analyzer, - mesh::ParallelP4estMesh{3}, equations, + mesh::Union{ParallelP4estMesh{3}, ParallelT8codeMesh{3}}, + equations, initial_condition, dg::DGSEM, cache, cache_analysis) @unpack vandermonde, weights = analyzer @unpack node_coordinates, inverse_jacobian = cache.elements @@ -64,7 +65,8 @@ function calc_error_norms(func, u, t, analyzer, end function integrate_via_indices(func::Func, u, - mesh::ParallelP4estMesh{3}, equations, + mesh::Union{ParallelP4estMesh{3}, ParallelT8codeMesh{3}}, + equations, dg::DGSEM, cache, args...; normalize = true) where {Func} @unpack weights = dg.basis diff --git a/src/callbacks_step/stepsize_dg2d.jl b/src/callbacks_step/stepsize_dg2d.jl index 673c3ba6aa6..c6d32c0f6dc 100644 --- a/src/callbacks_step/stepsize_dg2d.jl +++ b/src/callbacks_step/stepsize_dg2d.jl @@ -174,4 +174,36 @@ function max_dt(u, t, mesh::ParallelP4estMesh{2}, return dt end + +function max_dt(u, t, mesh::ParallelT8codeMesh{2}, + constant_speed::False, equations, dg::DG, cache) + # call the method accepting a general `mesh::T8codeMesh{2}` + # TODO: MPI, we should improve this; maybe we should dispatch on `u` + # and create some MPI array type, overloading broadcasting and mapreduce etc. + # Then, this specific array type should also work well with DiffEq etc. + dt = invoke(max_dt, + Tuple{typeof(u), typeof(t), T8codeMesh{2}, + typeof(constant_speed), typeof(equations), typeof(dg), + typeof(cache)}, + u, t, mesh, constant_speed, equations, dg, cache) + dt = MPI.Allreduce!(Ref(dt), min, mpi_comm())[] + + return dt +end + +function max_dt(u, t, mesh::ParallelT8codeMesh{2}, + constant_speed::True, equations, dg::DG, cache) + # call the method accepting a general `mesh::T8codeMesh{2}` + # TODO: MPI, we should improve this; maybe we should dispatch on `u` + # and create some MPI array type, overloading broadcasting and mapreduce etc. + # Then, this specific array type should also work well with DiffEq etc. + dt = invoke(max_dt, + Tuple{typeof(u), typeof(t), T8codeMesh{2}, + typeof(constant_speed), typeof(equations), typeof(dg), + typeof(cache)}, + u, t, mesh, constant_speed, equations, dg, cache) + dt = MPI.Allreduce!(Ref(dt), min, mpi_comm())[] + + return dt +end end # @muladd diff --git a/src/callbacks_step/stepsize_dg3d.jl b/src/callbacks_step/stepsize_dg3d.jl index 822ab2f87ec..664596f989e 100644 --- a/src/callbacks_step/stepsize_dg3d.jl +++ b/src/callbacks_step/stepsize_dg3d.jl @@ -150,4 +150,36 @@ function max_dt(u, t, mesh::ParallelP4estMesh{3}, return dt end + +function max_dt(u, t, mesh::ParallelT8codeMesh{3}, + constant_speed::False, equations, dg::DG, cache) + # call the method accepting a general `mesh::T8codeMesh{3}` + # TODO: MPI, we should improve this; maybe we should dispatch on `u` + # and create some MPI array type, overloading broadcasting and mapreduce etc. + # Then, this specific array type should also work well with DiffEq etc. + dt = invoke(max_dt, + Tuple{typeof(u), typeof(t), T8codeMesh{3}, + typeof(constant_speed), typeof(equations), typeof(dg), + typeof(cache)}, + u, t, mesh, constant_speed, equations, dg, cache) + dt = MPI.Allreduce!(Ref(dt), min, mpi_comm())[] + + return dt +end + +function max_dt(u, t, mesh::ParallelT8codeMesh{3}, + constant_speed::True, equations, dg::DG, cache) + # call the method accepting a general `mesh::T8codeMesh{3}` + # TODO: MPI, we should improve this; maybe we should dispatch on `u` + # and create some MPI array type, overloading broadcasting and mapreduce etc. + # Then, this specific array type should also work well with DiffEq etc. + dt = invoke(max_dt, + Tuple{typeof(u), typeof(t), T8codeMesh{3}, + typeof(constant_speed), typeof(equations), typeof(dg), + typeof(cache)}, + u, t, mesh, constant_speed, equations, dg, cache) + dt = MPI.Allreduce!(Ref(dt), min, mpi_comm())[] + + return dt +end end # @muladd diff --git a/src/meshes/p4est_mesh.jl b/src/meshes/p4est_mesh.jl index c5d39ef00c0..abe9d9345b5 100644 --- a/src/meshes/p4est_mesh.jl +++ b/src/meshes/p4est_mesh.jl @@ -1700,6 +1700,10 @@ function bilinear_interpolation!(coordinate, face_vertices, u, v) end end +function get_global_first_element_ids(mesh::P4estMesh) + return unsafe_wrap(Array, mesh.p4est.global_first_quadrant, mpi_nranks() + 1) +end + function balance!(mesh::P4estMesh{2}, init_fn = C_NULL) p4est_balance(mesh.p4est, P4EST_CONNECT_FACE, init_fn) # Due to a bug in `p4est`, the forest needs to be rebalanced twice sometimes diff --git a/src/meshes/t8code_mesh.jl b/src/meshes/t8code_mesh.jl index c9665a22af9..6fb4d861d10 100644 --- a/src/meshes/t8code_mesh.jl +++ b/src/meshes/t8code_mesh.jl @@ -7,8 +7,6 @@ to manage trees and mesh refinement. """ mutable struct T8codeMesh{NDIMS, RealT <: Real, IsParallel, NDIMSP2, NNODES} <: AbstractMesh{NDIMS} - cmesh :: Ptr{t8_cmesh} # cpointer to coarse mesh - scheme :: Ptr{t8_eclass_scheme} # cpointer to element scheme forest :: Ptr{t8_forest} # cpointer to forest is_parallel :: IsParallel @@ -25,14 +23,15 @@ mutable struct T8codeMesh{NDIMS, RealT <: Real, IsParallel, NDIMSP2, NNODES} <: nmortars :: Int nboundaries :: Int - function T8codeMesh{NDIMS}(cmesh, scheme, forest, tree_node_coordinates, nodes, + nmpiinterfaces :: Int + nmpimortars :: Int + + function T8codeMesh{NDIMS}(forest, tree_node_coordinates, nodes, boundary_names, current_filename) where {NDIMS} - is_parallel = False() + is_parallel = mpi_isparallel() ? True() : False() - mesh = new{NDIMS, Float64, typeof(is_parallel), NDIMS + 2, length(nodes)}(cmesh, - scheme, - forest, + mesh = new{NDIMS, Float64, typeof(is_parallel), NDIMS + 2, length(nodes)}(forest, is_parallel) mesh.nodes = nodes @@ -52,7 +51,7 @@ mutable struct T8codeMesh{NDIMS, RealT <: Real, IsParallel, NDIMSP2, NNODES} <: # further down. However, this might cause a pile-up of `mesh` # objects during long-running sessions. if !MPI.Finalized() - trixi_t8_unref_forest(mesh.forest) + t8_forest_unref(Ref(mesh.forest)) end end @@ -63,7 +62,7 @@ mutable struct T8codeMesh{NDIMS, RealT <: Real, IsParallel, NDIMSP2, NNODES} <: # more information. if haskey(ENV, "TRIXI_T8CODE_SC_FINALIZE") MPI.add_finalize_hook!() do - trixi_t8_unref_forest(mesh.forest) + t8_forest_unref(Ref(mesh.forest)) end end @@ -72,16 +71,15 @@ mutable struct T8codeMesh{NDIMS, RealT <: Real, IsParallel, NDIMSP2, NNODES} <: end const SerialT8codeMesh{NDIMS} = T8codeMesh{NDIMS, <:Real, <:False} +const ParallelT8codeMesh{NDIMS} = T8codeMesh{NDIMS, <:Real, <:True} @inline mpi_parallel(mesh::SerialT8codeMesh) = False() +@inline mpi_parallel(mesh::ParallelT8codeMesh) = True() @inline Base.ndims(::T8codeMesh{NDIMS}) where {NDIMS} = NDIMS @inline Base.real(::T8codeMesh{NDIMS, RealT}) where {NDIMS, RealT} = RealT -@inline ntrees(mesh::T8codeMesh) = Int(t8_forest_get_num_local_trees(mesh.forest)) +@inline ntrees(mesh::T8codeMesh) = size(mesh.tree_node_coordinates)[end] @inline ncells(mesh::T8codeMesh) = Int(t8_forest_get_local_num_elements(mesh.forest)) -@inline ninterfaces(mesh::T8codeMesh) = mesh.ninterfaces -@inline nmortars(mesh::T8codeMesh) = mesh.nmortars -@inline nboundaries(mesh::T8codeMesh) = mesh.nboundaries function Base.show(io::IO, mesh::T8codeMesh) print(io, "T8codeMesh{", ndims(mesh), ", ", real(mesh), "}") @@ -184,21 +182,23 @@ function T8codeMesh(trees_per_dimension; polydeg = 1, T8code.Libt8.p8est_connectivity_destroy(conn) end + do_face_ghost = mpi_isparallel() scheme = t8_scheme_new_default_cxx() - forest = t8_forest_new_uniform(cmesh, scheme, initial_refinement_level, 0, mpi_comm()) + forest = t8_forest_new_uniform(cmesh, scheme, initial_refinement_level, do_face_ghost, + mpi_comm()) basis = LobattoLegendreBasis(RealT, polydeg) nodes = basis.nodes + num_trees = t8_cmesh_get_num_trees(cmesh) + tree_node_coordinates = Array{RealT, NDIMS + 2}(undef, NDIMS, ntuple(_ -> length(nodes), NDIMS)..., - prod(trees_per_dimension)) + num_trees) # Get cell length in reference mesh: Omega_ref = [-1,1]^NDIMS. dx = [2 / n for n in trees_per_dimension] - num_local_trees = t8_cmesh_get_num_local_trees(cmesh) - # Non-periodic boundaries. boundary_names = fill(Symbol("---"), 2 * NDIMS, prod(trees_per_dimension)) @@ -208,7 +208,7 @@ function T8codeMesh(trees_per_dimension; polydeg = 1, mapping_ = mapping end - for itree in 1:num_local_trees + for itree in 1:num_trees veptr = t8_cmesh_get_tree_vertices(cmesh, itree - 1) verts = unsafe_wrap(Array, veptr, (3, 1 << NDIMS)) @@ -256,7 +256,7 @@ function T8codeMesh(trees_per_dimension; polydeg = 1, end end - return T8codeMesh{NDIMS}(cmesh, scheme, forest, tree_node_coordinates, nodes, + return T8codeMesh{NDIMS}(forest, tree_node_coordinates, nodes, boundary_names, "") end @@ -290,21 +290,25 @@ function T8codeMesh(cmesh::Ptr{t8_cmesh}; @assert (NDIMS == 2||NDIMS == 3) "NDIMS should be 2 or 3." + do_face_ghost = mpi_isparallel() scheme = t8_scheme_new_default_cxx() - forest = t8_forest_new_uniform(cmesh, scheme, initial_refinement_level, 0, mpi_comm()) + forest = t8_forest_new_uniform(cmesh, scheme, initial_refinement_level, do_face_ghost, + mpi_comm()) basis = LobattoLegendreBasis(RealT, polydeg) nodes = basis.nodes - num_local_trees = t8_cmesh_get_num_local_trees(cmesh) + num_trees = t8_cmesh_get_num_trees(cmesh) tree_node_coordinates = Array{RealT, NDIMS + 2}(undef, NDIMS, ntuple(_ -> length(nodes), NDIMS)..., - num_local_trees) + num_trees) nodes_in = [-1.0, 1.0] matrix = polynomial_interpolation_matrix(nodes_in, nodes) + num_local_trees = t8_cmesh_get_num_local_trees(cmesh) + if NDIMS == 2 data_in = Array{RealT, 3}(undef, 2, 2, 2) tmp1 = zeros(RealT, 2, length(nodes), length(nodes_in)) @@ -353,7 +357,7 @@ function T8codeMesh(cmesh::Ptr{t8_cmesh}; tmp1 = zeros(RealT, 3, length(nodes), length(nodes_in), length(nodes_in)) verts = zeros(3, 8) - for itree in 0:(num_local_trees - 1) + for itree in 0:(num_trees - 1) veptr = t8_cmesh_get_tree_vertices(cmesh, itree) # Note, `verts = unsafe_wrap(Array, veptr, (3, 1 << NDIMS))` @@ -387,9 +391,9 @@ function T8codeMesh(cmesh::Ptr{t8_cmesh}; map_node_coordinates!(tree_node_coordinates, mapping) # There's no simple and generic way to distinguish boundaries. Name all of them :all. - boundary_names = fill(:all, 2 * NDIMS, num_local_trees) + boundary_names = fill(:all, 2 * NDIMS, num_trees) - return T8codeMesh{NDIMS}(cmesh, scheme, forest, tree_node_coordinates, nodes, + return T8codeMesh{NDIMS}(forest, tree_node_coordinates, nodes, boundary_names, "") end @@ -442,7 +446,7 @@ function T8codeMesh(conn::Ptr{p8est_connectivity}; kwargs...) end """ - T8codeMesh{NDIMS}(meshfile::String; kwargs...) + T8codeMesh(meshfile::String, ndims; kwargs...) Main mesh constructor for the `T8codeMesh` that imports an unstructured, conforming mesh from a Gmsh mesh file (`.msh`). @@ -461,7 +465,6 @@ mesh from a Gmsh mesh file (`.msh`). - `initial_refinement_level::Integer`: refine the mesh uniformly to this level before the simulation starts. """ function T8codeMesh(meshfile::String, ndims; kwargs...) - # Prevent `t8code` from crashing Julia if the file doesn't exist. @assert isfile(meshfile) @@ -586,13 +589,525 @@ function adapt!(mesh::T8codeMesh, adapt_callback; recursive = true, balance = tr return nothing end -# TODO: Just a placeholder. Will be implemented later when MPI is supported. -function balance!(mesh::T8codeMesh, init_fn = C_NULL) +""" + Trixi.balance!(mesh::T8codeMesh) + +Balance a `T8codeMesh` to ensure 2^(NDIMS-1):1 face neighbors. +""" +function balance!(mesh::T8codeMesh) + new_forest_ref = Ref{t8_forest_t}() + t8_forest_init(new_forest_ref) + new_forest = new_forest_ref[] + + let set_from = mesh.forest, no_repartition = 1, do_ghost = 1 + t8_forest_set_balance(new_forest, set_from, no_repartition) + t8_forest_set_ghost(new_forest, do_ghost, T8_GHOST_FACES) + t8_forest_commit(new_forest) + end + + mesh.forest = new_forest + + return nothing +end + +""" + Trixi.partition!(mesh::T8codeMesh) + +Partition a `T8codeMesh` in order to redistribute elements evenly among MPI ranks. + +# Arguments +- `mesh::T8codeMesh`: Initialized mesh object. +""" +function partition!(mesh::T8codeMesh) + new_forest_ref = Ref{t8_forest_t}() + t8_forest_init(new_forest_ref) + new_forest = new_forest_ref[] + + let set_from = mesh.forest, do_ghost = 1, allow_for_coarsening = 1 + t8_forest_set_partition(new_forest, set_from, allow_for_coarsening) + t8_forest_set_ghost(new_forest, do_ghost, T8_GHOST_FACES) + t8_forest_commit(new_forest) + end + + mesh.forest = new_forest + return nothing end -# TODO: Just a placeholder. Will be implemented later when MPI is supported. -function partition!(mesh::T8codeMesh; allow_coarsening = true, weight_fn = C_NULL) +# Compute the global ids (zero-indexed) of first element in each MPI rank. +function get_global_first_element_ids(mesh::T8codeMesh) + n_elements_local = Int(t8_forest_get_local_num_elements(mesh.forest)) + n_elements_by_rank = Vector{Int}(undef, mpi_nranks()) + n_elements_by_rank[mpi_rank() + 1] = n_elements_local + MPI.Allgather!(MPI.UBuffer(n_elements_by_rank, 1), mpi_comm()) + return [sum(n_elements_by_rank[1:(rank - 1)]) for rank in 1:(mpi_nranks() + 1)] +end + +function count_interfaces(mesh::T8codeMesh) + @assert t8_forest_is_committed(mesh.forest) != 0 + + num_local_elements = t8_forest_get_local_num_elements(mesh.forest) + num_local_trees = t8_forest_get_num_local_trees(mesh.forest) + + current_index = t8_locidx_t(0) + + local_num_conform = 0 + local_num_mortars = 0 + local_num_boundary = 0 + + local_num_mpi_conform = 0 + local_num_mpi_mortars = 0 + + visited_global_mortar_ids = Set{UInt64}([]) + + max_level = t8_forest_get_maxlevel(mesh.forest) #UInt64 + max_tree_num_elements = UInt64(2^ndims(mesh))^max_level + + if mpi_isparallel() + ghost_num_trees = t8_forest_ghost_num_trees(mesh.forest) + + ghost_tree_element_offsets = [num_local_elements + + t8_forest_ghost_get_tree_element_offset(mesh.forest, + itree) + for itree in 0:(ghost_num_trees - 1)] + ghost_global_treeids = [t8_forest_ghost_get_global_treeid(mesh.forest, itree) + for itree in 0:(ghost_num_trees - 1)] + end + + for itree in 0:(num_local_trees - 1) + tree_class = t8_forest_get_tree_class(mesh.forest, itree) + eclass_scheme = t8_forest_get_eclass_scheme(mesh.forest, tree_class) + + num_elements_in_tree = t8_forest_get_tree_num_elements(mesh.forest, itree) + + global_itree = t8_forest_global_tree_id(mesh.forest, itree) + + for ielement in 0:(num_elements_in_tree - 1) + element = t8_forest_get_element_in_tree(mesh.forest, itree, ielement) + + level = t8_element_level(eclass_scheme, element) + + num_faces = t8_element_num_faces(eclass_scheme, element) + + # Note: This works only for forests of one element class. + current_linear_id = global_itree * max_tree_num_elements + + t8_element_get_linear_id(eclass_scheme, element, max_level) + + for iface in 0:(num_faces - 1) + pelement_indices_ref = Ref{Ptr{t8_locidx_t}}() + pneighbor_leafs_ref = Ref{Ptr{Ptr{t8_element}}}() + pneigh_scheme_ref = Ref{Ptr{t8_eclass_scheme}}() + + dual_faces_ref = Ref{Ptr{Cint}}() + num_neighbors_ref = Ref{Cint}() + + forest_is_balanced = Cint(1) + + t8_forest_leaf_face_neighbors(mesh.forest, itree, element, + pneighbor_leafs_ref, iface, dual_faces_ref, + num_neighbors_ref, + pelement_indices_ref, pneigh_scheme_ref, + forest_is_balanced) + + num_neighbors = num_neighbors_ref[] + dual_faces = unsafe_wrap(Array, dual_faces_ref[], num_neighbors) + neighbor_ielements = unsafe_wrap(Array, pelement_indices_ref[], + num_neighbors) + neighbor_leafs = unsafe_wrap(Array, pneighbor_leafs_ref[], num_neighbors) + neighbor_scheme = pneigh_scheme_ref[] + + if num_neighbors == 0 + local_num_boundary += 1 + else + neighbor_level = t8_element_level(neighbor_scheme, neighbor_leafs[1]) + + if all(neighbor_ielements .< num_local_elements) + # Conforming interface: The second condition ensures we + # only visit the interface once. + if level == neighbor_level && current_index <= neighbor_ielements[1] + local_num_conform += 1 + elseif level < neighbor_level + local_num_mortars += 1 + # `else level > neighbor_level` is ignored since we + # only want to count the mortar interface once. + end + else + if level == neighbor_level + local_num_mpi_conform += 1 + elseif level < neighbor_level + local_num_mpi_mortars += 1 + + global_mortar_id = 2 * ndims(mesh) * current_linear_id + iface + + else # level > neighbor_level + neighbor_global_ghost_itree = ghost_global_treeids[findlast(ghost_tree_element_offsets .<= + neighbor_ielements[1])] + neighbor_linear_id = neighbor_global_ghost_itree * + max_tree_num_elements + + t8_element_get_linear_id(neighbor_scheme, + neighbor_leafs[1], + max_level) + global_mortar_id = 2 * ndims(mesh) * neighbor_linear_id + + dual_faces[1] + + if !(global_mortar_id in visited_global_mortar_ids) + push!(visited_global_mortar_ids, global_mortar_id) + local_num_mpi_mortars += 1 + end + end + end + end + + t8_free(dual_faces_ref[]) + t8_free(pneighbor_leafs_ref[]) + t8_free(pelement_indices_ref[]) + end # for + + current_index += 1 + end # for + end # for + + return (interfaces = local_num_conform, + mortars = local_num_mortars, + boundaries = local_num_boundary, + mpi_interfaces = local_num_mpi_conform, + mpi_mortars = local_num_mpi_mortars) +end + +# I know this routine is an unmaintainable behemoth. However, I see no real +# and elegant way to refactor this into, for example, smaller parts. The +# `t8_forest_leaf_face_neighbors` routine is as of now rather costly and it +# makes sense to query it only once per face per element and extract all the +# information needed at once in order to fill the connectivity information. +# Instead, I opted for good documentation. +function fill_mesh_info!(mesh::T8codeMesh, interfaces, mortars, boundaries, + boundary_names; mpi_mesh_info = nothing) + @assert t8_forest_is_committed(mesh.forest) != 0 + + num_local_elements = t8_forest_get_local_num_elements(mesh.forest) + num_local_trees = t8_forest_get_num_local_trees(mesh.forest) + + if !isnothing(mpi_mesh_info) + #! format: off + remotes = t8_forest_ghost_get_remotes(mesh.forest) + ghost_num_trees = t8_forest_ghost_num_trees(mesh.forest) + + ghost_remote_first_elem = [num_local_elements + + t8_forest_ghost_remote_first_elem(mesh.forest, remote) + for remote in remotes] + + ghost_tree_element_offsets = [num_local_elements + + t8_forest_ghost_get_tree_element_offset(mesh.forest, itree) + for itree in 0:(ghost_num_trees - 1)] + + ghost_global_treeids = [t8_forest_ghost_get_global_treeid(mesh.forest, itree) + for itree in 0:(ghost_num_trees - 1)] + #! format: on + end + + # Process-local index of the current element in the space-filling curve. + current_index = t8_locidx_t(0) + + # Increment counters for the different interface/mortar/boundary types. + local_num_conform = 0 + local_num_mortars = 0 + local_num_boundary = 0 + + local_num_mpi_conform = 0 + local_num_mpi_mortars = 0 + + # Works for quads and hexs only. This mapping is needed in the MPI mortar + # sections below. + map_iface_to_ichild_to_position = [ + # 0 1 2 3 4 5 6 7 ichild/iface + [1, 0, 2, 0, 3, 0, 4, 0], # 0 + [0, 1, 0, 2, 0, 3, 0, 4], # 1 + [1, 2, 0, 0, 3, 4, 0, 0], # 2 + [0, 0, 1, 2, 0, 0, 3, 4], # 3 + [1, 2, 3, 4, 0, 0, 0, 0], # 4 + [0, 0, 0, 0, 1, 2, 3, 4], # 5 + ] + + # Helper variables to compute unique global MPI interface/mortar ids. + max_level = t8_forest_get_maxlevel(mesh.forest) #UInt64 + max_tree_num_elements = UInt64(2^ndims(mesh))^max_level + + # These two variables help to ensure that we count MPI mortars from smaller + # elements point of view only once. + visited_global_mortar_ids = Set{UInt64}([]) + global_mortar_id_to_local = Dict{UInt64, Int}([]) + + # Loop over all local trees. + for itree in 0:(num_local_trees - 1) + tree_class = t8_forest_get_tree_class(mesh.forest, itree) + eclass_scheme = t8_forest_get_eclass_scheme(mesh.forest, tree_class) + + num_elements_in_tree = t8_forest_get_tree_num_elements(mesh.forest, itree) + + global_itree = t8_forest_global_tree_id(mesh.forest, itree) + + # Loop over all local elements of the current local tree. + for ielement in 0:(num_elements_in_tree - 1) + element = t8_forest_get_element_in_tree(mesh.forest, itree, ielement) + + level = t8_element_level(eclass_scheme, element) + + num_faces = t8_element_num_faces(eclass_scheme, element) + + # Note: This works only for forests of one element class. + current_linear_id = global_itree * max_tree_num_elements + + t8_element_get_linear_id(eclass_scheme, element, max_level) + + # Loop over all faces of the current local element. + for iface in 0:(num_faces - 1) + # Compute the `orientation` of the touching faces. + if t8_element_is_root_boundary(eclass_scheme, element, iface) == 1 + cmesh = t8_forest_get_cmesh(mesh.forest) + itree_in_cmesh = t8_forest_ltreeid_to_cmesh_ltreeid(mesh.forest, itree) + iface_in_tree = t8_element_tree_face(eclass_scheme, element, iface) + orientation_ref = Ref{Cint}() + + t8_cmesh_get_face_neighbor(cmesh, itree_in_cmesh, iface_in_tree, C_NULL, + orientation_ref) + orientation = orientation_ref[] + else + orientation = zero(Cint) + end + + pelement_indices_ref = Ref{Ptr{t8_locidx_t}}() + pneighbor_leafs_ref = Ref{Ptr{Ptr{t8_element}}}() + pneigh_scheme_ref = Ref{Ptr{t8_eclass_scheme}}() + + dual_faces_ref = Ref{Ptr{Cint}}() + num_neighbors_ref = Ref{Cint}() + + forest_is_balanced = Cint(1) + + # Query neighbor information from t8code. + t8_forest_leaf_face_neighbors(mesh.forest, itree, element, + pneighbor_leafs_ref, iface, dual_faces_ref, + num_neighbors_ref, + pelement_indices_ref, pneigh_scheme_ref, + forest_is_balanced) + + num_neighbors = num_neighbors_ref[] + dual_faces = unsafe_wrap(Array, dual_faces_ref[], num_neighbors) + neighbor_ielements = unsafe_wrap(Array, pelement_indices_ref[], + num_neighbors) + neighbor_leafs = unsafe_wrap(Array, pneighbor_leafs_ref[], num_neighbors) + neighbor_scheme = pneigh_scheme_ref[] + + # Now we check for the different cases. The nested if-structure is as follows: + # + # if `boundary`: + # + # + # else: // It must be an interface or mortar. + # + # if `all neighbors are local elements`: + # + # if `local interface`: + # + # elseif `local mortar from larger element point of view`: + # + # else: // `local mortar from smaller elements point of view` + # // We only count local mortars once. + # + # else: // It must be either a MPI interface or a MPI mortar. + # + # if `MPI interface`: + # + # elseif `MPI mortar from larger element point of view`: + # + # else: // `MPI mortar from smaller elements point of view` + # + # + # // end + + # Domain boundary. + if num_neighbors == 0 + local_num_boundary += 1 + boundary_id = local_num_boundary + + boundaries.neighbor_ids[boundary_id] = current_index + 1 + + init_boundary_node_indices!(boundaries, iface, boundary_id) + + # One-based indexing. + boundaries.name[boundary_id] = boundary_names[iface + 1, itree + 1] + + # Interface or mortar. + else + neighbor_level = t8_element_level(neighbor_scheme, neighbor_leafs[1]) + + # Local interface or mortar. + if all(neighbor_ielements .< num_local_elements) + + # Local interface: The second condition ensures we only visit the interface once. + if level == neighbor_level && current_index <= neighbor_ielements[1] + local_num_conform += 1 + + interfaces.neighbor_ids[1, local_num_conform] = current_index + + 1 + interfaces.neighbor_ids[2, local_num_conform] = neighbor_ielements[1] + + 1 + + init_interface_node_indices!(interfaces, (iface, dual_faces[1]), + orientation, + local_num_conform) + # Local mortar. + elseif level < neighbor_level + local_num_mortars += 1 + + # Last entry is the large element. + mortars.neighbor_ids[end, local_num_mortars] = current_index + 1 + + init_mortar_neighbor_ids!(mortars, iface, dual_faces[1], + orientation, neighbor_ielements, + local_num_mortars) + + init_mortar_node_indices!(mortars, (dual_faces[1], iface), + orientation, local_num_mortars) + + # else: `level > neighbor_level` is skipped since we visit the mortar interface only once. + end + + # MPI interface or MPI mortar. + else + + # MPI interface. + if level == neighbor_level + local_num_mpi_conform += 1 + + neighbor_global_ghost_itree = ghost_global_treeids[findlast(ghost_tree_element_offsets .<= + neighbor_ielements[1])] + + neighbor_linear_id = neighbor_global_ghost_itree * + max_tree_num_elements + + t8_element_get_linear_id(neighbor_scheme, + neighbor_leafs[1], + max_level) + + if current_linear_id < neighbor_linear_id + local_side = 1 + smaller_iface = iface + smaller_linear_id = current_linear_id + faces = (iface, dual_faces[1]) + else + local_side = 2 + smaller_iface = dual_faces[1] + smaller_linear_id = neighbor_linear_id + faces = (dual_faces[1], iface) + end + + global_interface_id = 2 * ndims(mesh) * smaller_linear_id + + smaller_iface + + mpi_mesh_info.mpi_interfaces.local_neighbor_ids[local_num_mpi_conform] = current_index + + 1 + mpi_mesh_info.mpi_interfaces.local_sides[local_num_mpi_conform] = local_side + + init_mpi_interface_node_indices!(mpi_mesh_info.mpi_interfaces, + faces, local_side, orientation, + local_num_mpi_conform) + + neighbor_rank = remotes[findlast(ghost_remote_first_elem .<= + neighbor_ielements[1])] + mpi_mesh_info.neighbor_ranks_interface[local_num_mpi_conform] = neighbor_rank + + mpi_mesh_info.global_interface_ids[local_num_mpi_conform] = global_interface_id + + # MPI Mortar: from larger element point of view + elseif level < neighbor_level + local_num_mpi_mortars += 1 + + global_mortar_id = 2 * ndims(mesh) * current_linear_id + iface + + neighbor_ids = neighbor_ielements .+ 1 + + local_neighbor_positions = findall(neighbor_ids .<= + num_local_elements) + local_neighbor_ids = [neighbor_ids[i] + for i in local_neighbor_positions] + local_neighbor_positions = [map_iface_to_ichild_to_position[dual_faces[1] + 1][t8_element_child_id(neighbor_scheme, neighbor_leafs[i]) + 1] + for i in local_neighbor_positions] + + # Last entry is the large element. + push!(local_neighbor_ids, current_index + 1) + push!(local_neighbor_positions, 2^(ndims(mesh) - 1) + 1) + + mpi_mesh_info.mpi_mortars.local_neighbor_ids[local_num_mpi_mortars] = local_neighbor_ids + mpi_mesh_info.mpi_mortars.local_neighbor_positions[local_num_mpi_mortars] = local_neighbor_positions + + init_mortar_node_indices!(mpi_mesh_info.mpi_mortars, + (dual_faces[1], iface), orientation, + local_num_mpi_mortars) + + neighbor_ranks = [remotes[findlast(ghost_remote_first_elem .<= + ineighbor_ghost)] + for ineighbor_ghost in filter(x -> x >= + num_local_elements, + neighbor_ielements)] + mpi_mesh_info.neighbor_ranks_mortar[local_num_mpi_mortars] = neighbor_ranks + + mpi_mesh_info.global_mortar_ids[local_num_mpi_mortars] = global_mortar_id + + # MPI Mortar: from smaller elements point of view + else + neighbor_global_ghost_itree = ghost_global_treeids[findlast(ghost_tree_element_offsets .<= + neighbor_ielements[1])] + neighbor_linear_id = neighbor_global_ghost_itree * + max_tree_num_elements + + t8_element_get_linear_id(neighbor_scheme, + neighbor_leafs[1], + max_level) + global_mortar_id = 2 * ndims(mesh) * neighbor_linear_id + + dual_faces[1] + + if global_mortar_id in visited_global_mortar_ids + local_mpi_mortar_id = global_mortar_id_to_local[global_mortar_id] + + push!(mpi_mesh_info.mpi_mortars.local_neighbor_ids[local_mpi_mortar_id], + current_index + 1) + push!(mpi_mesh_info.mpi_mortars.local_neighbor_positions[local_mpi_mortar_id], + map_iface_to_ichild_to_position[iface + 1][t8_element_child_id(eclass_scheme, element) + 1]) + else + local_num_mpi_mortars += 1 + local_mpi_mortar_id = local_num_mpi_mortars + push!(visited_global_mortar_ids, global_mortar_id) + global_mortar_id_to_local[global_mortar_id] = local_mpi_mortar_id + + mpi_mesh_info.mpi_mortars.local_neighbor_ids[local_mpi_mortar_id] = [ + current_index + 1, + ] + mpi_mesh_info.mpi_mortars.local_neighbor_positions[local_mpi_mortar_id] = [ + map_iface_to_ichild_to_position[iface + 1][t8_element_child_id(eclass_scheme, element) + 1], + ] + init_mortar_node_indices!(mpi_mesh_info.mpi_mortars, + (iface, dual_faces[1]), + orientation, local_mpi_mortar_id) + + neighbor_ranks = [ + remotes[findlast(ghost_remote_first_elem .<= + neighbor_ielements[1])], + ] + mpi_mesh_info.neighbor_ranks_mortar[local_mpi_mortar_id] = neighbor_ranks + + mpi_mesh_info.global_mortar_ids[local_mpi_mortar_id] = global_mortar_id + end + end + end + end + + t8_free(dual_faces_ref[]) + t8_free(pneighbor_leafs_ref[]) + t8_free(pelement_indices_ref[]) + end # for iface + + current_index += 1 + end # for ielement + end # for itree + return nothing end diff --git a/src/solvers/dgsem_p4est/containers_parallel.jl b/src/solvers/dgsem_p4est/containers_parallel.jl index 7c7bd868457..fd2749155bb 100644 --- a/src/solvers/dgsem_p4est/containers_parallel.jl +++ b/src/solvers/dgsem_p4est/containers_parallel.jl @@ -43,7 +43,8 @@ function Base.resize!(mpi_interfaces::P4estMPIInterfaceContainer, capacity) end # Create MPI interface container and initialize interface data -function init_mpi_interfaces(mesh::ParallelP4estMesh, equations, basis, elements) +function init_mpi_interfaces(mesh::Union{ParallelP4estMesh, ParallelT8codeMesh}, + equations, basis, elements) NDIMS = ndims(elements) uEltype = eltype(elements) @@ -133,7 +134,8 @@ function Base.resize!(mpi_mortars::P4estMPIMortarContainer, capacity) end # Create MPI mortar container and initialize MPI mortar data -function init_mpi_mortars(mesh::ParallelP4estMesh, equations, basis, elements) +function init_mpi_mortars(mesh::Union{ParallelP4estMesh, ParallelT8codeMesh}, equations, + basis, elements) NDIMS = ndims(mesh) RealT = real(mesh) uEltype = eltype(elements) diff --git a/src/solvers/dgsem_p4est/dg_2d_parallel.jl b/src/solvers/dgsem_p4est/dg_2d_parallel.jl index a8887351c46..3bf0cd0cab5 100644 --- a/src/solvers/dgsem_p4est/dg_2d_parallel.jl +++ b/src/solvers/dgsem_p4est/dg_2d_parallel.jl @@ -6,7 +6,8 @@ #! format: noindent function prolong2mpiinterfaces!(cache, u, - mesh::ParallelP4estMesh{2}, + mesh::Union{ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, equations, surface_integral, dg::DG) @unpack mpi_interfaces = cache index_range = eachnode(dg) @@ -43,7 +44,8 @@ function prolong2mpiinterfaces!(cache, u, end function calc_mpi_interface_flux!(surface_flux_values, - mesh::ParallelP4estMesh{2}, + mesh::Union{ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, nonconservative_terms, equations, surface_integral, dg::DG, cache) @unpack local_neighbor_ids, node_indices, local_sides = cache.mpi_interfaces @@ -106,7 +108,8 @@ end # Inlined version of the interface flux computation for conservation laws @inline function calc_mpi_interface_flux!(surface_flux_values, - mesh::P4estMesh{2}, + mesh::Union{ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, nonconservative_terms::False, equations, surface_integral, dg::DG, cache, interface_index, normal_direction, @@ -131,7 +134,8 @@ end end function prolong2mpimortars!(cache, u, - mesh::ParallelP4estMesh{2}, equations, + mesh::Union{ParallelP4estMesh{2}, ParallelT8codeMesh{2}}, + equations, mortar_l2::LobattoLegendreMortarL2, surface_integral, dg::DGSEM) @unpack node_indices = cache.mpi_mortars @@ -199,7 +203,7 @@ function prolong2mpimortars!(cache, u, end function calc_mpi_mortar_flux!(surface_flux_values, - mesh::ParallelP4estMesh{2}, + mesh::Union{ParallelP4estMesh{2}, ParallelT8codeMesh{2}}, nonconservative_terms, equations, mortar_l2::LobattoLegendreMortarL2, surface_integral, dg::DG, cache) @@ -253,7 +257,8 @@ end # Inlined version of the mortar flux computation on small elements for conservation laws @inline function calc_mpi_mortar_flux!(fstar, - mesh::ParallelP4estMesh{2}, + mesh::Union{ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, nonconservative_terms::False, equations, surface_integral, dg::DG, cache, mortar_index, position_index, normal_direction, @@ -271,7 +276,9 @@ end end @inline function mpi_mortar_fluxes_to_elements!(surface_flux_values, - mesh::ParallelP4estMesh{2}, equations, + mesh::Union{ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, + equations, mortar_l2::LobattoLegendreMortarL2, dg::DGSEM, cache, mortar, fstar, u_buffer) diff --git a/src/solvers/dgsem_p4est/dg_3d_parallel.jl b/src/solvers/dgsem_p4est/dg_3d_parallel.jl index 13bf2a1a2eb..e504e06d2c4 100644 --- a/src/solvers/dgsem_p4est/dg_3d_parallel.jl +++ b/src/solvers/dgsem_p4est/dg_3d_parallel.jl @@ -6,7 +6,7 @@ #! format: noindent function rhs!(du, u, t, - mesh::ParallelP4estMesh{3}, equations, + mesh::Union{ParallelP4estMesh{3}, ParallelT8codeMesh{3}}, equations, initial_condition, boundary_conditions, source_terms::Source, dg::DG, cache) where {Source} # Start to receive MPI data @@ -113,7 +113,8 @@ function rhs!(du, u, t, end function prolong2mpiinterfaces!(cache, u, - mesh::ParallelP4estMesh{3}, + mesh::Union{ParallelP4estMesh{3}, + ParallelT8codeMesh{3}}, equations, surface_integral, dg::DG) @unpack mpi_interfaces = cache index_range = eachnode(dg) @@ -160,7 +161,8 @@ function prolong2mpiinterfaces!(cache, u, end function calc_mpi_interface_flux!(surface_flux_values, - mesh::ParallelP4estMesh{3}, + mesh::Union{ParallelP4estMesh{3}, + ParallelT8codeMesh{3}}, nonconservative_terms, equations, surface_integral, dg::DG, cache) @unpack local_neighbor_ids, node_indices, local_sides = cache.mpi_interfaces @@ -237,7 +239,8 @@ end # Inlined version of the interface flux computation for conservation laws @inline function calc_mpi_interface_flux!(surface_flux_values, - mesh::P4estMesh{3}, + mesh::Union{ParallelP4estMesh{3}, + ParallelT8codeMesh{3}}, nonconservative_terms::False, equations, surface_integral, dg::DG, cache, interface_index, normal_direction, @@ -265,7 +268,8 @@ end end function prolong2mpimortars!(cache, u, - mesh::ParallelP4estMesh{3}, equations, + mesh::Union{ParallelP4estMesh{3}, ParallelT8codeMesh{3}}, + equations, mortar_l2::LobattoLegendreMortarL2, surface_integral, dg::DGSEM) @unpack node_indices = cache.mpi_mortars @@ -374,7 +378,7 @@ function prolong2mpimortars!(cache, u, end function calc_mpi_mortar_flux!(surface_flux_values, - mesh::ParallelP4estMesh{3}, + mesh::Union{ParallelP4estMesh{3}, ParallelT8codeMesh{3}}, nonconservative_terms, equations, mortar_l2::LobattoLegendreMortarL2, surface_integral, dg::DG, cache) @@ -437,7 +441,8 @@ end # Inlined version of the mortar flux computation on small elements for conservation laws @inline function calc_mpi_mortar_flux!(fstar, - mesh::ParallelP4estMesh{3}, + mesh::Union{ParallelP4estMesh{3}, + ParallelT8codeMesh{3}}, nonconservative_terms::False, equations, surface_integral, dg::DG, cache, mortar_index, position_index, normal_direction, @@ -456,7 +461,9 @@ end end @inline function mpi_mortar_fluxes_to_elements!(surface_flux_values, - mesh::ParallelP4estMesh{3}, equations, + mesh::Union{ParallelP4estMesh{3}, + ParallelT8codeMesh{3}}, + equations, mortar_l2::LobattoLegendreMortarL2, dg::DGSEM, cache, mortar, fstar, u_buffer, fstar_tmp) diff --git a/src/solvers/dgsem_p4est/dg_parallel.jl b/src/solvers/dgsem_p4est/dg_parallel.jl index 712ede2bfce..eaa6ab5cee2 100644 --- a/src/solvers/dgsem_p4est/dg_parallel.jl +++ b/src/solvers/dgsem_p4est/dg_parallel.jl @@ -166,7 +166,8 @@ end # at `index_base`+1 in the MPI buffer. `data_size` is the data size associated with each small # position (i.e. position 1 or 2). The data corresponding to the large side (i.e. position 3) has # size `2 * data_size`. -@inline function buffer_mortar_indices(mesh::ParallelP4estMesh{2}, index_base, +@inline function buffer_mortar_indices(mesh::Union{ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, index_base, data_size) return ( # first, last for local element in position 1 (small element) @@ -185,7 +186,8 @@ end # at `index_base`+1 in the MPI buffer. `data_size` is the data size associated with each small # position (i.e. position 1 to 4). The data corresponding to the large side (i.e. position 5) has # size `4 * data_size`. -@inline function buffer_mortar_indices(mesh::ParallelP4estMesh{3}, index_base, +@inline function buffer_mortar_indices(mesh::Union{ParallelP4estMesh{3}, + ParallelT8codeMesh{3}}, index_base, data_size) return ( # first, last for local element in position 1 (small element) @@ -491,7 +493,8 @@ end # Exchange normal directions of small elements of the MPI mortars. They are needed on all involved # MPI ranks to calculate the mortar fluxes. -function exchange_normal_directions!(mpi_mortars, mpi_cache, mesh::ParallelP4estMesh, +function exchange_normal_directions!(mpi_mortars, mpi_cache, + mesh::Union{ParallelP4estMesh, ParallelT8codeMesh}, n_nodes) RealT = real(mesh) n_dims = ndims(mesh) diff --git a/src/solvers/dgsem_t8code/containers.jl b/src/solvers/dgsem_t8code/containers.jl index 093feb2985a..d7ff79fbf2f 100644 --- a/src/solvers/dgsem_t8code/containers.jl +++ b/src/solvers/dgsem_t8code/containers.jl @@ -18,19 +18,22 @@ function reinitialize_containers!(mesh::T8codeMesh, equations, dg::DGSEM, cache) @unpack boundaries = cache resize!(boundaries, mesh.nboundaries) - trixi_t8_fill_mesh_info(mesh.forest, elements, interfaces, mortars, boundaries, - mesh.boundary_names) + fill_mesh_info!(mesh, interfaces, mortars, boundaries, + mesh.boundary_names) return nothing end function count_required_surfaces!(mesh::T8codeMesh) - counts = trixi_t8_count_interfaces(mesh.forest) + counts = count_interfaces(mesh) mesh.nmortars = counts.mortars mesh.ninterfaces = counts.interfaces mesh.nboundaries = counts.boundaries + mesh.nmpimortars = counts.mpi_mortars + mesh.nmpiinterfaces = counts.mpi_interfaces + return counts end @@ -38,7 +41,9 @@ end function count_required_surfaces(mesh::T8codeMesh) return (interfaces = mesh.ninterfaces, mortars = mesh.nmortars, - boundaries = mesh.nboundaries) + boundaries = mesh.nboundaries, + mpi_interfaces = mesh.nmpiinterfaces, + mpi_mortars = mesh.nmpimortars) end # Compatibility to `dgsem_p4est/containers.jl`. diff --git a/src/solvers/dgsem_t8code/containers_2d.jl b/src/solvers/dgsem_t8code/containers_2d.jl index bf77826a34b..ce525bfdf65 100644 --- a/src/solvers/dgsem_t8code/containers_2d.jl +++ b/src/solvers/dgsem_t8code/containers_2d.jl @@ -26,6 +26,7 @@ function calc_node_coordinates!(node_coordinates, tree_class = t8_forest_get_tree_class(mesh.forest, itree) eclass_scheme = t8_forest_get_eclass_scheme(mesh.forest, tree_class) num_elements_in_tree = t8_forest_get_tree_num_elements(mesh.forest, itree) + global_itree = t8_forest_global_tree_id(mesh.forest, itree) for ielement in 0:(num_elements_in_tree - 1) element = t8_forest_get_element_in_tree(mesh.forest, itree, ielement) @@ -55,7 +56,7 @@ function calc_node_coordinates!(node_coordinates, multiply_dimensionwise!(view(node_coordinates, :, :, :, current_index += 1), matrix1, matrix2, view(mesh.tree_node_coordinates, :, :, :, - itree + 1), + global_itree + 1), tmp1) end end diff --git a/src/solvers/dgsem_t8code/containers_3d.jl b/src/solvers/dgsem_t8code/containers_3d.jl index f2d54ff07da..4d56bc734aa 100644 --- a/src/solvers/dgsem_t8code/containers_3d.jl +++ b/src/solvers/dgsem_t8code/containers_3d.jl @@ -28,6 +28,7 @@ function calc_node_coordinates!(node_coordinates, tree_class = t8_forest_get_tree_class(mesh.forest, itree) eclass_scheme = t8_forest_get_eclass_scheme(mesh.forest, tree_class) num_elements_in_tree = t8_forest_get_tree_num_elements(mesh.forest, itree) + global_itree = t8_forest_global_tree_id(mesh.forest, itree) for ielement in 0:(num_elements_in_tree - 1) element = t8_forest_get_element_in_tree(mesh.forest, itree, ielement) @@ -63,7 +64,7 @@ function calc_node_coordinates!(node_coordinates, current_index += 1), matrix1, matrix2, matrix3, view(mesh.tree_node_coordinates, :, :, :, :, - itree + 1), + global_itree + 1), tmp1) end end diff --git a/src/solvers/dgsem_t8code/containers_parallel.jl b/src/solvers/dgsem_t8code/containers_parallel.jl new file mode 100644 index 00000000000..0cb3f5887a0 --- /dev/null +++ b/src/solvers/dgsem_t8code/containers_parallel.jl @@ -0,0 +1,65 @@ +function reinitialize_containers!(mesh::ParallelT8codeMesh, equations, dg::DGSEM, cache) + @unpack elements, interfaces, boundaries, mortars, mpi_interfaces, mpi_mortars, + mpi_cache = cache + resize!(elements, ncells(mesh)) + init_elements!(elements, mesh, dg.basis) + + count_required_surfaces!(mesh) + required = count_required_surfaces(mesh) + + resize!(interfaces, required.interfaces) + + resize!(boundaries, required.boundaries) + + resize!(mortars, required.mortars) + + resize!(mpi_interfaces, required.mpi_interfaces) + + resize!(mpi_mortars, required.mpi_mortars) + + mpi_mesh_info = (mpi_mortars = mpi_mortars, + mpi_interfaces = mpi_interfaces, + + # Temporary arrays for updating `mpi_cache`. + global_mortar_ids = fill(UInt64(0), nmpimortars(mpi_mortars)), + global_interface_ids = fill(UInt64(0), nmpiinterfaces(mpi_interfaces)), + neighbor_ranks_mortar = Vector{Vector{Int}}(undef, + nmpimortars(mpi_mortars)), + neighbor_ranks_interface = fill(-1, nmpiinterfaces(mpi_interfaces))) + + fill_mesh_info!(mesh, interfaces, mortars, boundaries, + mesh.boundary_names; mpi_mesh_info = mpi_mesh_info) + + init_mpi_cache!(mpi_cache, mesh, mpi_mesh_info, nvariables(equations), nnodes(dg), + eltype(elements)) + + empty!(mpi_mesh_info.global_mortar_ids) + empty!(mpi_mesh_info.global_interface_ids) + empty!(mpi_mesh_info.neighbor_ranks_mortar) + empty!(mpi_mesh_info.neighbor_ranks_interface) + + # Re-initialize and distribute normal directions of MPI mortars; requires + # MPI communication, so the MPI cache must be re-initialized beforehand. + init_normal_directions!(mpi_mortars, dg.basis, elements) + exchange_normal_directions!(mpi_mortars, mpi_cache, mesh, nnodes(dg)) + + return nothing +end + +# Compatibility to `dgsem_p4est/containers.jl`. +function init_mpi_interfaces!(interfaces, mesh::ParallelT8codeMesh) + # Do nothing. + return nothing +end + +# Compatibility to `dgsem_p4est/containers.jl`. +function init_mpi_mortars!(mortars, mesh::ParallelT8codeMesh) + # Do nothing. + return nothing +end + +# Compatibility to `dgsem_p4est/containers_parallel.jl`. +function init_mpi_mortars!(mpi_mortars, mesh::ParallelT8codeMesh, basis, elements) + # Do nothing. + return nothing +end diff --git a/src/solvers/dgsem_t8code/dg.jl b/src/solvers/dgsem_t8code/dg.jl index 6e9660c917d..e01b12e0f80 100644 --- a/src/solvers/dgsem_t8code/dg.jl +++ b/src/solvers/dgsem_t8code/dg.jl @@ -13,8 +13,8 @@ function create_cache(mesh::T8codeMesh, equations::AbstractEquations, dg::DG, :: boundaries = init_boundaries(mesh, equations, dg.basis, elements) mortars = init_mortars(mesh, equations, dg.basis, elements) - trixi_t8_fill_mesh_info(mesh.forest, elements, interfaces, mortars, boundaries, - mesh.boundary_names) + fill_mesh_info!(mesh, interfaces, mortars, boundaries, + mesh.boundary_names) cache = (; elements, interfaces, boundaries, mortars) @@ -29,4 +29,7 @@ end include("containers.jl") include("containers_2d.jl") include("containers_3d.jl") + +include("containers_parallel.jl") +include("dg_parallel.jl") end # @muladd diff --git a/src/solvers/dgsem_t8code/dg_parallel.jl b/src/solvers/dgsem_t8code/dg_parallel.jl new file mode 100644 index 00000000000..ece614b7d75 --- /dev/null +++ b/src/solvers/dgsem_t8code/dg_parallel.jl @@ -0,0 +1,135 @@ +@muladd begin +#! format: noindent + +# This method is called when a `SemidiscretizationHyperbolic` is constructed. +# It constructs the basic `cache` used throughout the simulation to compute +# the RHS etc. +function create_cache(mesh::ParallelT8codeMesh, equations::AbstractEquations, dg::DG, + ::Any, + ::Type{uEltype}) where {uEltype <: Real} + # Make sure to balance and partition the forest before creating any + # containers in case someone has tampered with forest after creating the + # mesh. + balance!(mesh) + partition!(mesh) + + count_required_surfaces!(mesh) + + elements = init_elements(mesh, equations, dg.basis, uEltype) + mortars = init_mortars(mesh, equations, dg.basis, elements) + interfaces = init_interfaces(mesh, equations, dg.basis, elements) + boundaries = init_boundaries(mesh, equations, dg.basis, elements) + + mpi_mortars = init_mpi_mortars(mesh, equations, dg.basis, elements) + mpi_interfaces = init_mpi_interfaces(mesh, equations, dg.basis, elements) + + mpi_mesh_info = (mpi_mortars = mpi_mortars, + mpi_interfaces = mpi_interfaces, + global_mortar_ids = fill(UInt64(0), nmpimortars(mpi_mortars)), + global_interface_ids = fill(UInt64(0), + nmpiinterfaces(mpi_interfaces)), + neighbor_ranks_mortar = Vector{Vector{Int}}(undef, + nmpimortars(mpi_mortars)), + neighbor_ranks_interface = fill(-1, + nmpiinterfaces(mpi_interfaces))) + + fill_mesh_info!(mesh, interfaces, mortars, boundaries, + mesh.boundary_names; mpi_mesh_info = mpi_mesh_info) + + mpi_cache = init_mpi_cache(mesh, mpi_mesh_info, nvariables(equations), nnodes(dg), + uEltype) + + empty!(mpi_mesh_info.global_mortar_ids) + empty!(mpi_mesh_info.global_interface_ids) + empty!(mpi_mesh_info.neighbor_ranks_mortar) + empty!(mpi_mesh_info.neighbor_ranks_interface) + + init_normal_directions!(mpi_mortars, dg.basis, elements) + exchange_normal_directions!(mpi_mortars, mpi_cache, mesh, nnodes(dg)) + + cache = (; elements, interfaces, mpi_interfaces, boundaries, mortars, mpi_mortars, + mpi_cache) + + # Add specialized parts of the cache required to compute the volume integral etc. + cache = (; cache..., + create_cache(mesh, equations, dg.volume_integral, dg, uEltype)...) + cache = (; cache..., create_cache(mesh, equations, dg.mortar, uEltype)...) + + return cache +end + +function init_mpi_cache(mesh::ParallelT8codeMesh, mpi_mesh_info, nvars, nnodes, uEltype) + mpi_cache = P4estMPICache(uEltype) + init_mpi_cache!(mpi_cache, mesh, mpi_mesh_info, nvars, nnodes, uEltype) + return mpi_cache +end + +function init_mpi_cache!(mpi_cache::P4estMPICache, mesh::ParallelT8codeMesh, + mpi_mesh_info, nvars, nnodes, uEltype) + mpi_neighbor_ranks, mpi_neighbor_interfaces, mpi_neighbor_mortars = init_mpi_neighbor_connectivity(mpi_mesh_info, + mesh) + + mpi_send_buffers, mpi_recv_buffers, mpi_send_requests, mpi_recv_requests = init_mpi_data_structures(mpi_neighbor_interfaces, + mpi_neighbor_mortars, + ndims(mesh), + nvars, + nnodes, + uEltype) + + n_elements_global = Int(t8_forest_get_global_num_elements(mesh.forest)) + n_elements_local = Int(t8_forest_get_local_num_elements(mesh.forest)) + + n_elements_by_rank = Vector{Int}(undef, mpi_nranks()) + n_elements_by_rank[mpi_rank() + 1] = n_elements_local + + MPI.Allgather!(MPI.UBuffer(n_elements_by_rank, 1), mpi_comm()) + + n_elements_by_rank = OffsetArray(n_elements_by_rank, 0:(mpi_nranks() - 1)) + + # Account for 1-based indexing in Julia. + first_element_global_id = sum(n_elements_by_rank[0:(mpi_rank() - 1)]) + 1 + + @assert n_elements_global==sum(n_elements_by_rank) "error in total number of elements" + + @pack! mpi_cache = mpi_neighbor_ranks, mpi_neighbor_interfaces, + mpi_neighbor_mortars, + mpi_send_buffers, mpi_recv_buffers, + mpi_send_requests, mpi_recv_requests, + n_elements_by_rank, n_elements_global, + first_element_global_id + + return mpi_cache +end + +function init_mpi_neighbor_connectivity(mpi_mesh_info, mesh::ParallelT8codeMesh) + @unpack mpi_interfaces, mpi_mortars, global_interface_ids, neighbor_ranks_interface, global_mortar_ids, neighbor_ranks_mortar = mpi_mesh_info + + mpi_neighbor_ranks = vcat(neighbor_ranks_interface, neighbor_ranks_mortar...) |> + sort |> unique + + p = sortperm(global_interface_ids) + + neighbor_ranks_interface .= neighbor_ranks_interface[p] + interface_ids = collect(1:nmpiinterfaces(mpi_interfaces))[p] + + p = sortperm(global_mortar_ids) + neighbor_ranks_mortar .= neighbor_ranks_mortar[p] + mortar_ids = collect(1:nmpimortars(mpi_mortars))[p] + + # For each neighbor rank, init connectivity data structures + mpi_neighbor_interfaces = Vector{Vector{Int}}(undef, length(mpi_neighbor_ranks)) + mpi_neighbor_mortars = Vector{Vector{Int}}(undef, length(mpi_neighbor_ranks)) + for (index, d) in enumerate(mpi_neighbor_ranks) + mpi_neighbor_interfaces[index] = interface_ids[findall(==(d), + neighbor_ranks_interface)] + mpi_neighbor_mortars[index] = mortar_ids[findall(x -> (d in x), + neighbor_ranks_mortar)] + end + + # Check that all interfaces were counted exactly once + @assert mapreduce(length, +, mpi_neighbor_interfaces; init = 0) == + nmpiinterfaces(mpi_interfaces) + + return mpi_neighbor_ranks, mpi_neighbor_interfaces, mpi_neighbor_mortars +end +end # @muladd diff --git a/src/solvers/dgsem_tree/dg_2d_parallel.jl b/src/solvers/dgsem_tree/dg_2d_parallel.jl index 8095dae123a..157d462aa2f 100644 --- a/src/solvers/dgsem_tree/dg_2d_parallel.jl +++ b/src/solvers/dgsem_tree/dg_2d_parallel.jl @@ -446,7 +446,8 @@ function init_mpi_neighbor_connectivity(elements, mpi_interfaces, mpi_mortars, end function rhs!(du, u, t, - mesh::Union{ParallelTreeMesh{2}, ParallelP4estMesh{2}}, equations, + mesh::Union{ParallelTreeMesh{2}, ParallelP4estMesh{2}, + ParallelT8codeMesh{2}}, equations, initial_condition, boundary_conditions, source_terms::Source, dg::DG, cache) where {Source} # Start to receive MPI data diff --git a/test/test_mpi.jl b/test/test_mpi.jl index ad1ba4e835d..1ab1282b891 100644 --- a/test/test_mpi.jl +++ b/test/test_mpi.jl @@ -19,10 +19,12 @@ CI_ON_WINDOWS = (get(ENV, "GITHUB_ACTIONS", false) == "true") && Sys.iswindows() # TreeMesh tests include("test_mpi_tree.jl") - # P4estMesh tests + # P4estMesh and T8codeMesh tests include("test_mpi_p4est_2d.jl") + include("test_mpi_t8code_2d.jl") if !CI_ON_WINDOWS # see comment on `CI_ON_WINDOWS` above include("test_mpi_p4est_3d.jl") + include("test_mpi_t8code_3d.jl") end end # MPI diff --git a/test/test_mpi_p4est_2d.jl b/test/test_mpi_p4est_2d.jl index da90537fcfd..6d66bc68a26 100644 --- a/test/test_mpi_p4est_2d.jl +++ b/test/test_mpi_p4est_2d.jl @@ -33,6 +33,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") @test errors.linf≈[0.00011787417954578494] rtol=1.0e-4 end end + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_nonconforming_flag.jl" begin @@ -40,6 +49,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") "elixir_advection_nonconforming_flag.jl"), l2=[3.198940059144588e-5], linf=[0.00030636069494005547]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_unstructured_flag.jl" begin @@ -47,6 +65,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") "elixir_advection_unstructured_flag.jl"), l2=[0.0005379687442422346], linf=[0.007438525029884735]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_amr_solution_independent.jl" begin @@ -56,6 +83,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") l2=[4.949660644033807e-5], linf=[0.0004867846262313763], coverage_override=(maxiters = 6,)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_amr_unstructured_flag.jl" begin @@ -64,6 +100,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") l2=[0.0012766060609964525], linf=[0.01750280631586159], coverage_override=(maxiters = 6,)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_restart.jl" begin @@ -73,6 +118,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") # With the default `maxiters = 1` in coverage tests, # there would be no time steps after the restart. coverage_override=(maxiters = 100_000,)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_euler_source_terms_nonconforming_unstructured_flag.jl" begin @@ -90,6 +144,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_2d_dgsem") 0.03759938693042297, 0.08039824959535657, ]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end end end # P4estMesh MPI diff --git a/test/test_mpi_p4est_3d.jl b/test/test_mpi_p4est_3d.jl index 75f43650082..cca9093ec51 100644 --- a/test/test_mpi_p4est_3d.jl +++ b/test/test_mpi_p4est_3d.jl @@ -33,6 +33,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") @test errors.linf≈[0.0014548839020096516] rtol=1.0e-4 end end + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_amr.jl" begin @@ -46,6 +55,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") initial_refinement_level = 2, base_level = 2, med_level = 3, max_level = 4)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_amr_unstructured_curved.jl" begin @@ -58,6 +76,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") initial_refinement_level = 0, base_level = 0, med_level = 1, max_level = 2)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_restart.jl" begin @@ -67,12 +94,30 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") # With the default `maxiters = 1` in coverage tests, # there would be no time steps after the restart. coverage_override=(maxiters = 100_000,)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_advection_cubed_sphere.jl" begin @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_advection_cubed_sphere.jl"), l2=[0.002006918015656413], linf=[0.027655117058380085]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end # Compressible Euler @@ -94,6 +139,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") 0.008526972236273522, ], tspan=(0.0, 0.01)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_euler_source_terms_nonperiodic.jl" begin @@ -114,6 +168,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") 0.01562861968368434, ], tspan=(0.0, 1.0)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_euler_ec.jl" begin @@ -134,6 +197,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") ], tspan=(0.0, 0.2), coverage_override=(polydeg = 3,)) # Prevent long compile time in CI + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end @trixi_testset "elixir_euler_source_terms_nonperiodic_hohqmesh.jl" begin @@ -153,6 +225,15 @@ const EXAMPLES_DIR = pkgdir(Trixi, "examples", "p4est_3d_dgsem") 0.048396544302230504, 0.1154589758186293, ]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end end end end # P4estMesh MPI diff --git a/test/test_mpi_t8code_2d.jl b/test/test_mpi_t8code_2d.jl new file mode 100644 index 00000000000..7c7fc03898c --- /dev/null +++ b/test/test_mpi_t8code_2d.jl @@ -0,0 +1,142 @@ +module TestExamplesMPIT8codeMesh2D + +using Test +using Trixi + +include("test_trixi.jl") + +const EXAMPLES_DIR = pkgdir(Trixi, "examples", "t8code_2d_dgsem") + +@testset "T8codeMesh MPI 2D" begin +#! format: noindent + +# Run basic tests +@testset "Examples 2D" begin + # Linear scalar advection + @trixi_testset "elixir_advection_basic.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_advection_basic.jl"), + # Expected errors are exactly the same as with TreeMesh! + l2=[8.311947673061856e-6], + linf=[6.627000273229378e-5]) + + @testset "error-based step size control" begin + Trixi.mpi_isroot() && println("-"^100) + Trixi.mpi_isroot() && + println("elixir_advection_basic.jl with error-based step size control") + + sol = solve(ode, RDPK3SpFSAL35(); abstol = 1.0e-4, reltol = 1.0e-4, + ode_default_options()..., callback = callbacks) + summary_callback() + errors = analysis_callback(sol) + if Trixi.mpi_isroot() + @test errors.l2≈[3.3022040342579066e-5] rtol=1.0e-4 + @test errors.linf≈[0.00011787417954578494] rtol=1.0e-4 + end + end + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_advection_nonconforming_flag.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_advection_nonconforming_flag.jl"), + l2=[3.198940059144588e-5], + linf=[0.00030636069494005547]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_advection_unstructured_flag.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_advection_unstructured_flag.jl"), + l2=[0.0005379687442422346], + linf=[0.007438525029884735]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_advection_amr_solution_independent.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_advection_amr_solution_independent.jl"), + # Expected errors are exactly the same as with TreeMesh! + l2=[4.933027431215839e-5], + linf=[0.00048678461161243136], + coverage_override=(maxiters = 6,)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_advection_amr_unstructured_flag.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_advection_amr_unstructured_flag.jl"), + l2=[0.001980652042312077], + linf=[0.0328882442132265], + coverage_override=(maxiters = 6,)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_euler_source_terms_nonconforming_unstructured_flag.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_euler_source_terms_nonconforming_unstructured_flag.jl"), + l2=[ + 0.0034516244508588046, + 0.0023420334036925493, + 0.0024261923964557187, + 0.004731710454271893, + ], + linf=[ + 0.04155789011775046, + 0.024772109862748914, + 0.03759938693042297, + 0.08039824959535657, + ]) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end +end +end # T8codeMesh MPI + +end # module diff --git a/test/test_mpi_t8code_3d.jl b/test/test_mpi_t8code_3d.jl new file mode 100644 index 00000000000..a15690a7629 --- /dev/null +++ b/test/test_mpi_t8code_3d.jl @@ -0,0 +1,180 @@ +module TestExamplesMPIT8codeMesh3D + +using Test +using Trixi + +include("test_trixi.jl") + +const EXAMPLES_DIR = pkgdir(Trixi, "examples", "t8code_3d_dgsem") + +@testset "T8codeMesh MPI 3D" begin +#! format: noindent + +# Run basic tests +@testset "Examples 3D" begin + # Linear scalar advection + @trixi_testset "elixir_advection_basic.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_advection_basic.jl"), + # Expected errors are exactly the same as with TreeMesh! + l2=[0.00016263963870641478], + linf=[0.0014537194925779984]) + + @testset "error-based step size control" begin + Trixi.mpi_isroot() && println("-"^100) + Trixi.mpi_isroot() && + println("elixir_advection_basic.jl with error-based step size control") + + sol = solve(ode, RDPK3SpFSAL35(); abstol = 1.0e-4, reltol = 1.0e-4, + ode_default_options()..., callback = callbacks) + summary_callback() + errors = analysis_callback(sol) + if Trixi.mpi_isroot() + @test errors.l2≈[0.00016800412839949264] rtol=1.0e-4 + @test errors.linf≈[0.0014548839020096516] rtol=1.0e-4 + end + end + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_advection_amr.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_advection_amr.jl"), + # Expected errors are exactly the same as with TreeMesh! + l2=[1.1302812803902801e-5], + linf=[0.0007889950196294793], + # override values are different from the serial tests to ensure each process holds at least + # one element, otherwise OrdinaryDiffEq fails during initialization + coverage_override=(maxiters = 6, + initial_refinement_level = 2, + base_level = 2, med_level = 3, + max_level = 4)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_advection_amr_unstructured_curved.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_advection_amr_unstructured_curved.jl"), + l2=[2.0556575425846923e-5], + linf=[0.00105682693484822], + tspan=(0.0, 1.0), + coverage_override=(maxiters = 6, + initial_refinement_level = 0, + base_level = 0, med_level = 1, + max_level = 2)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + # Compressible Euler + @trixi_testset "elixir_euler_source_terms_nonconforming_unstructured_curved.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_euler_source_terms_nonconforming_unstructured_curved.jl"), + l2=[ + 4.070355207909268e-5, + 4.4993257426833716e-5, + 5.10588457841744e-5, + 5.102840924036687e-5, + 0.00019986264001630542, + ], + linf=[ + 0.0016987332417202072, + 0.003622956808262634, + 0.002029576258317789, + 0.0024206977281964193, + 0.008526972236273522, + ], + tspan=(0.0, 0.01)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_euler_source_terms_nonperiodic.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, + "elixir_euler_source_terms_nonperiodic.jl"), + l2=[ + 0.0015106060984283647, + 0.0014733349038567685, + 0.00147333490385685, + 0.001473334903856929, + 0.0028149479453087093, + ], + linf=[ + 0.008070806335238156, + 0.009007245083113125, + 0.009007245083121784, + 0.009007245083102688, + 0.01562861968368434, + ], + tspan=(0.0, 1.0)) + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end + + @trixi_testset "elixir_euler_ec.jl" begin + @test_trixi_include(joinpath(EXAMPLES_DIR, "elixir_euler_ec.jl"), + l2=[ + 0.010380390326164493, + 0.006192950051354618, + 0.005970674274073704, + 0.005965831290564327, + 0.02628875593094754, + ], + linf=[ + 0.3326911600075694, + 0.2824952141320467, + 0.41401037398065543, + 0.45574161423218573, + 0.8099577682187109, + ], + tspan=(0.0, 0.2), + coverage_override=(polydeg = 3,)) # Prevent long compile time in CI + + # Ensure that we do not have excessive memory allocations + # (e.g., from type instabilities) + let + t = sol.t[end] + u_ode = sol.u[end] + du_ode = similar(u_ode) + @test (@allocated Trixi.rhs!(du_ode, u_ode, semi, t)) < 1000 + end + end +end +end # T8codeMesh MPI + +end # module diff --git a/test/test_t8code_2d.jl b/test/test_t8code_2d.jl index ab95e068d02..d536a6dd73a 100644 --- a/test/test_t8code_2d.jl +++ b/test/test_t8code_2d.jl @@ -33,10 +33,9 @@ end @trixi_testset "test check_for_negative_volumes" begin @test_warn "Discovered negative volumes" begin # Unstructured mesh with six cells which have left-handed node ordering. - mesh_file = joinpath(EXAMPLES_DIR, "rectangle_with_negative_volumes.msh") - isfile(mesh_file) || - download("https://gist.githubusercontent.com/jmark/bfe0d45f8e369298d6cc637733819013/raw/cecf86edecc736e8b3e06e354c494b2052d41f7a/rectangle_with_negative_volumes.msh", - mesh_file) + mesh_file = Trixi.download("https://gist.githubusercontent.com/jmark/bfe0d45f8e369298d6cc637733819013/raw/cecf86edecc736e8b3e06e354c494b2052d41f7a/rectangle_with_negative_volumes.msh", + joinpath(EXAMPLES_DIR, + "rectangle_with_negative_volumes.msh")) # This call should throw a warning about negative volumes detected. mesh = T8codeMesh(mesh_file, 2)