diff --git a/.github/workflows/idefix-ci-doc.yml b/.github/workflows/idefix-ci-doc.yml new file mode 100644 index 00000000..03640465 --- /dev/null +++ b/.github/workflows/idefix-ci-doc.yml @@ -0,0 +1,24 @@ +name: Build the docs +on: + workflow_dispatch: + push: + branches: + - master + - develop + pull_request: + paths-ignore: + - '.github/ISSUE_TEMPLATE/*' + + +jobs: + ReadTheDocs: + runs-on: ubuntu-latest + steps: + - name: Check out repo + uses: actions/checkout@v3 + - name: install doxygen + run: sudo apt-get install -y doxygen + - name: install python dependencies + run: python -m pip install --exists-action=w --no-cache-dir -r doc/python_requirements.txt + - name: compile documentation + run: python -m sphinx -T -b html -d _build/doctrees -D language=en doc/source doc/html diff --git a/.github/workflows/idefix-ci-jobs.yml b/.github/workflows/idefix-ci-jobs.yml new file mode 100644 index 00000000..c0e62296 --- /dev/null +++ b/.github/workflows/idefix-ci-jobs.yml @@ -0,0 +1,203 @@ +on: + workflow_call: + inputs: + TESTME_OPTIONS: + required: true + type: string + IDEFIX_COMPILER: + required: true + type: string + +# concurrency: +# # auto-cancel any concurrent job *in the same context* +# # see https://docs.github.com/en/actions/using-workflows/workflow-syntax-for-github-actions#concurrency +# # see https://docs.github.com/en/actions/learn-github-actions/contexts#github-context +# group: ${{ github.workflow }}-${{ github.ref }} +# cancel-in-progress: true + +env: + IDEFIX_COMPILER: ${{ inputs.IDEFIX_COMPILER }} + TESTME_OPTIONS: ${{ inputs.TESTME_OPTIONS }} + PYTHONPATH: ${{ github.workspace }} + IDEFIX_DIR: ${{ github.workspace }} + +jobs: + ShocksHydro: + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Sod test + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/sod -all $TESTME_OPTIONS + - name: Isothermal Sod test + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/sod-iso -all $TESTME_OPTIONS + - name: Mach reflection test + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD//MachReflection -all $TESTME_OPTIONS + + ParabolicHydro: + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Viscous flow past cylinder + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/ViscousFlowPastCylinder -all $TESTME_OPTIONS + - name: Viscous disk + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/ViscousDisk -all $TESTME_OPTIONS + - name: Thermal diffusion + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/thermalDiffusion -all $TESTME_OPTIONS + + ShocksMHD: + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: MHD Sod test + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/sod -all $TESTME_OPTIONS + - name: MHD Isothermal Sod test + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/sod-iso -all $TESTME_OPTIONS + - name: Orszag Tang + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/OrszagTang -all $TESTME_OPTIONS + - name: Orszag Tang 3D+restart tests + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/OrszagTang3D -all $TESTME_OPTIONS + + ParabolicMHD: + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Ambipolar C Shock + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/AmbipolarCshock -all $TESTME_OPTIONS + - name: Ambipolar C Shock 3D + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/AmbipolarCshock3D -all $TESTME_OPTIONS + - name: Resistive Alfvén wave + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/ResistiveAlfvenWave -all $TESTME_OPTIONS + - name: Grid coarsening diffusion + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/Coarsening -all $TESTME_OPTIONS + - name: Hall whistler waves + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/HallWhistler -all $TESTME_OPTIONS + + Fargo: + needs: [ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Fargo + planet + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/FargoPlanet -all $TESTME_OPTIONS + - name: Fargo MHD spherical + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/FargoMHDSpherical -all $TESTME_OPTIONS + + ShearingBox: + needs: [ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Hydro shearing box + run: scripts/ci/run-tests $IDEFIX_DIR/test/HD/ShearingBox -all $TESTME_OPTIONS + - name: MHD shearing box + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/ShearingBox -all $TESTME_OPTIONS + + SelfGravity: + needs: [ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Jeans Instability + run: scripts/ci/run-tests $IDEFIX_DIR/test/SelfGravity/JeansInstability -all $TESTME_OPTIONS + - name: Random sphere spherical + run: scripts/ci/run-tests $IDEFIX_DIR/test/SelfGravity/RandomSphere -all $TESTME_OPTIONS + - name: Random sphere cartesian + run: scripts/ci/run-tests $IDEFIX_DIR/test/SelfGravity/RandomSphereCartesian -all $TESTME_OPTIONS + - name: Uniform spherical collapse + run: scripts/ci/run-tests $IDEFIX_DIR/test/SelfGravity/UniformCollapse -all $TESTME_OPTIONS + - name: Dusty spherical collapse + run: scripts/ci/run-tests $IDEFIX_DIR/test/SelfGravity/DustyCollapse -all $TESTME_OPTIONS + + Planet: + needs: [ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: 3 body + run: scripts/ci/run-tests $IDEFIX_DIR/test/Planet/Planet3Body -all $TESTME_OPTIONS + - name: migration + run: scripts/ci/run-tests $IDEFIX_DIR/test/Planet/PlanetMigration2D -all $TESTME_OPTIONS + - name: planet-planet + run: scripts/ci/run-tests $IDEFIX_DIR/test/Planet/PlanetPlanetRK42D -all $TESTME_OPTIONS + - name: spiral wake + run: scripts/ci/run-tests $IDEFIX_DIR/test/Planet/PlanetSpiral2D -all $TESTME_OPTIONS + - name: torques + run: scripts/ci/run-tests $IDEFIX_DIR/test/Planet/PlanetTorque3D -all $TESTME_OPTIONS + - name: RK5 + run: scripts/ci/run-tests $IDEFIX_DIR/test/Planet/PlanetsIsActiveRK52D -all $TESTME_OPTIONS + + Dust: + needs: [ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Energy conservation + run: scripts/ci/run-tests $IDEFIX_DIR/test/Dust/DustEnergy -all $TESTME_OPTIONS + - name: Dusty wave + run: scripts/ci/run-tests $IDEFIX_DIR/test/Dust/DustyWave -all $TESTME_OPTIONS + + Braginskii: + needs: [ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: MTI + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/MTI -all $TESTME_OPTIONS + - name: Spherical anisotropic diffusion + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/sphBragTDiffusion -all $TESTME_OPTIONS + - name: Spherical anisotropic viscosity + run: scripts/ci/run-tests $IDEFIX_DIR/test/MHD/sphBragViscosity -all $TESTME_OPTIONS + + Examples: + needs: [Fargo, Dust, Planet, ShearingBox, SelfGravity] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Run examples test + run: cd test && ./checks_examples.sh $TEST_OPTIONS + + Utils: + needs: [Fargo, Dust, Planet, ShearingBox, SelfGravity] + runs-on: self-hosted + steps: + - name: Check out repo + uses: actions/checkout@v3 + with: + submodules: recursive + - name: Lookup table + run: scripts/ci/run-tests $IDEFIX_DIR/test/utils/lookupTable -all $TESTME_OPTIONS + - name: Dump Image + run: scripts/ci/run-tests $IDEFIX_DIR/test/utils/dumpImage -all $TESTME_OPTIONS diff --git a/.github/workflows/idefix-ci.yml b/.github/workflows/idefix-ci.yml index e8201828..cb946143 100644 --- a/.github/workflows/idefix-ci.yml +++ b/.github/workflows/idefix-ci.yml @@ -4,23 +4,11 @@ on: push: branches: - master - - v2.0 + - develop pull_request: paths-ignore: - '.github/ISSUE_TEMPLATE/*' -concurrency: - # auto-cancel any concurrent job *in the same context* - # see https://docs.github.com/en/actions/using-workflows/workflow-syntax-for-github-actions#concurrency - # see https://docs.github.com/en/actions/learn-github-actions/contexts#github-context - group: ${{ github.workflow }}-${{ github.ref }} - cancel-in-progress: true - -env: - TESTME_OPTIONS: -cuda -Werror - PYTHONPATH: ${{ github.workspace }} - IDEFIX_DIR: ${{ github.workspace }} - jobs: Linter: # Don't do this in forks @@ -35,261 +23,26 @@ jobs: - uses: pre-commit-ci/lite-action@v1.0.0 if: always() - ShocksHydro: + icc-jobs: needs: Linter - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Sod test - run: | - cd $IDEFIX_DIR/test/HD/sod - ./testme.py -all $TESTME_OPTIONS - - name: Isothermal Sod test - run: | - cd $IDEFIX_DIR/test/HD/sod-iso - ./testme.py -all $TESTME_OPTIONS - - name: Mach reflection test - run: | - cd $IDEFIX_DIR/test/HD//MachReflection - ./testme.py -all $TESTME_OPTIONS + name: CPU Jobs (intel OneApi) + uses: ./.github/workflows/idefix-ci-jobs.yml + with: + TESTME_OPTIONS: -intel -Werror + IDEFIX_COMPILER: icc - ParabolicHydro: + gcc-jobs: needs: Linter - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Viscous flow past cylinder - run: | - cd $IDEFIX_DIR/test/HD/ViscousFlowPastCylinder - ./testme.py -all $TESTME_OPTIONS - - name: Viscous disk - run: | - cd $IDEFIX_DIR/test/HD/ViscousDisk - ./testme.py -all $TESTME_OPTIONS - - name: Thermal diffusion - run: | - cd $IDEFIX_DIR/test/HD/thermalDiffusion - ./testme.py -all $TESTME_OPTIONS + name: CPU Jobs (gcc) + uses: ./.github/workflows/idefix-ci-jobs.yml + with: + TESTME_OPTIONS: -Werror + IDEFIX_COMPILER: gcc - ShocksMHD: + cuda-jobs: needs: Linter - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: MHD Sod test - run: | - cd $IDEFIX_DIR/test/MHD/sod - ./testme.py -all $TESTME_OPTIONS - - name: MHD Isothermal Sod test - run: | - cd $IDEFIX_DIR/test/MHD/sod-iso - ./testme.py -all $TESTME_OPTIONS - - name: Orszag Tang - run: | - cd $IDEFIX_DIR/test/MHD/OrszagTang - ./testme.py -all $TESTME_OPTIONS - - name: Orszag Tang 3D+restart tests - run: | - cd $IDEFIX_DIR/test/MHD/OrszagTang3D - ./testme.py -all $TESTME_OPTIONS - - - ParabolicMHD: - needs: Linter - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Ambipolar C Shock - run: | - cd $IDEFIX_DIR/test/MHD/AmbipolarCshock - ./testme.py -all $TESTME_OPTIONS - - name: Ambipolar C Shock 3D - run: | - cd $IDEFIX_DIR/test/MHD/AmbipolarCshock3D - ./testme.py -all $TESTME_OPTIONS - - name: Resistive Alfvén wave - run: | - cd $IDEFIX_DIR/test/MHD/ResistiveAlfvenWave - ./testme.py -all $TESTME_OPTIONS - - name: Grid coarsening diffusion - run: | - cd $IDEFIX_DIR/test/MHD/Coarsening - ./testme.py -all $TESTME_OPTIONS - - name: Hall whistler waves - run: | - cd $IDEFIX_DIR/test/MHD/HallWhistler - ./testme.py -all $TESTME_OPTIONS - - Fargo: - needs: [Linter, ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Fargo + planet - run: | - cd $IDEFIX_DIR/test/HD/FargoPlanet - ./testme.py -all $TESTME_OPTIONS - - name: Fargo MHD spherical - run: | - cd $IDEFIX_DIR/test/MHD/FargoMHDSpherical - ./testme.py -all $TESTME_OPTIONS - - ShearingBox: - needs: [Linter, ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Hydro shearing box - run: | - cd $IDEFIX_DIR/test/HD/ShearingBox - ./testme.py -all $TESTME_OPTIONS - - name: MHD shearing box - run: | - cd $IDEFIX_DIR/test/MHD/ShearingBox - ./testme.py -all $TESTME_OPTIONS - - SelfGravity: - needs: [Linter, ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Jeans Instability - run: | - cd $IDEFIX_DIR/test/SelfGravity/JeansInstability - ./testme.py -all $TESTME_OPTIONS - - name: Random sphere spherical - run: | - cd $IDEFIX_DIR/test/SelfGravity/RandomSphere - ./testme.py -all $TESTME_OPTIONS - - name: Random sphere cartesian - run: | - cd $IDEFIX_DIR/test/SelfGravity/RandomSphereCartesian - ./testme.py -all $TESTME_OPTIONS - - name: Uniform spherical collapse - run: | - cd $IDEFIX_DIR/test/SelfGravity/UniformCollapse - ./testme.py -all $TESTME_OPTIONS - - name: Dusty spherical collapse - run: | - cd $IDEFIX_DIR/test/SelfGravity/DustyCollapse - ./testme.py -all $TESTME_OPTIONS - - Planet: - needs: [Linter, ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: 3 body - run: | - cd $IDEFIX_DIR/test/Planet/Planet3Body - ./testme.py -all $TESTME_OPTIONS - - name: migration - run: | - cd $IDEFIX_DIR/test/Planet/PlanetMigration2D - ./testme.py -all $TESTME_OPTIONS - - name: planet-planet - run: | - cd $IDEFIX_DIR/test/Planet/PlanetPlanetRK42D - ./testme.py -all $TESTME_OPTIONS - - name: spiral wake - run: | - cd $IDEFIX_DIR/test/Planet/PlanetSpiral2D - ./testme.py -all $TESTME_OPTIONS - - name: torques - run: | - cd $IDEFIX_DIR/test/Planet/PlanetTorque3D - ./testme.py -all $TESTME_OPTIONS - - name: RK5 - run: | - cd $IDEFIX_DIR/test/Planet/PlanetsIsActiveRK52D - ./testme.py -all $TESTME_OPTIONS - - Dust: - needs: [Linter, ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Energy conservation - run: | - cd $IDEFIX_DIR/test/Dust/DustEnergy - ./testme.py -all $TESTME_OPTIONS - - name: Dusty wave - run: | - cd $IDEFIX_DIR/test/Dust/DustyWave - ./testme.py -all $TESTME_OPTIONS - - Braginskii: - needs: [Linter, ShocksHydro, ParabolicHydro, ShocksMHD, ParabolicMHD] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: MTI - run: | - cd $IDEFIX_DIR/test/MHD/MTI - ./testme.py -all $TESTME_OPTIONS - - name: Spherical anisotropic diffusion - run: | - cd $IDEFIX_DIR/test/MHD/sphBragTDiffusion - ./testme.py -all $TESTME_OPTIONS - - name: Spherical anisotropic viscosity - run: | - cd $IDEFIX_DIR/test/MHD/sphBragViscosity - ./testme.py -all $TESTME_OPTIONS - - Examples: - needs: [Fargo, Dust, Planet, ShearingBox, SelfGravity] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Run examples test - run: cd test && ./checks_examples.sh $TEST_OPTIONS - - Utils: - needs: [Fargo, Dust, Planet, ShearingBox, SelfGravity] - runs-on: self-hosted - steps: - - name: Check out repo - uses: actions/checkout@v3 - with: - submodules: recursive - - name: Lookup table - run: | - cd $IDEFIX_DIR/test/utils/lookupTable - ./testme.py -all $TESTME_OPTIONS - - name: Dump Image - run: | - cd $IDEFIX_DIR/test/utils/dumpImage - ./testme.py -all $TESTME_OPTIONS + name: CUDA Jobs + uses: ./.github/workflows/idefix-ci-jobs.yml + with: + TESTME_OPTIONS: -cuda -Werror + IDEFIX_COMPILER: nvcc diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 726f0045..72360011 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -34,6 +34,10 @@ repos: - B # flake8-bugbear - I # isort - NPY # numpy-specific rules + - --ignore + - F405 + - --ignore + - F403 # ignore import * - repo: https://github.com/neutrinoceros/inifix rev: v5.0.2 diff --git a/CHANGELOG.md b/CHANGELOG.md index 8d791f8f..66f01f7b 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,6 +4,21 @@ All notable changes to this project will be documented in this file. The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/), and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). +## [2.2.00] 2025-01-17 +### Changed + +- Fix a bug that could lead to a segmentation fault when MHD+Fargo+MPI (with X3 domain decomposition) are all enabled (#295) +- Fix a bug that could result in an incorrect magnetic field when initialising B from the vector potential in non-axisymmetric spherical geometry (#293) +- Fix a bug that could result in Idefix believing the MPI library is not Cuda aware for some versions of OpenMPI (#310) +- Ensure that the behaviour in 1D Spherical geometry is identical to Pluto (#291) + +### Added + +- Add the python interface "pydefix", allowing users to initialise and analyse Idefix simulations live from Python without writing any file (#277) +- Add the native Idefix coordinates in VTK file to simplify postprocessing (#292) +- Add code testing on CPU targets using gcc and Intel Oneapi (#300) + + ## [2.1.02] 2024-10-24 ### Changed diff --git a/CMakeLists.txt b/CMakeLists.txt index 8df0de73..a8c7ff94 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -5,16 +5,17 @@ endif() set (CMAKE_CXX_STANDARD 17) set(Idefix_VERSION_MAJOR 2) -set(Idefix_VERSION_MINOR 1) -set(Idefix_VERSION_PATCH 02) +set(Idefix_VERSION_MINOR 2) +set(Idefix_VERSION_PATCH 00) -project (idefix VERSION 2.1.02) +project (idefix VERSION 2.2.00) option(Idefix_MHD "enable MHD" OFF) option(Idefix_MPI "enable Message Passing Interface parallelisation" OFF) option(Idefix_HIGH_ORDER_FARGO "Force Fargo to use a PPM reconstruction scheme" OFF) option(Idefix_DEBUG "Enable Idefix debug features (makes the code very slow)" OFF) option(Idefix_RUNTIME_CHECKS "Enable runtime sanity checks" OFF) option(Idefix_WERROR "Treat compiler warnings as errors" OFF) +option(Idefix_PYTHON "Enable python bindings (requires pybind11)" OFF) set(Idefix_CXX_FLAGS "" CACHE STRING "Additional compiler/linker flag") set(Idefix_DEFS "definitions.hpp" CACHE FILEPATH "Problem definition header file") option(Idefix_CUSTOM_EOS "Use custom equation of state" OFF) @@ -45,7 +46,8 @@ include(SetRequiredBuildSettingsForGCC8) #Idefix requires Cuda Lambdas (experimental) if(Kokkos_ENABLE_CUDA) set(Kokkos_ENABLE_CUDA_LAMBDA ON CACHE BOOL "Idefix requires lambdas on Cuda" FORCE) - set(Kokkos_ENABLE_IMPL_CUDA_MALLOC_ASYNC OFF CACHE BOOL "Disable Async malloc to avoid bugs on PSM2" FORCE) + # CUDA_MALLOC_ASYNC disbaled by default in Kokkos 4.5, so not required here + #set(Kokkos_ENABLE_IMPL_CUDA_MALLOC_ASYNC OFF CACHE BOOL "Disable Async malloc to avoid bugs on PSM2" FORCE) endif() # Add kokkos CMAKE files (required early since these set compiler options) @@ -76,7 +78,7 @@ add_subdirectory(src build) if(EXISTS ${PROJECT_BINARY_DIR}/setup.cpp) target_sources(idefix PUBLIC ${PROJECT_BINARY_DIR}/setup.cpp) else() - message(WARNING "No specific setup.cpp found in the problem directory") + message(WARNING "No specific setup.cpp found in the problem directory (this message can be ignored if using python to define your problem)") endif() # If a CMakeLists.txt is in the problem dir (for problem-specific source files) @@ -114,12 +116,26 @@ if(Idefix_HDF5) ) find_package(HDF5 REQUIRED) target_link_libraries(idefix "${HDF5_LIBRARIES}") - target_include_directories(idefix PUBLIC "${HDF5_INCLUDE_DIRS}") + target_include_directories(idefix "${HDF5_INCLUDE_DIRS}") message(STATUS "XDMF (hdf5+xmf) dumps enabled") else() set(Idefix_HDF5 OFF) endif() +if(Idefix_PYTHON) + add_compile_definitions("WITH_PYTHON") + if (NOT DEFINED Python_FIND_FRAMEWORK) + set(Python_FIND_FRAMEWORK "LAST") # Use Apple's python only at last resort on Macos + endif () + set(PYBIND11_FINDPYTHON ON CACHE BOOL "Idefix requires python" FORCE) + find_package(pybind11 REQUIRED) + target_link_libraries(idefix pybind11::embed) + target_sources(idefix + PUBLIC src/pydefix.cpp + PUBLIC src/pydefix.hpp + ) +endif() + if(Idefix_DEBUG) add_compile_definitions("DEBUG") set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -O0 -g") @@ -231,6 +247,7 @@ else() endif() message(STATUS " MPI: ${Idefix_MPI}") message(STATUS " HDF5: ${Idefix_HDF5}") +message(STATUS " Python: ${Idefix_PYTHON}") message(STATUS " Reconstruction: ${Idefix_RECONSTRUCTION}") message(STATUS " Precision: ${Idefix_PRECISION}") message(STATUS " Version: ${Idefix_VERSION}") diff --git a/doc/source/conf.py b/doc/source/conf.py index b712d026..24f9709b 100644 --- a/doc/source/conf.py +++ b/doc/source/conf.py @@ -23,7 +23,7 @@ author = 'Geoffroy Lesur' # The full version, including alpha/beta/rc tags -release = '2.1.02' +release = '2.2.00' diff --git a/doc/source/index.rst b/doc/source/index.rst index fa8dadb3..128d591a 100644 --- a/doc/source/index.rst +++ b/doc/source/index.rst @@ -33,6 +33,10 @@ MPI library When using MPI parallelisation, *Idefix* relies on an external MPI library. *Idefix* has been tested successfully with OpenMPI and IntelMPI libraries. When used on GPU architectures, *Idefix* assumes that the MPI library is GPU-Aware. If unsure, check this last point with your system administrator. +Python + When using *Idefix* with its python interface through the module `Pydefix`, *Idefix* relies on an external python>=3.8 interpreter with the module `pybind11 `_ + installed. + ================ Features ================ diff --git a/doc/source/modules.rst b/doc/source/modules.rst index 13287ece..2f1f3787 100644 --- a/doc/source/modules.rst +++ b/doc/source/modules.rst @@ -28,6 +28,9 @@ In this section, you will find a more detailed documentation about each module t :ref:`gridCoarseningModule` The grid coarsening module, that allows to derefine the grid in particular locations to speed up the computation. +:ref:`pydefixModule` + The Python interaction module, that allows Idefix to interact directly with a python interpreter. + .. toctree:: :maxdepth: 2 @@ -40,3 +43,4 @@ In this section, you will find a more detailed documentation about each module t modules/selfGravity.rst modules/braginskii.rst modules/gridCoarsening.rst + modules/pydefix.rst diff --git a/doc/source/modules/pydefix.rst b/doc/source/modules/pydefix.rst new file mode 100644 index 00000000..f6828420 --- /dev/null +++ b/doc/source/modules/pydefix.rst @@ -0,0 +1,126 @@ +.. _pydefixModule: + +Pydefix module +============== + +The Pydefix module allows Idefix to talk directly to a Python interpreter while running. It can be used to create your own initial condition +and/or for custom outputs produced from Python. Pydefix relies on the pybind11 python package + +The essence of Pydefix is to allows the user to have a direct access to Idefix data structure from Python without writing/accessing any file. In particular, IdefixArrays are viewed as numpy arrays in Python. +Note however that to keep things simple, Pydefix works on the host memory space only, and hence sync data to/from the GPU (if used) before invoking Python functions. Hence, using Pydefix for outputs induces +a significant loss of performances. + + +Before you start +---------------- +Pybind11 installation ++++++++++++++++++++++ + +In order to use pydefix, you need to be working in a python>=3.8 environement that includes `pybind11 `_. Follow the instruction of your package manager to install pybind11>=2.12. + +Pydefix usage +------------- +Idefix Configuration +++++++++++++++++++++ + +In order to use Pydefix, you need to switch on ``Idefix_PYTHON`` in cmake. This will auto-detect Python and check that pybind11 can be used effectively. + + +Run Idefix with Pydefix ++++++++++++++++++++++++ + +Pydefix is typically enabled from your input file `idefix.ini` in the block ``[Python]``. The following parameters are available: + ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ +| Entry name | Parameter type | Comment | ++========================+=======================+===========================================================================================================+ +| script | string | | (Mandatory) Filename (*without ".py"!*) of the python script that Idefix should use. | +| | | | The script should be in the same location as the Idefix executable file. | ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ +| output_function | string | | (Optional) Name of the function that will be called for each output event (the function should be | +| | | | defined in the python script above). When ommited, pydefix output functions are disabled. | +| | | | The periodicity of the pydefix output routine is set in the block:entry `[Output]:python` | ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ +| initflow_function | string | | (Optional) Name of the python function that will be called to initialize the flow in place of the C++ | +| | | | function ``Setup::InitFlow``. Revert to ``Setup::Initflow`` when ommited. | ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ + +Python script ++++++++++++++ + +When using Pydefix, idefix expects a python script to be specified in the input file (see ``script`` parameter above). To be fully functionnal, you should import the ``pydefix`` module at the beginning +of your python script (you can also import other python module, as any python script). + +Your python script should define functions that will be called while Idefix is running: +* The signature of the ``initflow`` function should be ``(data)`` where ``data`` is a python structure matching Idefix's ``DataBlockHost`` class. +* The signature of the ``output`` function should be ``(data,grid,n)`` where ``data`` is a python structure matching Idefix's ``DataBlockHost`` class, ``grid`` is Idefix's ``GridHost`` class, and ``n`` is an integer representing the current number of the output + +MPI parallelism ++++++++++++++++ +When Idefix runs with MPI parallelism enabled, a python interpreter and script is launched by each MPI process. Each of these script is independent +and have access to its local ``dataBlockHost``. The `pydefix` module however gives access to the local rank ``prank`` and total MPI size ``psize``. In addition, +pydefix provides the function ``GatherIdefixArray`` to gather the data distributed among each process without invoking MPI directly in python. This function +expects a 3D distributed IdefixArray in entry following the signature + +.. code-block:: c++ + + GatherIdefixArray(IdefixHostArray3D in, // 3D distributed array + DataBlockHost dataHost, // dataBlock structure + bool keepBoundaries = true, // Whether we keep the ghost zones in the returned array + bool broadcast = true) // Whether the returned array is available only in proc #0 or in every proc (caution! possibly requires lots of memory) + +This function is used as follows: + +.. code-block:: python + + import pydefix as pdfx # mandatory + import numpy as np + import matplotlib.pyplot as plt + + def output(data,grid,n): + # Gather pressure field from every process in process #0 (set broadcast to True to distribute the result to all processes) + prs = pdfx.GatherIdefixArray(data.Vc[pdfx.PRS,:,:,:],data,broadcast=False) + + # Only root process performs this + if pdfx.prank==0: + x=grid.x[pdfx.IDIR] # The grid contains the full domain size, while the datablock contains only local information + y=grid.x[pdfx.JDIR] + plt.figure() + plt.pcolormesh(x,y,prs[0,:,:],cmap='plasma') + + +.. note:: + For more advanced usage, it is also possibly to directly call MPI routines from python using the `Mpi4py `_ module. + +Example ++++++++ + +An example is provided in the directory `test/IO/python`. This example shows how to use Idefix with pure python initial conditions and outputs. +It reproduces the 2D OrszagTang vortex available in MHD/OrszagTang without requiring any single line of C++ code from the user. + +The python script `pydefix_example.py` initializes the initial condition of the OT test (``initflow``) and produces a series of PNG files through matplotlib (`output`). + +Troubleshooting +--------------- + +It during configuration stage, you get:: + + CMake Error at CMakeLists.txt:122 (find_package): + By not providing "Findpybind11.cmake" in CMAKE_MODULE_PATH this project has + asked CMake to find a package configuration file provided by "pybind11", + but CMake did not find one. + +It means that cmake cannot find the location of pybind11 (this typically happens on MacOs). In order to locate pybind11, open a python interpreter, and get pybind11 install dir through: + +.. code-block:: python + + import pybind11 + print(pybind11.__file__) + +You can then exit the interpreter and set the pybind11_DIR environement variable to the right path: + +.. code-block:: bash + + export pybind11_DIR=env/lib/python3.10/site-packages/pybind11 + +you can then run cmake which should be able to find pybind11, and compile the code. diff --git a/doc/source/reference/idefix.ini.rst b/doc/source/reference/idefix.ini.rst index a7309f10..c9f9c578 100644 --- a/doc/source/reference/idefix.ini.rst +++ b/doc/source/reference/idefix.ini.rst @@ -358,6 +358,23 @@ and ``X1-end``, ``X2-end``, ``X3-end`` for the right boundaries. ``X2`` boundari | | | (see :ref:`userdefBoundaries`) | +----------------+------------------------------------------------------------------------------------------------------------------+ +``Python`` section +------------------ + +This section describes the python script and function that can interact with Idefix while running using the Pydefix module (see :ref:`pydefixModule`) + ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ +| Entry name | Parameter type | Comment | ++========================+=======================+===========================================================================================================+ +| script | string | | (Mandatory) Filename (*without ".py"!*) of the python script that Idefix should use. | +| | | | The script should be in location of Idefix executable file | ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ +| output_function | string | | (Optional) Name of the function that will be called for each output event (the function should be | +| | | | defined in the python script above). When ommited, pydefix output functions are disabled. | ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ +| initflow_function | string | | (optional) Name of the python function that will be called to initialize the flow in place of the C++ | +| | | | function `Setup::InitFlow`. Revert to `Setup::Initflow`` when ommited. | ++------------------------+-----------------------+-----------------------------------------------------------------------------------------------------------+ .. _outputSection: @@ -411,6 +428,9 @@ This section describes the outputs *Idefix* produces. For more details about eac | | | | (see :ref:`functionEnrollment`). The user-defined variables defined by this function | | | | | are then written as new variables in vtk and/or xdmf outputs. | +----------------+-------------------------+--------------------------------------------------------------------------------------------------+ +| python | float | | Time interval between pydefix outputs, in code units. | +| | | | If negative, periodic pydefix outputs are disabled. | ++----------------+-------------------------+--------------------------------------------------------------------------------------------------+ .. note:: Even if dumps are not mentionned in your input file (and are therefore disabled), dump files are still produced when *Idefix* captures a signal diff --git a/doc/source/reference/makefile.rst b/doc/source/reference/makefile.rst index 4f948e5d..aadc205c 100644 --- a/doc/source/reference/makefile.rst +++ b/doc/source/reference/makefile.rst @@ -108,18 +108,11 @@ We recommend the following modules and environement variables on AdAstra: .. code-block:: bash - module load cpe/23.12 + module load cpe/24.07 module load craype-accel-amd-gfx90a craype-x86-trento module load PrgEnv-cray - module load amd-mixed/5.7.1 - module load rocm/5.7.1 # nécessaire a cause d'un bug de path pas encore fix.. - export HIPCC_COMPILE_FLAGS_APPEND="-isystem ${CRAY_MPICH_PREFIX}/include" - export HIPCC_LINK_FLAGS_APPEND="-L${CRAY_MPICH_PREFIX}/lib -lmpi ${PE_MPICH_GTL_DIR_amd_gfx90a} ${PE_MPICH_GTL_LIBS_amd_gfx90a} -lstdc++fs" - export CXX=hipcc - export CC=hipcc - -The `-lstdc++fs` option being there to guarantee the link to the HIP library and the access to specific -C++17 functions. + module load amd-mixed + module load cray-python/3.11.7 Finally, *Idefix* can be configured to run on Mi250 by enabling HIP and the desired architecture with the following options to ccmake: diff --git a/doc/source/reference/outputs.rst b/doc/source/reference/outputs.rst index d531b497..2604d9ab 100644 --- a/doc/source/reference/outputs.rst +++ b/doc/source/reference/outputs.rst @@ -7,7 +7,7 @@ Output formats -------------- *Idefix* uses several types of outputs you may want for your setup. By default, *Idefix* allows -for 4 kinds of outputs: +for 5 kinds of outputs: * logs which essentially tells the user what *Idefix* is currently doing. When running in serial, logs are sent to stdout, but when MPI is enabled, only the logs of the rank 0 process is sent to stdout, and each process (including rank 0) simultaneously writes a @@ -21,17 +21,20 @@ for 4 kinds of outputs: or `Visit `_. The XDMF format relies on the HDF5 format and therefore requires *Idefix* to be configured with HDF5 support. * user-defined analysis files. These are totally left to the user. They usually consist of ascii tables defined by the user, but they can be anything. +* python script, that relies on the :ref:`pydefixModule`. This launches a user-defined python function fed with Idefix data. One can then directly plot or interact with Idefix outputs from python. The output periodicity and the userdef variables should all be declared in the input file, as described in :ref:`outputSection`. Defining your own outputs ------------------------- -*Idefix* provides two ways to define your own outputs: analysis, which are used to make your -own output file (e.g. an ascii-tabulated file); and user variables, which are written by *Idefix* output routines. +*Idefix* provides three ways to define your own outputs: analysis, which are used to make your +own output file (e.g. an ascii-tabulated file); user variables, which are written by *Idefix* output routines, and python user-defined +functions that process Idefix's data. Both analysis and uservar requires the definition of a user function which needs to be enrolled following the procedure described -in :ref:`functionEnrollment` and using the function signatures declared in `output.hpp`. +in :ref:`functionEnrollment` and using the function signatures declared in `output.hpp`. The python outputs are described +in the :ref:`pydefixModule`. We provide below an example of a setup using both analysis outputs and uservar outputs diff --git a/pytools/idfx_test.py b/pytools/idfx_test.py index dc2da0a0..25bac76b 100644 --- a/pytools/idfx_test.py +++ b/pytools/idfx_test.py @@ -58,6 +58,10 @@ def __init__ (self): help="Test on Nvidia GPU using CUDA", action="store_true") + parser.add_argument("-intel", + help="Test compiling with Intel OneAPI", + action="store_true") + parser.add_argument("-hip", help="Test on AMD GPU using HIP", action="store_true") @@ -119,6 +123,12 @@ def configure(self,definitionFile=""): # disable Async cuda malloc for tests performed on old UCX implementations comm.append("-DKokkos_ENABLE_IMPL_CUDA_MALLOC_ASYNC=OFF") + if self.intel: + # disable fmad operations on Cuda to make it compatible with CPU arithmetics + comm.append("-DIdefix_CXX_FLAGS=-fp-model=strict") + comm.append("-DCMAKE_CXX_COMPILER=icpx") + comm.append("-DCMAKE_C_COMPILER=icx") + if self.hip: comm.append("-DKokkos_ENABLE_HIP=ON") # disable fmad operations on HIP to make it compatible with CPU arithmetics diff --git a/pytools/vtk_io.py b/pytools/vtk_io.py index 1f693318..eb5cf412 100644 --- a/pytools/vtk_io.py +++ b/pytools/vtk_io.py @@ -7,7 +7,7 @@ import warnings import numpy as np import os - +import re # restrict what's included with `import *` to public API __all__ = [ @@ -21,6 +21,8 @@ dt = np.dtype(">f") # Big endian single precision floats dint = np.dtype(">i4") # Big endian integer +NATIVE_COORDINATE_REGEXP = re.compile(r"X(1|2|3)(L|C)_NATIVE_COORDINATES") + KNOWN_GEOMETRIES = { 0: "cartesian", 1: "polar", @@ -33,10 +35,14 @@ class VTKDataset(object): def __init__(self, filename, geometry=None): self.filename = os.path.abspath(filename) self.data = {} + self.native_coordinates = {} with open(filename, "rb") as fh: self._load_header(fh, geometry=geometry) self._load(fh) + if self.native_coordinates: + self._setup_coordinates_from_native() + def _load(self, fh): if self._dataset_type in ("RECTILINEAR_GRID", "STRUCTURED_GRID"): self._load_hydro(fh) @@ -88,6 +94,9 @@ def _load_header(self, fh, geometry=None): self.t = np.fromfile(fh, dt, 1) elif d.startswith("PERIODICITY"): self.periodicity = np.fromfile(fh, dtype=dint, count=3).astype(bool) + elif NATIVE_COORDINATE_REGEXP.match(d): + native_name, _ncomp, native_dim, _dtype = d.split() + self.native_coordinates[native_name] = np.fromfile(fh, dtype=dt, count=int(native_dim)) else: warnings.warn("Found unknown field %s" % d) fh.readline() # skip extra linefeed (empty line) @@ -341,6 +350,29 @@ def _load_hydro(self, fh): def _load_particles(self, fh): raise NotImplementedError("Particles vtk are not supported yet !") + def _setup_coordinates_from_native(self): + if self.geometry == "spherical": + native2attr = { + "X1L_NATIVE_COORDINATES": "rl", + "X1C_NATIVE_COORDINATES": "r", + "X2L_NATIVE_COORDINATES": "thetal", + "X2C_NATIVE_COORDINATES": "theta", + "X3L_NATIVE_COORDINATES": "phil", + "X3C_NATIVE_COORDINATES": "phi", + } + elif self.geometry in ("cartesian", "cylindrical", "polar"): + native2attr = { + "X1L_NATIVE_COORDINATES": "xl", + "X1C_NATIVE_COORDINATES": "x", + "X2L_NATIVE_COORDINATES": "yl", + "X2C_NATIVE_COORDINATES": "y", + "X3L_NATIVE_COORDINATES": "zl", + "X3C_NATIVE_COORDINATES": "z", + } + + for native_field, attr in native2attr.items(): + setattr(self, attr, self.native_coordinates[native_field]) + def __repr__(self): return "VTKDataset('%s')" % self.filename diff --git a/reference b/reference index 609a6016..b675bcea 160000 --- a/reference +++ b/reference @@ -1 +1 @@ -Subproject commit 609a60164461fb2d4fbeffa89ce148acf1191525 +Subproject commit b675bceaa6aabc01dded346e2d631857f698dc76 diff --git a/scripts/ci/run-tests b/scripts/ci/run-tests new file mode 100755 index 00000000..32983102 --- /dev/null +++ b/scripts/ci/run-tests @@ -0,0 +1,7 @@ +#!/usr/bin/env bash +if [ "$IDEFIX_COMPILER" == icc ]; then + source /opt/intel/oneapi/setvars.sh +fi +set -ue +cd "$1" +./testme.py "${@:2}" diff --git a/src/dataBlock/dataBlockHost.cpp b/src/dataBlock/dataBlockHost.cpp index b5dff915..d570648f 100644 --- a/src/dataBlock/dataBlockHost.cpp +++ b/src/dataBlock/dataBlockHost.cpp @@ -34,6 +34,9 @@ DataBlockHost::DataBlockHost(DataBlock& datain) { nghost = data->nghost; + lbound = data->lbound; + rbound = data->rbound; + xbeg = data->xbeg; xend = data->xend; beg = data->beg; @@ -95,6 +98,9 @@ DataBlockHost::DataBlockHost(DataBlock& datain) { this->haveplanetarySystem = data->haveplanetarySystem; this->planetarySystem = data->planetarySystem.get(); + this->t = data->t; + this->dt = data->dt; + idfx::popRegion(); } @@ -102,6 +108,8 @@ DataBlockHost::DataBlockHost(DataBlock& datain) { void DataBlockHost::SyncToDevice() { idfx::pushRegion("DataBlockHost::SyncToDevice()"); + data->t = this->t; + data->dt = this->dt; Kokkos::deep_copy(data->hydro->Vc,Vc); Kokkos::deep_copy(data->hydro->InvDt,InvDt); @@ -138,6 +146,9 @@ void DataBlockHost::SyncToDevice() { void DataBlockHost::SyncFromDevice() { idfx::pushRegion("DataBlockHost::SyncFromDevice()"); + this->t = data->t; + this->dt = data->dt; + Kokkos::deep_copy(Vc,data->hydro->Vc); Kokkos::deep_copy(InvDt,data->hydro->InvDt); @@ -223,8 +234,8 @@ void DataBlockHost::MakeVsFromAmag(IdefixHostArray4D &Ain) { + 1/(x1m(i)*(cos(x2m(j)) - cos(x2m(j+1)))) * (sin(x2m(j+1))*Ain(KDIR,k,j+1,i) - sin(x2m(j))*Ain(KDIR,k,j,i) ) , - - 1/(x1m(i)*sin(x2(j))*dx3(k)) * (Ain(JDIR,k+1,j,i) - - Ain(JDIR,k,j,i) ) ); + - dx2(j)/(x1m(i)*((cos(x2m(j))- cos(x2m(j+1))))*dx3(k)) + * (Ain(JDIR,k+1,j,i) - Ain(JDIR,k,j,i) ) ); real Ax2m = fabs(sin(x2m(j))); // Regularisation along the axis diff --git a/src/dataBlock/dataBlockHost.hpp b/src/dataBlock/dataBlockHost.hpp index 86c8c570..f7dc4041 100644 --- a/src/dataBlock/dataBlockHost.hpp +++ b/src/dataBlock/dataBlockHost.hpp @@ -25,33 +25,33 @@ class DataBlockHost { public: // Local grid information - std::array::HostMirror,3> x; ///< geometrical central points - std::array::HostMirror,3> xr; ///< cell right interface - std::array::HostMirror,3> xl; ///< cell left interface - std::array::HostMirror,3> dx; ///< cell width + std::array,3> x; ///< geometrical central points + std::array,3> xr; ///< cell right interface + std::array,3> xl; ///< cell left interface + std::array,3> dx; ///< cell width - IdefixArray3D::HostMirror dV; ///< cell volume - std::array::HostMirror,3> A; ///< cell right interface area + IdefixHostArray3D dV; ///< cell volume + std::array,3> A; ///< cell right interface area - IdefixArray4D::HostMirror Vc; ///< Main cell-centered primitive variables index + IdefixHostArray4D Vc; ///< Main cell-centered primitive variables index bool haveDust{false}; std::vector> dustVc; ///< Cell-centered primitive variables index for dust #if MHD == YES - IdefixArray4D::HostMirror Vs; ///< Main face-centered primitive variables index - IdefixArray4D::HostMirror Ve; ///< Main edge-centered primitive variables index - IdefixArray4D::HostMirror J; ///< Current (only when haveCurrent is enabled) + IdefixHostArray4D Vs; ///< Main face-centered primitive variables index + IdefixHostArray4D Ve; ///< Main edge-centered primitive variables index + IdefixHostArray4D J; ///< Current (only when haveCurrent is enabled) - IdefixArray3D::HostMirror Ex1; ///< x1 electric field - IdefixArray3D::HostMirror Ex2; ///< x2 electric field - IdefixArray3D::HostMirror Ex3; ///< x3 electric field + IdefixHostArray3D Ex1; ///< x1 electric field + IdefixHostArray3D Ex2; ///< x2 electric field + IdefixHostArray3D Ex3; ///< x3 electric field #endif - IdefixArray4D::HostMirror Uc; ///< Main cell-centered conservative variables - IdefixArray3D::HostMirror InvDt; ///< Inverse of maximum timestep in each cell + IdefixHostArray4D Uc; ///< Main cell-centered conservative variables + IdefixHostArray3D InvDt; ///< Inverse of maximum timestep in each cell - std::array::HostMirror,3> coarseningLevel; ///< Grid coarsening level + std::array,3> coarseningLevel; ///< Grid coarsening level ///< (only defined when coarsening ///< is enabled) @@ -75,6 +75,9 @@ class DataBlockHost { std::array gbeg; ///< Begining of local block in the grid (internal) std::array gend; ///< End of local block in the grid (internal) + real dt; ///< Current timestep + real t; ///< Current time + explicit DataBlockHost(DataBlock &); ///< Constructor from a device datablock ///< (NB: does not sync any data) @@ -95,7 +98,6 @@ class DataBlockHost { GridCoarsening haveGridCoarsening{GridCoarsening::disabled}; ///< Is grid coarsening enabled? - private: // Data object to which we are the mirror DataBlock *data; }; diff --git a/src/dataBlock/dumpToFile.cpp b/src/dataBlock/dumpToFile.cpp index 98b4a474..35b60b81 100644 --- a/src/dataBlock/dumpToFile.cpp +++ b/src/dataBlock/dumpToFile.cpp @@ -92,6 +92,16 @@ void DataBlock::DumpToFile(std::string filebase) { WriteVariable(fileHdl, 4, dims, fieldName, locVc.data()); + if (this->gravity->haveSelfGravityPotential) { + IdefixArray3D::HostMirror locPot = Kokkos::create_mirror_view(this->gravity->phiP); + Kokkos::deep_copy(locPot, this->gravity->phiP); + + dims[3] = 1; + std::snprintf(fieldName,NAMESIZE,"Pot"); + + WriteVariable(fileHdl, 4, dims, fieldName, locPot.data()); + } + // Write Flux /* nx1=this->np_tot[IDIR]; diff --git a/src/dataBlock/fargo.hpp b/src/dataBlock/fargo.hpp index 3c3a9401..099c9935 100644 --- a/src/dataBlock/fargo.hpp +++ b/src/dataBlock/fargo.hpp @@ -508,7 +508,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s-m ; ss < s ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -518,7 +518,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s ; ss < s-m ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -544,7 +544,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s-m ; ss < s ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -554,7 +554,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s ; ss < s-m ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -637,7 +637,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s-m ; ss < s ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -647,7 +647,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s ; ss < s-m ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -672,7 +672,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s-m ; ss < s ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } @@ -682,7 +682,7 @@ void Fargo::ShiftFluid(const real t, const real dt, Fluid* hydro) { for(int ss = s ; ss < s-m ; ss++) { int sc; if(haveDomainDecomposition) { - sc = ss; + sc = ss + maxShift; } else { sc = sbeg + modPositive(ss-sbeg,n); } diff --git a/src/fluid/constrainedTransport/calcRiemannEmf.hpp b/src/fluid/constrainedTransport/calcRiemannEmf.hpp index 3c9479c2..2779ffba 100644 --- a/src/fluid/constrainedTransport/calcRiemannEmf.hpp +++ b/src/fluid/constrainedTransport/calcRiemannEmf.hpp @@ -22,8 +22,6 @@ template void ConstrainedTransport::CalcRiemannAverage() { idfx::pushRegion("ConstrainedTransport::calcRiemannAverage"); -#if EMF_AVERAGE == UCT_HLLD || EMF_AVERAGE == UCT_HLL - // Corned EMFs IdefixArray3D ez = this->ez; #if DIMENSIONS == 3 @@ -428,7 +426,6 @@ void ConstrainedTransport::CalcRiemannAverage() { #endif } ); -#endif // EMF_AVERAGE idfx::popRegion(); } diff --git a/src/grid.cpp b/src/grid.cpp index e4c6b7b4..ace5eb0c 100644 --- a/src/grid.cpp +++ b/src/grid.cpp @@ -158,8 +158,16 @@ Grid::Grid(Input &input) { for(int i=0 ; i < 3 ; i++) period[i] = 0; - // Check if the dec option has been passed when number of procs > 1 + // Check that number of procs > 1 if(idfx::psize>1) { + int ngridtot=1; + for(int dir=0 ; dir < DIMENSIONS; dir++) { + ngridtot *= np_int[dir]; + } + // Check that the total grid dimension is effectively divisible by number of procs + if(ngridtot % idfx::psize) + IDEFIX_ERROR("Total grid size must be a multiple of the number of mpi process"); + // Check that dec option has been passed if(input.CheckEntry("CommandLine","dec") != DIMENSIONS) { // No command line decomposition, make auto-decomposition if possible // (only when nproc and dimensions are powers of 2, and in 1D) @@ -185,7 +193,7 @@ Grid::Grid(Input &input) { int ntot=1; for(int dir=0 ; dir < DIMENSIONS; dir++) { nproc[dir] = input.Get("CommandLine","dec",dir); - // Check that the dimension is effectively divisible by number of procs + // Check that the dimension is effectively divisible by number of procs along each direction if(np_int[dir] % nproc[dir]) IDEFIX_ERROR("Grid size must be a multiple of the domain decomposition"); // Count the total number of procs we'll need for the specified domain decomposition diff --git a/src/gridHost.cpp b/src/gridHost.cpp index d31a232c..ad61f043 100644 --- a/src/gridHost.cpp +++ b/src/gridHost.cpp @@ -201,15 +201,25 @@ void GridHost::MakeGrid(Input &input) { // dir >= DIMENSIONS/ Init simple uniform grid for(int i = 0 ; i < np_tot[dir] ; i++) { // Initialize to default values - xstart[dir] = -0.5; - xend[dir] = 0.5; - + #if GEOMETRY != SPHERICAL + xstart[dir] = -0.5; + xend[dir] = 0.5; + #elif GEOMETRY == SPHERICAL + if(dir == JDIR) { + xstart[dir] = 0; + xend[dir] = M_PI; + } + if(dir == KDIR) { + xstart[dir] = -0.5; + xend[dir] = 0.5; + } + #endif this->xbeg[dir] = xstart[dir]; this->xend[dir] = xend[dir]; - dx[dir](i) = 1.0; - x[dir](i)=0.0; - xl[dir](i)=-0.5; - xr[dir](i)=0.5; + dx[dir](i) = xend[dir] - xstart[dir]; + x[dir](i)=0.5*(xstart[dir]+xend[dir]); + xl[dir](i)=xstart[dir]; + xr[dir](i)=xend[dir]; } } } diff --git a/src/gridHost.hpp b/src/gridHost.hpp index 065c65af..55283a5c 100644 --- a/src/gridHost.hpp +++ b/src/gridHost.hpp @@ -22,10 +22,10 @@ class GridHost { public: - std::array::HostMirror,3> x; ///< geometrical central points - std::array::HostMirror,3> xr; ///< cell right interface - std::array::HostMirror,3> xl; ///< cell left interface - std::array::HostMirror,3> dx; ///< cell width + std::array,3> x; ///< geometrical central points + std::array,3> xr; ///< cell right interface + std::array,3> xl; ///< cell left interface + std::array,3> dx; ///< cell width std::array xbeg; ///< Beginning of grid std::array xend; ///< End of grid diff --git a/src/kokkos b/src/kokkos index 15dc143e..09e775bf 160000 --- a/src/kokkos +++ b/src/kokkos @@ -1 +1 @@ -Subproject commit 15dc143e5f39949eece972a798e175c4b463d4b8 +Subproject commit 09e775bfc585840bb9ab1156cbd8d7d1c9e0cc6d diff --git a/src/main.cpp b/src/main.cpp index c242ccef..6b65c4f1 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -9,7 +9,7 @@ //@HEADER // ************************************************************************ // -// IDEFIX v 2.1.00 +// IDEFIX v 2.2.00 // // ************************************************************************ //@HEADER @@ -91,8 +91,12 @@ int main( int argc, char* argv[] ) { // instantiate required objects. DataBlock data(grid, input); TimeIntegrator Tint(input,data); + #ifdef WITH_PYTHON + Pydefix pydefix(input); + #endif Output output(input, data); Setup mysetup(input, grid, data, output); + idfx::cout << "Main: initialisation finished." << std::endl; char host[1024]; @@ -111,6 +115,9 @@ int main( int argc, char* argv[] ) { grid.ShowConfig(); data.ShowConfig(); Tint.ShowConfig(); + #ifdef WITH_PYTHON + pydefix.ShowConfig(); + #endif /////////////////////////////// // Initial conditions (or restart) @@ -118,10 +125,22 @@ int main( int argc, char* argv[] ) { // Are we restarting? if(input.restartRequested) { if(input.forceInitRequested) { - idfx::pushRegion("Setup::Initflow"); - mysetup.InitFlow(data); - data.DeriveVectorPotential(); - idfx::popRegion(); + #ifdef WITH_PYTHON + if(pydefix.haveInitflow) { + idfx::pushRegion("Pydefix::Initflow"); + pydefix.InitFlow(data); + } else { + idfx::pushRegion("Setup::Initflow"); + mysetup.InitFlow(data); + } + data.DeriveVectorPotential(); + idfx::popRegion(); + #else + idfx::pushRegion("Setup::Initflow"); + mysetup.InitFlow(data); + data.DeriveVectorPotential(); + idfx::popRegion(); + #endif } idfx::cout << "Main: Restarting from dump file." << std::endl; bool restartSuccess = output.RestartFromDump(data,input.restartFileNumber); @@ -134,8 +153,18 @@ int main( int argc, char* argv[] ) { } if(!input.restartRequested) { idfx::cout << "Main: Creating initial conditions." << std::endl; - idfx::pushRegion("Setup::Initflow"); - mysetup.InitFlow(data); + #ifdef WITH_PYTHON + if(pydefix.haveInitflow) { + idfx::pushRegion("Pydefix::Initflow"); + pydefix.InitFlow(data); + } else { + idfx::pushRegion("Setup::Initflow"); + mysetup.InitFlow(data); + } + #else + idfx::pushRegion("Setup::Initflow"); + mysetup.InitFlow(data); + #endif idfx::popRegion(); data.DeriveVectorPotential(); // This does something only when evolveVectorPotential is on data.SetBoundaries(); diff --git a/src/mpi.cpp b/src/mpi.cpp index d25056fa..4de7739c 100644 --- a/src/mpi.cpp +++ b/src/mpi.cpp @@ -860,54 +860,81 @@ void Mpi::CheckConfig() { IdefixArray1D src("MPIChecksrc",1); IdefixArray1D::HostMirror srcHost = Kokkos::create_mirror_view(src); - srcHost(0) = idfx::prank; - Kokkos::deep_copy(src, srcHost); - - IdefixArray1D dst("MPICheckdst",1); - IdefixArray1D::HostMirror dstHost = Kokkos::create_mirror_view(dst); - - MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_RETURN); - - // Capture segfaults - struct sigaction newHandler; - struct sigaction oldHandler; - memset(&newHandler, 0, sizeof(newHandler)); - newHandler.sa_flags = SA_SIGINFO; - newHandler.sa_sigaction = Mpi::SigErrorHandler; - sigaction(SIGSEGV, &newHandler, &oldHandler); - - try { - int ierr = MPI_Allreduce(src.data(), dst.data(), 1, MPI_INT64_T, MPI_SUM, MPI_COMM_WORLD); - if(ierr != 0) { - char MPImsg[MPI_MAX_ERROR_STRING]; - int MPImsgLen; - MPI_Error_string(ierr, MPImsg, &MPImsgLen); - throw std::runtime_error(std::string(MPImsg, MPImsgLen)); + if(idfx::prank == 0) { + srcHost(0) = 0; + Kokkos::deep_copy(src, srcHost); + } + + if(idfx::psize > 1) { + MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_RETURN); + + // Capture segfaults + struct sigaction newHandler; + struct sigaction oldHandler; + memset(&newHandler, 0, sizeof(newHandler)); + newHandler.sa_flags = SA_SIGINFO; + newHandler.sa_sigaction = Mpi::SigErrorHandler; + sigaction(SIGSEGV, &newHandler, &oldHandler); + try { + // We next circulate the info round-robin accross all the nodes to check that + // MPI can exchange buffers in idefix arrays + + MPI_Status status; + int ierrSend, ierrRecv; + if(idfx::prank == 0) { + ierrSend = MPI_Send(src.data(), 1, MPI_INT64_T, idfx::prank+1, 1, MPI_COMM_WORLD); + ierrRecv = MPI_Recv(src.data(), 1, MPI_INT64_T, idfx::psize-1, 1, MPI_COMM_WORLD, &status); + } else { + ierrRecv = MPI_Recv(src.data(), 1, MPI_INT64_T, idfx::prank-1, 1, MPI_COMM_WORLD, &status); + // Add our own rank to the data + Kokkos::deep_copy(srcHost, src); + srcHost(0) += idfx::prank; + Kokkos::deep_copy(src, srcHost); + ierrSend = MPI_Send(src.data(), 1, MPI_INT64_T, (idfx::prank+1)%idfx::psize, 1, + MPI_COMM_WORLD); + } + + if(ierrSend != 0) { + char MPImsg[MPI_MAX_ERROR_STRING]; + int MPImsgLen; + MPI_Error_string(ierrSend, MPImsg, &MPImsgLen); + throw std::runtime_error(std::string(MPImsg, MPImsgLen)); + } + if(ierrRecv != 0) { + char MPImsg[MPI_MAX_ERROR_STRING]; + int MPImsgLen; + MPI_Error_string(ierrSend, MPImsg, &MPImsgLen); + throw std::runtime_error(std::string(MPImsg, MPImsgLen)); + } + } catch(std::exception &e) { + std::stringstream errmsg; + errmsg << "Your MPI library is unable to perform Send/Recv on Idefix arrays."; + errmsg << std::endl; + #ifdef KOKKOS_ENABLE_CUDA + errmsg << "Check that your MPI library is CUDA aware." << std::endl; + #elif defined(KOKKOS_ENABLE_HIP) + errmsg << "Check that your MPI library is RocM aware." << std::endl; + #else + errmsg << "Check your MPI library configuration." << std::endl; + #endif + errmsg << "Error: " << e.what() << std::endl; + IDEFIX_ERROR(errmsg); } - } catch(std::exception &e) { - std::stringstream errmsg; - errmsg << "Your MPI library is unable to perform reductions on Idefix arrays."; - errmsg << std::endl; - #ifdef KOKKOS_ENABLE_CUDA - errmsg << "Check that your MPI library is CUDA aware." << std::endl; - #elif defined(KOKKOS_ENABLE_HIP) - errmsg << "Check that your MPI library is RocM aware." << std::endl; - #else - errmsg << "Check your MPI library configuration." << std::endl; - #endif - errmsg << "Error: " << e.what() << std::endl; - IDEFIX_ERROR(errmsg); + // Restore old handlers + sigaction(SIGSEGV, &oldHandler, NULL ); + MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_ARE_FATAL); } - // Restore old handlers - sigaction(SIGSEGV, &oldHandler, NULL ); - MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_ARE_FATAL); - // Check that we have the proper end result - Kokkos::deep_copy(dstHost, dst); + // Check that we have the proper end result + Kokkos::deep_copy(srcHost, src); int64_t size = static_cast(idfx::psize); - if(dstHost(0) != size*(size-1)/2) { + int64_t rank = static_cast(idfx::prank); + int64_t result = rank == 0 ? size*(size-1)/2 : rank*(rank+1)/2; + + if(srcHost(0) != result) { + idfx::cout << "got " << srcHost(0) << " expected " << result << std::endl; std::stringstream errmsg; - errmsg << "Your MPI library managed to perform reduction on Idefix Arrays, but the result "; + errmsg << "Your MPI library managed to perform MPI exchanges on Idefix Arrays, but the result "; errmsg << "is incorrect. " << std::endl; errmsg << "Check your MPI library configuration." << std::endl; IDEFIX_ERROR(errmsg); diff --git a/src/output/output.cpp b/src/output/output.cpp index 658e1cdd..a67cb166 100644 --- a/src/output/output.cpp +++ b/src/output/output.cpp @@ -10,8 +10,17 @@ #include "dataBlock.hpp" #include "fluid.hpp" #include "slice.hpp" +#include "dataBlockHost.hpp" -Output::Output(Input &input, DataBlock &data) { +#ifdef WITH_PYTHON +#include "pydefix.hpp" +#endif + +Output::Output(Input &input, DataBlock &data) +#ifdef WITH_PYTHON +:pydefix(input) +#endif +{ idfx::pushRegion("Output::Output"); // initialise the output objects for each format @@ -116,6 +125,20 @@ Output::Output(Input &input, DataBlock &data) { } } + // Initialise python script outputs + if(input.CheckEntry("Output","python")>0) { + #ifndef WITH_PYTHON + IDEFIX_ERROR("Usage of python outputs requires Idefix to be compiled with Python capabilities"); + #else + pythonPeriod = input.Get("Output","python",0); + if(pythonPeriod>=0.0) { // backward compatibility (negative value means no file) + pythonLast = data.t - pythonPeriod; // write something in the next CheckForWrite() + pythonEnabled = true; + pythonNumber = 0; + } + #endif + } + // Register variables that are needed in restart dumps data.dump->RegisterVariable(&dumpLast, "dumpLast"); data.dump->RegisterVariable(&analysisLast, "analysisLast"); @@ -123,6 +146,10 @@ Output::Output(Input &input, DataBlock &data) { #ifdef WITH_HDF5 data.dump->RegisterVariable(&xdmfLast, "xdmfLast"); #endif + #ifdef WITH_PYTHON + data.dump->RegisterVariable(&pythonLast, "pythonLast"); + data.dump->RegisterVariable(&pythonNumber,"pythonNumber"); + #endif idfx::popRegion(); } @@ -228,6 +255,24 @@ int Output::CheckForWrites(DataBlock &data) { slices[i]->CheckForWrite(data); } } + + #ifdef WITH_PYTHON + if(pythonEnabled) { + if(data.t >= pythonLast + pythonPeriod) { + elapsedTime -= timer.seconds(); + pydefix.Output(data,pythonNumber); + pythonNumber++; + elapsedTime += timer.seconds(); + // Check if our next predicted output should already have happened + if((pythonLast+pythonPeriod <= data.t) && pythonPeriod>0.0) { + // Move forward analysisLast + while(pythonLast <= data.t - pythonPeriod) { + pythonLast += pythonPeriod; + } + } + } + } + #endif // Do we need a restart dump? if(dumpEnabled) { bool haveClockDump = false; diff --git a/src/output/output.hpp b/src/output/output.hpp index 54cb97af..039ab157 100644 --- a/src/output/output.hpp +++ b/src/output/output.hpp @@ -18,6 +18,9 @@ #ifdef WITH_HDF5 #include "xdmf.hpp" #endif +#ifdef WITH_PYTHON +#include "pydefix.hpp" +#endif #include "dump.hpp" #include "slice.hpp" @@ -79,6 +82,14 @@ class Output { bool haveSlices = false; std::vector> slices; + #ifdef WITH_PYTHON + Pydefix pydefix; + bool pythonEnabled = false; + real pythonPeriod; + real pythonLast; + real pythonNumber; + #endif + Kokkos::Timer timer; double elapsedTime{0.0}; }; diff --git a/src/output/vtk.cpp b/src/output/vtk.cpp index ae4a4c60..98e4c257 100644 --- a/src/output/vtk.cpp +++ b/src/output/vtk.cpp @@ -124,24 +124,51 @@ Vtk::Vtk(Input &input, DataBlock *datain, std::string filebase) { this->xnode = new float[nx1+ioffset]; this->ynode = new float[nx2+joffset]; this->znode = new float[nx3+koffset]; + this->xcenter = new float[nx1]; + this->ycenter = new float[nx2]; + this->zcenter = new float[nx3]; for (int32_t i = 0; i < nx1 + ioffset; i++) { - if(grid.np_tot[IDIR] == 1) // only one dimension in this direction + if(grid.np_tot[IDIR] == 1) { // only one dimension in this direction xnode[i] = bigEndian(static_cast(grid.x[IDIR](i))); - else + } else { xnode[i] = bigEndian(static_cast(grid.xl[IDIR](i + grid.nghost[IDIR]))); + } } for (int32_t j = 0; j < nx2 + joffset; j++) { - if(grid.np_tot[JDIR] == 1) // only one dimension in this direction + if(grid.np_tot[JDIR] == 1) { // only one dimension in this direction ynode[j] = bigEndian(static_cast(grid.x[JDIR](j))); - else + } else { ynode[j] = bigEndian(static_cast(grid.xl[JDIR](j + grid.nghost[JDIR]))); + } } for (int32_t k = 0; k < nx3 + koffset; k++) { - if(grid.np_tot[KDIR] == 1) + if(grid.np_tot[KDIR] == 1) { znode[k] = bigEndian(static_cast(grid.x[KDIR](k))); - else + } else { znode[k] = bigEndian(static_cast(grid.xl[KDIR](k + grid.nghost[KDIR]))); + } + } + for (int32_t i = 0; i < nx1; i++) { + if(grid.np_tot[IDIR] == 1) { // only one dimension in this direction + xcenter[i] = xnode[i]; + } else { + xcenter[i] = bigEndian(static_cast(grid.x[IDIR](i + grid.nghost[IDIR]))); + } + } + for (int32_t j = 0; j < nx2; j++) { + if(grid.np_tot[JDIR] == 1) { // only one dimension in this direction + ycenter[j] = ynode[j]; + } else { + ycenter[j] = bigEndian(static_cast(grid.x[JDIR](j + grid.nghost[JDIR]))); + } + } + for (int32_t k = 0; k < nx3; k++) { + if(grid.np_tot[KDIR] == 1) { + zcenter[k] = znode[k]; + } else { + zcenter[k] = bigEndian(static_cast(grid.x[KDIR](k + grid.nghost[KDIR]))); + } } #if VTK_FORMAT == VTK_STRUCTURED_GRID // VTK_FORMAT /* -- Allocate memory for node_coord which is later used -- */ @@ -375,8 +402,8 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) { #elif VTK_FORMAT == VTK_STRUCTURED_GRID ssheader << "DATASET STRUCTURED_GRID" << std::endl; #endif - // fields: geometry, periodicity, time - int nfields = 3; + // fields: geometry, periodicity, time, 6 NativeCoordinates (x1l, x2l, x3l, x1c, x2c, x3c) + int nfields = 9; // Write grid geometry in the VTK file ssheader << "FIELD FieldData " << nfields << std::endl; @@ -426,7 +453,71 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) { // Done, add cariage return for next ascii write ssheader << std::endl; + // write x1l native coordinates + ssheader << "X1L_NATIVE_COORDINATES 1 " << (nx1 + ioffset) << " float" << std::endl; + // Flush the ascii header + header = ssheader.str(); + WriteHeaderString(header.c_str(), fvtk); + // reset the string stream + ssheader.str(std::string()); + WriteHeaderBinary(this->xnode, nx1 + ioffset, fvtk); + // Done, add cariage return for next ascii write + ssheader << std::endl; + + // write x2l native coordinates + ssheader << "X2L_NATIVE_COORDINATES 1 " << (nx2 + joffset) << " float" << std::endl; + // Flush the ascii header + header = ssheader.str(); + WriteHeaderString(header.c_str(), fvtk); + // reset the string stream + ssheader.str(std::string()); + WriteHeaderBinary(this->ynode, nx2 + joffset, fvtk); + // Done, add cariage return for next ascii write + ssheader << std::endl; + + // write x3l native coordinates + ssheader << "X3L_NATIVE_COORDINATES 1 " << (nx3 + koffset) << " float" << std::endl; + // Flush the ascii header + header = ssheader.str(); + WriteHeaderString(header.c_str(), fvtk); + // reset the string stream + ssheader.str(std::string()); + WriteHeaderBinary(this->znode, nx3 + koffset, fvtk); + // Done, add cariage return for next ascii write + ssheader << std::endl; + + // write x1 native coordinates + ssheader << "X1C_NATIVE_COORDINATES 1 " << nx1 << " float" << std::endl; + // Flush the ascii header + header = ssheader.str(); + WriteHeaderString(header.c_str(), fvtk); + // reset the string stream + ssheader.str(std::string()); + WriteHeaderBinary(this->xcenter, nx1, fvtk); + // Done, add cariage return for next ascii write + ssheader << std::endl; + + // write x2 native coordinates + ssheader << "X2C_NATIVE_COORDINATES 1 " << nx2 << " float" << std::endl; + // Flush the ascii header + header = ssheader.str(); + WriteHeaderString(header.c_str(), fvtk); + // reset the string stream + ssheader.str(std::string()); + WriteHeaderBinary(this->ycenter, nx2, fvtk); + // Done, add cariage return for next ascii write + ssheader << std::endl; + // write x3 native coordinates + ssheader << "X3C_NATIVE_COORDINATES 1 " << nx3 << " float" << std::endl; + // Flush the ascii header + header = ssheader.str(); + WriteHeaderString(header.c_str(), fvtk); + // reset the string stream + ssheader.str(std::string()); + WriteHeaderBinary(this->zcenter, nx3, fvtk); + // Done, add cariage return for next ascii write + ssheader << std::endl; ssheader << "DIMENSIONS " << nx1 + ioffset << " " << nx2 + joffset << " " << nx3 + koffset << std::endl; @@ -486,7 +577,6 @@ void Vtk::WriteHeader(IdfxFileHandler fvtk, real time) { #undef VTK_STRUCTURED_GRID #undef VTK_RECTILINEAR_GRID - /* ********************************************************************* */ void Vtk::WriteScalar(IdfxFileHandler fvtk, float* Vin, const std::string &var_name) { /*! diff --git a/src/output/vtk.hpp b/src/output/vtk.hpp index a9503a5a..1a55f67f 100644 --- a/src/output/vtk.hpp +++ b/src/output/vtk.hpp @@ -129,6 +129,7 @@ class Vtk : public BaseVtk { // Coordinates needed by VTK outputs float *xnode, *ynode, *znode; + float *xcenter, *ycenter, *zcenter; IdefixHostArray4D node_coord; diff --git a/src/pydefix.cpp b/src/pydefix.cpp new file mode 100644 index 00000000..f0d8bc78 --- /dev/null +++ b/src/pydefix.cpp @@ -0,0 +1,365 @@ +// *********************************************************************************** +// Idefix MHD astrophysical code +// Copyright(C) Geoffroy R. J. Lesur +// and other code contributors +// Licensed under CeCILL 2.1 License, see COPYING for more information +// *********************************************************************************** + +#define PYBIND11_DETAILED_ERROR_MESSAGES + +#include "pydefix.hpp" +#include // everything needed for embedding +#include // for numpy arrays +#include // For STL vectors and containers +#include +#include +#include "idefix.hpp" +#include "dataBlock.hpp" +#include "dataBlockHost.hpp" + + +namespace py = pybind11; + +int Pydefix::ninstance = 0; + + +namespace PydefixTools { +// Functions provided by Idefix in Pydefix for user convenience +py::array_t GatherIdefixArray(IdefixHostArray3D in, + DataBlockHost dataHost, + bool keepBoundaries = true, + bool broadcast = true) { + idfx::pushRegion("PydefixTools::GatherIdefixArray"); + Grid *grid = dataHost.data->mygrid; + IdefixHostArray3D out; + py::array_t pyOut; + if(broadcast || idfx::prank==0) { + if(keepBoundaries) { + // Create a python-managed array, with memory accessible from Kokkos + pyOut = py::array_t({grid->np_tot[KDIR], + grid->np_tot[JDIR], + grid->np_tot[IDIR]}); + out = Kokkos::View> + (reinterpret_cast(pyOut.request().ptr), + grid->np_tot[KDIR], + grid->np_tot[JDIR], + grid->np_tot[IDIR]); + + } else { + pyOut = py::array_t({grid->np_int[KDIR], + grid->np_int[JDIR], + grid->np_int[IDIR]}); + out = Kokkos::View> + (reinterpret_cast(pyOut.request().ptr), + grid->np_int[KDIR], + grid->np_int[JDIR], + grid->np_int[IDIR]); + } + } + if(idfx::prank == 0) { + for(int rank = 0 ; rank < idfx::psize ; rank++) { + // np_tot: total size of the incoming array + // np_int: size that should be copied into global + // beg: offset in the incoming array where copy should begin + // gbeg: offset in the global array where copy should be begin + std::array np_int,np_tot, beg, gbeg; + IdefixHostArray3D buf; + + if(rank==0) { + np_int = dataHost.np_int; + np_tot = dataHost.np_tot; + gbeg = dataHost.gbeg; + beg = dataHost.beg; + + // Add back boundaries + if(keepBoundaries) { + for(int dir = 0 ; dir < DIMENSIONS ; dir++) { + if(dataHost.lbound[dir] != internal) { + np_int[dir] += dataHost.nghost[dir]; + gbeg[dir] -= dataHost.nghost[dir]; + beg[dir] -= dataHost.nghost[dir]; + } + if(dataHost.rbound[dir] != internal) { + np_int[dir] += dataHost.nghost[dir]; + } + } + } + buf = in; + } else { // target rank >0 + #ifdef WITH_MPI + MPI_Status status; + // Fetch the local array size + MPI_Recv(np_int.data(), 3, MPI_INT, rank, 010, MPI_COMM_WORLD, &status); + MPI_Recv(np_tot.data(), 3, MPI_INT, rank, 011, MPI_COMM_WORLD, &status); + MPI_Recv(beg.data(), 3, MPI_INT, rank, 012, MPI_COMM_WORLD, &status); + MPI_Recv(gbeg.data(), 3, MPI_INT, rank, 013, MPI_COMM_WORLD, &status); + + buf = IdefixHostArray3D("pydefix::tempArray", + np_tot[KDIR],np_tot[JDIR],np_tot[IDIR]); + // Fetch data + MPI_Recv(buf.data(), np_tot[IDIR]*np_tot[JDIR]*np_tot[KDIR], + realMPI, rank, 014, MPI_COMM_WORLD,&status); + #else + IDEFIX_ERROR("Can't deal with psize>1 without MPI."); + #endif + } // target rank + // Copy data from the buffer + + for(int k = 0 ; k < np_int[KDIR] ; k++) { + int kt = k+gbeg[KDIR]; + if(!keepBoundaries) kt -= dataHost.nghost[KDIR]; + for(int j = 0 ; j < np_int[JDIR] ; j++) { + int jt = j+gbeg[JDIR]; + if(!keepBoundaries) jt -= dataHost.nghost[JDIR]; + for(int i = 0 ; i < np_int[IDIR] ; i++) { + int it = i+gbeg[IDIR]; + if(!keepBoundaries) it -= dataHost.nghost[IDIR]; + out(kt,jt,it) = buf(k+beg[KDIR],j+beg[JDIR],i+beg[IDIR]); + } + } + }// End for + }// End loop on target rank for root process + } else { // MPI prank >0 + std::array np_int = dataHost.np_int; + std::array np_tot = dataHost.np_tot; + std::array gbeg = dataHost.gbeg; + std::array beg = dataHost.beg; + + // Add back boundaries + if(keepBoundaries) { + for(int dir = 0 ; dir < DIMENSIONS ; dir++) { + if(dataHost.lbound[dir] != internal) { + np_int[dir] += dataHost.nghost[dir]; + gbeg[dir] -= dataHost.nghost[dir]; + beg[dir] -= dataHost.nghost[dir]; + } + if(dataHost.rbound[dir] != internal) { + np_int[dir] += dataHost.nghost[dir]; + } + } + } + #ifdef WITH_MPI + // send the local array size + MPI_Send(np_int.data(), 3, MPI_INT, 0, 010, MPI_COMM_WORLD); + MPI_Send(np_tot.data(), 3, MPI_INT, 0, 011, MPI_COMM_WORLD); + MPI_Send(beg.data(), 3, MPI_INT, 0, 012, MPI_COMM_WORLD); + MPI_Send(gbeg.data(), 3, MPI_INT, 0, 013, MPI_COMM_WORLD); + MPI_Send(in.data(), np_tot[IDIR]*np_tot[JDIR]*np_tot[KDIR], realMPI, 0, 014, MPI_COMM_WORLD); + #else + IDEFIX_ERROR("Can't deal with psize>1 without MPI."); + #endif + } + // All is transfered + #ifdef WITH_MPI + if(broadcast) { + MPI_Bcast(out.data(), out.extent(0)*out.extent(1)*out.extent(2), realMPI, 0, MPI_COMM_WORLD); + } + #endif + + idfx::popRegion(); + return pyOut; +} +}// namespace PydefixTools + +/************************************ + * DataBlockHost Python binding + * **********************************/ + +PYBIND11_EMBEDDED_MODULE(pydefix, m) { + py::class_(m, "DataBlockHost") + .def(py::init<>()) + .def_readwrite("x", &DataBlockHost::x) + .def_readwrite("xr", &DataBlockHost::xr) + .def_readwrite("xl", &DataBlockHost::xl) + .def_readwrite("dx", &DataBlockHost::dx) + .def_readwrite("dV", &DataBlockHost::dV) + .def_readwrite("A", &DataBlockHost::A) + .def_readwrite("Vc", &DataBlockHost::Vc) + .def_readwrite("dustVc", &DataBlockHost::dustVc) + #if MHD == YES + .def_readwrite("Vs", &DataBlockHost::Vs) + .def_readwrite("Ve", &DataBlockHost::Ve) + .def_readwrite("J", &DataBlockHost::J) + .def_readwrite("Ex1", &DataBlockHost::Ex1) + .def_readwrite("Ex2", &DataBlockHost::Ex2) + .def_readwrite("Ex3", &DataBlockHost::Ex3) + #endif + .def_readwrite("xbeg", &DataBlockHost::xbeg) + .def_readwrite("xend", &DataBlockHost::xend) + .def_readwrite("gbeg", &DataBlockHost::gbeg) + .def_readwrite("gend", &DataBlockHost::gend) + .def_readwrite("beg", &DataBlockHost::beg) + .def_readwrite("end", &DataBlockHost::end) + .def_readwrite("np_tot", &DataBlockHost::np_tot) + .def_readwrite("np_int", &DataBlockHost::np_int) + .def_readwrite("nghost", &DataBlockHost::nghost) + .def_readwrite("InvDt", &DataBlockHost::InvDt) + .def_readwrite("t",&DataBlockHost::t) + .def_readwrite("dt",&DataBlockHost::dt); + + py::class_(m, "GridHost") + .def(py::init<>()) + .def_readwrite("x", &GridHost::x) + .def_readwrite("xr", &GridHost::xr) + .def_readwrite("xl", &GridHost::xl) + .def_readwrite("dx", &GridHost::dx) + .def_readwrite("xbeg", &GridHost::xbeg) + .def_readwrite("xend", &GridHost::xend) + .def_readwrite("np_tot", &GridHost::np_tot) + .def_readwrite("np_int", &GridHost::np_int) + .def_readwrite("nghost", &GridHost::nghost); + + m.attr("RHO") = RHO; + m.attr("VX1") = VX1; + m.attr("VX2") = VX2; + m.attr("VX3") = VX3; + m.attr("PRS") = PRS; + #if MHD == YES + m.attr("BX1") = BX1; + m.attr("BX2") = BX2; + m.attr("BX3") = BX3; + D_EXPAND( + m.attr("BX1s") = BX1s; , + m.attr("BX2s") = BX2s; , + m.attr("BX3s") = BX3s; ) + #endif + m.attr("IDIR") = IDIR; + m.attr("JDIR") = JDIR; + m.attr("KDIR") = KDIR; + + m.attr("prank") = idfx::prank; + m.attr("psize") = idfx::psize; + + m.def("GatherIdefixArray",&PydefixTools::GatherIdefixArray, + py::arg("in"), + py::arg("data"), + py::arg("keepBoundaries") = true, + py::arg("broadcast") = true, + "Gather arrays from MPI domain decomposition"); +} + + +template +void Pydefix::CallScript(std::string scriptName, std::string funcName, Ts... args) { + idfx::pushRegion("Pydefix::CallScript"); + try { + py::module_ script = py::module_::import(scriptName.c_str()); + py::object result = script.attr(funcName.c_str())(&args...); + } catch(std::exception &e) { + std::stringstream message; + message << "An exception occured while calling the Python interpreter" << std::endl + << "in file \"" << scriptName << ".py\" function \"" << funcName << "\":" + << std::endl + << e.what() << std::endl; + IDEFIX_ERROR(message); + } + idfx::popRegion(); +} + + +Pydefix::Pydefix(Input &input) { + // Check that the input has a [python] block + if(input.CheckBlock("Python")) { + this->isActive = true; + ninstance++; + // Check whether we need to start an interpreter + if(ninstance==1) { + idfx::cout << "Pydefix: start Python interpreter." << std::endl; + + py::initialize_interpreter(); + } + this->scriptFilename = input.Get("Python","script",0); + if(scriptFilename.substr(scriptFilename.length() - 3, 3).compare(".py")==0) { + IDEFIX_ERROR("You should not include the python script .py extension in your input file"); + } + if(input.CheckEntry("Python","output_function")>0) { + this->haveOutput = true; + this->outputFunctionName = input.Get("Python","output_function",0); + } + if(input.CheckEntry("Python","initflow_function")>0) { + this->haveInitflow = true; + this->initflowFunctionName = input.Get("Python","initflow_function",0); + } + } +} + +void Pydefix::Output(DataBlock &data, int n) { + idfx::pushRegion("Pydefix::Output"); + if(!this->isActive) { + IDEFIX_ERROR("Python Outputs requires the [python] block to be defined in your input file."); + } + if(!this->haveOutput) { + IDEFIX_ERROR("No python output function has been defined " + "in your input file [python]:output_function"); + } + DataBlockHost dataHost(data); + GridHost gridHost(*data.mygrid); + gridHost.SyncFromDevice(); + dataHost.SyncFromDevice(); + this->CallScript(this->scriptFilename,this->outputFunctionName,dataHost, gridHost, n); + idfx::popRegion(); +} + +void Pydefix::InitFlow(DataBlock &data) { + idfx::pushRegion("Pydefix::InitFlow"); + if(!this->isActive) { + IDEFIX_ERROR("Python Initflow requires the [python] block to be defined in your input file."); + } + if(!this->haveInitflow) { + IDEFIX_ERROR("No python initflow function has been defined " + "in your input file [python]:initflow_function"); + } + DataBlockHost dataHost(data); + dataHost.SyncFromDevice(); + this->CallScript(this->scriptFilename,this->initflowFunctionName,dataHost); + dataHost.SyncToDevice(); + idfx::popRegion(); +} + +void Pydefix::ShowConfig() { + if(isActive == false) { + idfx::cout << "Pydefix: DISABLED." << std::endl; + } else { + idfx::cout << "Pydefix: ENABLED." << std::endl; + if(haveOutput) { + idfx::cout << "Pydefix: output function ENABLED." << std::endl; + } else { + idfx::cout << "Pydefix: output function DISABLED." << std::endl; + } + if(haveInitflow) { + idfx::cout << "Pydefix: initflow function ENABLED." << std::endl; + } else { + idfx::cout << "Pydefix: initflow function DISABLED." << std::endl; + } + } +} + +Pydefix::~Pydefix() { + if(isActive) { + if(ninstance == 1) { + py::finalize_interpreter(); + idfx::cout << "Pydefix: shutdown Python interpreter." << std::endl; + } + ninstance--; + isActive = false; + } +} + + +/* +py::array_t Pydefix::toNumpyArray(const IdefixHostArray3D& in) { + py::array_t array({in.extent(0),in.extent(1),in.extent(2)},in.data()); + return array; +} + +py::array_t Pydefix::toNumpyArray(const IdefixHostArray4D& in) { + py::array_t array({in.extent(0),in.extent(1),in.extent(2),in.extent(3)},in.data()); + return array; +} +*/ diff --git a/src/pydefix.hpp b/src/pydefix.hpp new file mode 100644 index 00000000..a12086f8 --- /dev/null +++ b/src/pydefix.hpp @@ -0,0 +1,227 @@ +// *********************************************************************************** +// Idefix MHD astrophysical code +// Copyright(C) Geoffroy R. J. Lesur +// and other code contributors +// Licensed under CeCILL 2.1 License, see COPYING for more information +// *********************************************************************************** + +#ifndef PYDEFIX_HPP_ +#define PYDEFIX_HPP_ + + +#define PYBIND11_DETAILED_ERROR_MESSAGES +#include +#include +#include +#include +#include +#include "idefix.hpp" +#include "input.hpp" + +namespace py = pybind11; + +class DataBlock; +class DataBlockHost; + +class Pydefix { + public: + explicit Pydefix(Input&); + ~Pydefix(); + void Output(DataBlock &, int); + void InitFlow(DataBlock &); + void ShowConfig(); + bool isActive{false}; + bool haveOutput{false}; + bool haveInitflow{false}; + private: + template + void CallScript(std::string, std::string, Ts...); + static int ninstance; + std::string scriptFilename; + std::string outputFunctionName; + std::string initflowFunctionName; +}; + + +namespace pybind11 { namespace detail { +// Caster for IdefixArray4D +template struct type_caster> { + public: + PYBIND11_TYPE_CASTER(IdefixHostArray4D, _("IdefixHostArray4D")); + + // Conversion part 1 (Python -> C++) + bool load(py::handle src, bool convert) { + if ( !convert && !py::array_t::check_(src) ) + return false; + + auto array = py::array_t::ensure(src); + if ( !array ) + return false; + + auto dims = array.ndim(); + if ( dims != 4 ) + return false; + + auto buf = array.request(); + + + + value = Kokkos::View> (reinterpret_cast(buf.ptr), + array.shape()[0], + array.shape()[1], + array.shape()[2], + array.shape()[3]); + + + return true; + } + + //Conversion part 2 (C++ -> Python) + static py::handle cast(const IdefixHostArray4D& src, + py::return_value_policy policy, + py::handle parent) { + py::none dummyDataOwner; + py::array_t a({src.extent(0), + src.extent(1), + src.extent(2), + src.extent(3)}, + src.data(), dummyDataOwner); + + return a.release(); + } +}; +// Caster for IdefixArray3D +template struct type_caster> { + public: + PYBIND11_TYPE_CASTER(IdefixHostArray3D, _("IdefixHostArray3D")); + + // Conversion part 1 (Python -> C++) + bool load(py::handle src, bool convert) { + idfx::pushRegion("Pydefix::TypeCaster3D Python->C"); + if ( !convert && !py::array_t::check_(src) ) + return false; + + auto array = py::array_t::ensure(src); + if ( !array ) + return false; + + auto dims = array.ndim(); + if ( dims != 3 ) + return false; + + auto buf = array.request(); + + value = Kokkos::View> (reinterpret_cast(buf.ptr), + array.shape()[0], + array.shape()[1], + array.shape()[2]); + + idfx::popRegion(); + return true; + } + + //Conversion part 2 (C++ -> Python) + static py::handle cast(const IdefixHostArray3D& src, + py::return_value_policy policy, + py::handle parent) { + idfx::pushRegion("Pydefix::TypeCaster3D C->Python"); + // Keep a local reference + auto arr = src; + py::none dummyDataOwner; + py::array_t a({src.extent(0), + src.extent(1), + src.extent(2)}, + src.data(),dummyDataOwner); + idfx::popRegion(); + return a.release(); + } +}; +// Caster for IdefixArray2D +template struct type_caster> { + public: + PYBIND11_TYPE_CASTER(IdefixHostArray2D, _("IdefixHostArray2D")); + + // Conversion part 1 (Python -> C++) + bool load(py::handle src, bool convert) { + if ( !convert && !py::array_t::check_(src) ) + return false; + + auto array = py::array_t::ensure(src); + if ( !array ) + return false; + + auto dims = array.ndim(); + if ( dims != 2 ) + return false; + + auto buf = array.request(); + + value = Kokkos::View> (reinterpret_cast(buf.ptr), + array.shape()[0], + array.shape()[1]); + + idfx::popRegion(); + return true; + } + + //Conversion part 2 (C++ -> Python) + static py::handle cast(const IdefixHostArray2D& src, + py::return_value_policy policy, + py::handle parent) { + py::none dummyOwner; + py::array_t a({src.extent(0),src.extent(1)},src.data(),dummyOwner); + return a.release(); + } +}; +// Caster for IdefixArray1D +template struct type_caster> { + public: + PYBIND11_TYPE_CASTER(IdefixHostArray1D, _("IdefixHostArray1D")); + + // Conversion part 1 (Python -> C++) + bool load(py::handle src, bool convert) { + if ( !convert && !py::array_t::check_(src) ) + return false; + + auto array = py::array_t::ensure(src); + if ( !array ) + return false; + + auto dims = array.ndim(); + if ( dims != 1 ) + return false; + + auto buf = array.request(); + + value = Kokkos::View> (reinterpret_cast(buf.ptr), + array.shape()[0]); + + idfx::popRegion(); + return true; + } + + //Conversion part 2 (C++ -> Python) + static py::handle cast(const IdefixHostArray1D& src, + py::return_value_policy policy, + py::handle parent) { + py::none dummyDataOwner; + py::array_t a(src.extent(0),src.data(),dummyDataOwner); + return a.release(); + } +}; +} // namespace detail +} // namespace pybind11 + +#endif // PYDEFIX_HPP_ diff --git a/test/HD/sod/idefix_py.ini b/test/HD/sod/idefix_py.ini new file mode 100644 index 00000000..59c58f64 --- /dev/null +++ b/test/HD/sod/idefix_py.ini @@ -0,0 +1,32 @@ +[Grid] +X1-grid 1 0.0 500 u 1.0 +X2-grid 1 0.0 1 u 1.0 +X3-grid 1 0.0 1 u 1.0 + +[TimeIntegrator] +CFL 0.8 +tstop 0.2 +first_dt 1.e-4 +nstages 2 + +[Hydro] +solver roe +gamma 1.4 + +[Boundary] +X1-beg outflow +X1-end outflow +X2-beg outflow +X2-end outflow +X3-beg outflow +X3-end outflow + +[Python] +script pydefix_example +output_function output +initflow_function initflow + +[Output] +vtk 0.1 +dmp 0.2 +python 0.01 diff --git a/test/HD/sod/pydefix_example.py b/test/HD/sod/pydefix_example.py new file mode 100644 index 00000000..58d73642 --- /dev/null +++ b/test/HD/sod/pydefix_example.py @@ -0,0 +1,19 @@ +import pydefix as pdfx +import numpy as np +import matplotlib.pyplot as plt + +def output(data): + + plt.figure() + plt.plot(data.x[pdfx.IDIR],data.Vc[pdfx.VX1,0,0,:],label='VX1') + plt.plot(data.x[pdfx.IDIR],data.Vc[pdfx.RHO,0,0,:],label='RHO') + plt.legend() + plt.show() + #data.Vc[0,0,0,10] = 2.0 + +def initflow(data): + + # Initialize the flow + data.Vc[pdfx.RHO,0,0,:] = 1.0 + data.Vc[pdfx.PRS,0,0,:] = 2.0 + data.Vc[pdfx.VX1,0,0,:] = np.sin(2.0*np.pi*data.x[0][:]) diff --git a/test/IO/pydefix/CMakeLists.txt b/test/IO/pydefix/CMakeLists.txt new file mode 100644 index 00000000..83d16a52 --- /dev/null +++ b/test/IO/pydefix/CMakeLists.txt @@ -0,0 +1,2 @@ +enable_idefix_property(Idefix_MHD) +enable_idefix_property(Idefix_PYTHON) diff --git a/test/IO/pydefix/README.md b/test/IO/pydefix/README.md new file mode 100644 index 00000000..820a0958 --- /dev/null +++ b/test/IO/pydefix/README.md @@ -0,0 +1,44 @@ +This directory shows how to use Idefix with pure python initial conditions and outputs. It reproduces the 2D OrszagTang vortex available in MHD/OrszagTang without requiring any single line of C++ code from the user. + +The python script `pydefix_example.py` initializes the initial condition of the OT test (`initflow`) and produces a series of PNG files through matplotlib (`output`). + +# Python modules installation + +In order to use pydefix, you need to be working in a python environement that includes pybind11. The simplest way is to install the suggested list of packages (that you can deploy in dedicated venv environement) + +```bash +pip install -r python_requirements.txt +``` + +# Configuration + +In order to use Pydefix, you need to switch the `Idefix_PYTHON` to ON in cmake. In this particular setup, it is automatically done for you in CMakeLists.txt to avoid mistakes. + +# Running + +Just run idefix as usual. + +# Troubleshooting + +It during configuration stage, you get: + +` +CMake Error at CMakeLists.txt:122 (find_package): + By not providing "Findpybind11.cmake" in CMAKE_MODULE_PATH this project has + asked CMake to find a package configuration file provided by "pybind11", + but CMake did not find one.` + +It means that cmake cannot find the location of pybind11 (this typically happens on MacOs). In order to locate pybind11, open a python interpreter, and get pybind11 install dir through: + +```python +import pybind11 +print(pybind11.__file__) +``` + +You can then exit the interpreter and set the pybind11_DIR environement variable to the right path: + +```bash +export pybind11_DIR=env/lib/python3.10/site-packages/pybind11 +``` + +you can then run cmake which should be able to find pybind11, and compile the code. diff --git a/test/IO/pydefix/definitions.hpp b/test/IO/pydefix/definitions.hpp new file mode 100644 index 00000000..69f0c49e --- /dev/null +++ b/test/IO/pydefix/definitions.hpp @@ -0,0 +1,4 @@ +#define COMPONENTS 2 +#define DIMENSIONS 2 + +#define GEOMETRY CARTESIAN diff --git a/test/IO/pydefix/idefix.ini b/test/IO/pydefix/idefix.ini new file mode 100644 index 00000000..57fdcfc4 --- /dev/null +++ b/test/IO/pydefix/idefix.ini @@ -0,0 +1,30 @@ +[Grid] +X1-grid 1 0.0 256 u 1.0 +X2-grid 1 0.0 256 u 1.0 +X3-grid 1 0.0 1 u 1.0 + +[TimeIntegrator] +CFL 0.6 +tstop 0.5 +first_dt 1.e-4 +nstages 2 + +[Hydro] +solver roe + +[Python] +script pydefix_example +output_function output +initflow_function initflow + +[Boundary] +X1-beg periodic +X1-end periodic +X2-beg periodic +X2-end periodic +X3-beg outflow +X3-end outflow + +[Output] +log 10 +python 0.02 diff --git a/test/IO/pydefix/pydefix_example.py b/test/IO/pydefix/pydefix_example.py new file mode 100644 index 00000000..a08e8646 --- /dev/null +++ b/test/IO/pydefix/pydefix_example.py @@ -0,0 +1,33 @@ +from pydefix import * +import numpy as np +import matplotlib.pyplot as plt + +# The output function +# the only argument is dataBlockHost python object, wrapping a dataBlockHost Idefix object +def output(data,grid,n): + # Note: if MPI is not enabled, GatherIdefixArray still works (but merely does a local copy) + pressure = GatherIdefixArray(data.Vc[PRS,:,:,:],data,broadcast=False,keepBoundaries=True) + # only process #0 performs the output + if prank==0: + plt.close() + plt.figure() + plt.pcolormesh(grid.x[IDIR],grid.x[JDIR],pressure[0,:,:],label='PRS',vmin=0.02,vmax=0.5,cmap='plasma') + plt.title("t=%.2f"%data.t) + plt.colorbar() + plt.savefig("PRS.%.4d.png"%n) + + +def initflow(data): + # Field amplitude + B0 = 1/np.sqrt(4*np.pi) + + [z,y,x] = np.meshgrid(data.x[KDIR], data.x[JDIR], data.x[IDIR], indexing='ij') + + # Initialize the flow + data.Vc[RHO,:,:,:] = 25/(36*np.pi) + data.Vc[PRS,:,:,:] = 5/(12*np.pi) + data.Vc[VX1,:,:,:] = -np.sin(2*np.pi*y) + data.Vc[VX2,:,:,:] = np.sin(2*np.pi*x) + + data.Vs[BX1s,:,:-1,:-1] = -B0*np.sin(2*np.pi*y) + data.Vs[BX2s,:,:-1,:-1] = B0*np.sin(4*np.pi*x) diff --git a/test/IO/pydefix/python_requirements.txt b/test/IO/pydefix/python_requirements.txt new file mode 100644 index 00000000..6afd3019 --- /dev/null +++ b/test/IO/pydefix/python_requirements.txt @@ -0,0 +1,6 @@ +# note: version requirements are indicative and tests # should be able to run with +# older versions of our dependencies, though it is not guaranteed. + +numpy>=1.16.6 +matplotlib>=2.2.5 +pybind11>=2.12.0 diff --git a/test/MHD/AmbipolarCshock/python/testidefix.py b/test/MHD/AmbipolarCshock/python/testidefix.py index ed622fdf..8c95dbb1 100755 --- a/test/MHD/AmbipolarCshock/python/testidefix.py +++ b/test/MHD/AmbipolarCshock/python/testidefix.py @@ -72,7 +72,7 @@ def f(x,D): Dth=np.zeros(x.shape) for i in range(x.size): r.set_initial_value(DSim[iref],x[iref]) - Dth[i]=r.integrate(x[i]) + Dth[i]=r.integrate(x[i]).item() #print("Dth[%d]=%g"%(i,Dth[i])) bTh=np.sqrt(b0**2+2*A**2*(Dth-1)*(1/Dth-1/M**2)) diff --git a/test/MHD/AmbipolarCshock3D/python/testidefix.py b/test/MHD/AmbipolarCshock3D/python/testidefix.py index 52f10b63..a7354fba 100755 --- a/test/MHD/AmbipolarCshock3D/python/testidefix.py +++ b/test/MHD/AmbipolarCshock3D/python/testidefix.py @@ -75,7 +75,7 @@ def f(x,D): Dth=np.zeros(x.shape) for i in range(x.size): r.set_initial_value(DSim[iref],x[iref]) - Dth[i]=r.integrate(x[i]) + Dth[i]=r.integrate(x[i]).item() #print("Dth[%d]=%g"%(i,Dth[i])) bTh=np.sqrt(b0**2+2*A**2*(Dth-1)*(1/Dth-1/M**2)) diff --git a/test/MHD/FargoMHDSpherical/idefix.ini b/test/MHD/FargoMHDSpherical/idefix.ini index 4165ebc7..06bfcdf7 100644 --- a/test/MHD/FargoMHDSpherical/idefix.ini +++ b/test/MHD/FargoMHDSpherical/idefix.ini @@ -1,7 +1,7 @@ [Grid] -X1-grid 1 1.0 16 l 2.0 -X2-grid 1 1.2707963267948965 64 u 1.8707963267948966 -X3-grid 1 0.0 64 u 6.283185307179586 +X1-grid 1 1.0 16 l 2.0 +X2-grid 1 1.2707963267948965 16 u 1.8707963267948966 +X3-grid 1 0.0 128 u 6.283185307179586 [TimeIntegrator] CFL 0.5 diff --git a/test/MHD/FargoMHDSpherical/testme.py b/test/MHD/FargoMHDSpherical/testme.py index a09d3eeb..7ee3100b 100755 --- a/test/MHD/FargoMHDSpherical/testme.py +++ b/test/MHD/FargoMHDSpherical/testme.py @@ -53,3 +53,7 @@ def testMe(test): test.vectPot=False test.mpi=True testMe(test) + + test.vectPot=True + test.mpi=True + testMe(test)