diff --git a/.github/ISSUE_TEMPLATE/1-bug.yml b/.github/ISSUE_TEMPLATE/1-bug.yml
new file mode 100644
index 0000000000..72a46551b4
--- /dev/null
+++ b/.github/ISSUE_TEMPLATE/1-bug.yml
@@ -0,0 +1,119 @@
+name: 🐞 Bug Report
+description: File a bug report
+labels: ["triage"]
+
+body:
+
+ - type: markdown
+ attributes:
+ value: |
+ Thanks for taking the time to fill out this bug report!
+
+ This issue tracker is for reporting bugs and issues found in the GAMER project.
+ If you are looking for help with GAMER, please check out the [GAMER Slack](https://join.slack.com/t/gamer-project/shared_invite/enQtNTUwMDA5ODAwMTMzLTc3ZWY2MWE2YTlmMDI0MTQ4M2JjOTg2NmU4OWVkOGY1ZTI3MmY5NjUxOTk1ZjM5ZjNjOGViMGY3ZGExMDdiYzU).
+
+ Please fill out the following form with as much detail as possible to help us diagnose and resolve the issue efficiently.
+
+ - type: textarea
+ id: issue_content
+ attributes:
+ label: 🔎 What happened?
+ description: Explain the problem and expected behavior.
+ validations:
+ required: true
+
+ - type: textarea
+ id: reproduce_step
+ attributes:
+ label: 📃 Steps to reproduce
+ description: Please provide detailed steps to help us reproduce the issue.
+ placeholder: |
+ e.g.,
+ 1. Copy the shocktube test problem from `example/Hydro/Riemann/*`.
+ 2. Build GAMER with CXXFLAG -g -O0 in the machine configuration file.
+ 3. Run `./gamer`, and a segmentation fault occurs.
+ validations:
+ required: true
+
+ - type: input
+ id: commit_hash
+ attributes:
+ label: ⌚ Commit hash
+ description: Please provide the commit hash of the version you are using by running `git rev-parse HEAD`. If you are not using the public `main` branch, also include the repository link and branch name.
+ placeholder: e.g., 2f1ceb7ceb6249f0252d58cdc0269383631bdd68 or 2f1ceb7
+ validations:
+ required: true
+
+ - type: input
+ id: config_cmd
+ attributes:
+ label: 🔧 Configuration command
+ description: Please provide the configuration command you used to generate the `Makefile`. Alternatively, copy and paste the contents of `Makefile.log`.
+ placeholder: e.g., python configure.py --fftw=FFTW2 --gravity=true --gpu=true
+
+ - type: textarea
+ id: files_modified
+ attributes:
+ label: 🔨 Source files modified
+ description: Please provide a list of source files you have modified, if any.
+ placeholder: e.g., src/Hydro/Gravity/Init_TestProb_Hydro_Gravity.cpp
+
+ - type: dropdown
+ id: operation_system
+ attributes:
+ label: 💻 Operating system
+ description: Please specify the OS of your machine. If your OS is not listed, please select "Other" and specify it in the "Additional information" section.
+ multiple: false
+ options:
+ - linux (x86)
+ - linux (ARM)
+ - macOS (Intel)
+ - macOS (Apple silicon)
+ - Windows (x86)
+ - Windows (ARM)
+ - Other (Please specify)
+
+ - type: textarea
+ id: machine
+ attributes:
+ label: 💾 Machine configuration file
+ description: Please provide the contents of the machine configuration file you used under `gamer/configs`. If you used a built-in configuration, simply specify its name (e.g., `spock_intel.config`).
+ placeholder: |
+ e.g.,
+ # NTU-spock
+ CUDA_PATH /software/cuda/12.1
+ FFTW2_PATH /software/fftw/2.1.5-intel-2023.1.0-openmpi-4.1.5-ucx_mt
+ FFTW3_PATH /software/fftw/3.3.10-intel-2023.1.0-openmpi-4.1.5-ucx_mt
+ MPI_PATH /software/openmpi/4.1.5-ucx_mt-intel-2023.1.0
+ ...
+
+ # compilers
+ CXX icpc
+ CXX_MPI mpicxx
+ ...
+
+ - type: checkboxes
+ id: labels
+ attributes:
+ label: 🔖 Related topics
+ description: Select relevant topics (you may choose more than one).
+ options:
+ - label: Hydro
+ - label: MHD
+ - label: FDM
+ - label: AMR
+ - label: Gravity
+ - label: Particle
+ - label: Parallelization
+ - label: GPU
+ - label: Memory
+ - label: YT
+ - label: Tool
+ - label: Docs
+ - label: Other
+
+ - type: textarea
+ id: additional_info
+ attributes:
+ label: 💬 Additional information
+ description: Please provide any additional information that may help diagnose the issue (e.g., screenshots, stdout/stderr, and `Record__Note`).
diff --git a/.github/ISSUE_TEMPLATE/2-feature.yml b/.github/ISSUE_TEMPLATE/2-feature.yml
new file mode 100644
index 0000000000..b56db8c472
--- /dev/null
+++ b/.github/ISSUE_TEMPLATE/2-feature.yml
@@ -0,0 +1,62 @@
+name: 💪 Feature Request
+description: Suggest a new feature
+labels: ["triage"]
+
+body:
+
+ - type: markdown
+ attributes:
+ value: |
+ Thanks for suggesting a new feature!
+
+ Please fill out the following form to help us better understand your request.
+
+ - type: textarea
+ id: feature_description
+ attributes:
+ label: 💭 Proposed feature
+ description: Provide a clear and concise description of the feature you are suggesting.
+ validations:
+ required: true
+
+ - type: textarea
+ id: use_case
+ attributes:
+ label: 🔆 Use case
+ description: Describe where this feature would be applied.
+ validations:
+ required: true
+
+ - type: textarea
+ id: motivation
+ attributes:
+ label: 💡 Motivation
+ description: Explain why this feature is needed and how it benefits the project.
+ validations:
+ required: true
+
+ - type: textarea
+ id: additional_info
+ attributes:
+ label: 💬 Additional information
+ description: Provide any additional information or context that may be helpful.
+
+ - type: checkboxes
+ id: labels
+ attributes:
+ label: 🔖 Related topics
+ description: Select relevant topics (you may choose more than one).
+ options:
+ - label: Hydro
+ - label: MHD
+ - label: FDM
+ - label: AMR
+ - label: Gravity
+ - label: Particle
+ - label: Parallelization
+ - label: GPU
+ - label: Memory
+ - label: YT
+ - label: Tool
+ - label: Docs
+ - label: Other
diff --git a/.github/ISSUE_TEMPLATE/3-assign-task.yml b/.github/ISSUE_TEMPLATE/3-assign-task.yml
new file mode 100644
index 0000000000..5cb039b190
--- /dev/null
+++ b/.github/ISSUE_TEMPLATE/3-assign-task.yml
@@ -0,0 +1,48 @@
+name: 🐟 Assign Task
+description: Assign a new set of tasks (Dev Team Only)
+
+body:
+
+ - type: markdown
+ attributes:
+ value: |
+ 🐟🐟🐟🐟🐟🐟🐟🐟🐟|| DEV TEAM ONLY ||🐟🐟🐟🐟🐟🐟🐟🐟🐟
+ This is a template for the development team to assign a new set of tasks.
+ If you are not a member of the team, please do not use this template.
+ 🐟🐟🐟🐟🐟🐟🐟🐟🐟|| DEV TEAM ONLY ||🐟🐟🐟🐟🐟🐟🐟🐟🐟
+
+ - type: textarea
+ id: goal
+ attributes:
+ label: 🎯 Goals
+ validations:
+ required: true
+
+ - type: textarea
+ id: task
+ attributes:
+ label: 📋 Mandatory tasks
+ value: |
+ - [ ] Task 1
+ - [ ] Task 2
+ - [ ] Task 3
+ validations:
+ required: true
+
+ - type: textarea
+ id: opt_task
+ attributes:
+ label: Optional tasks
+ value: |
+ - [ ] Task A
+ - [ ] Task B
+
+ - type: textarea
+ id: additional
+ attributes:
+ label: References
+
+ - type: textarea
+ id: note
+ attributes:
+ label: Supplementary info
diff --git a/.github/ISSUE_TEMPLATE/config.yml b/.github/ISSUE_TEMPLATE/config.yml
new file mode 100644
index 0000000000..7c074225b1
--- /dev/null
+++ b/.github/ISSUE_TEMPLATE/config.yml
@@ -0,0 +1,5 @@
+blank_issues_enabled: false
+contact_links:
+ - name: ❓ General help or Question
+ url: https://join.slack.com/t/gamer-project/shared_invite/enQtNTUwMDA5ODAwMTMzLTc3ZWY2MWE2YTlmMDI0MTQ4M2JjOTg2NmU4OWVkOGY1ZTI3MmY5NjUxOTk1ZjM5ZjNjOGViMGY3ZGExMDdiYzU
+ about: Join our Slack for discussion
diff --git a/configs/eureka_gnu.config b/configs/eureka_gnu.config
index 1896846d45..b8ba966717 100644
--- a/configs/eureka_gnu.config
+++ b/configs/eureka_gnu.config
@@ -8,6 +8,7 @@ GRACKLE_PATH
GSL_PATH /software/gsl/default
LIBYT_PATH
CUFFTDX_PATH /software/cuFFTDx/default
+HYPRE_PATH
# compilers
CXX g++
diff --git a/configs/eureka_intel.config b/configs/eureka_intel.config
index 2647d6be68..69265adec7 100644
--- a/configs/eureka_intel.config
+++ b/configs/eureka_intel.config
@@ -8,6 +8,7 @@ GRACKLE_PATH
GSL_PATH /software/gsl/default
LIBYT_PATH
CUFFTDX_PATH /software/cuFFTDx/default
+HYPRE_PATH
# compilers
CXX icpc
diff --git a/configs/spock_intel.config b/configs/spock_intel.config
index 83e807b775..832dc39303 100644
--- a/configs/spock_intel.config
+++ b/configs/spock_intel.config
@@ -8,6 +8,7 @@ GRACKLE_PATH
GSL_PATH /software/gsl/2.6-intel-2023.1.0
LIBYT_PATH
CUFFTDX_PATH /software/cuFFTDx/22.11
+HYPRE_PATH
# compilers
CXX icpc
diff --git a/configs/template.config b/configs/template.config
index 93bddd6578..3ff39365f7 100644
--- a/configs/template.config
+++ b/configs/template.config
@@ -10,6 +10,7 @@ GRACKLE_PATH /path/to/grackle
GSL_PATH /path/to/gsl
LIBYT_PATH /path/to/libyt
CUFFTDX_PATH /path/to/cufftdx
+HYPRE_PATH
# 2. Compiler type
CXX icpc # Serial compiler
diff --git a/doc/wiki/Developing-with-VS-Code.md b/doc/wiki/Developing-with-VS-Code.md
new file mode 100644
index 0000000000..499cedc0f7
--- /dev/null
+++ b/doc/wiki/Developing-with-VS-Code.md
@@ -0,0 +1,66 @@
+This guide provides step-by-step instructions for setting up and using Visual Studio Code (VS Code) to develop the GAMER codebase.
+
+## Setup
+
+### Prerequisites
+
+- **Quick Start**: Follow the [[Quick Start | Quick-Start]] guide to download GAMER and complete at least one [[demo | Quick-Start:-1D-Shock-Tube]].
+- **Visual Studio Code**: Download and install VS Code from [https://code.visualstudio.com/](https://code.visualstudio.com/).
+- **C/C++ Extension**: Install the "C/C++" extension from the VS Code Marketplace.
+
+### Setting Up the Workspace
+
+1. **Launch VS Code**.
+2. **Open the GAMER Project Folder**:
+ - Go to `File` > `Open Folder...`.
+ - Select your GAMER project directory.
+
+> [!TIP]
+> When using remote-SSH, open the directory as an absolute path to avoid [this issue](https://github.com/microsoft/vscode-cpptools/issues/4818).
+
+### Configuring VS Code for GAMER
+
+Run the following script from the root directory of the GAMER project:
+```bash
+sh tool/vscode/copy_to_vscode.sh
+```
+This script copies the necessary configuration files to the `.vscode` directory, integrating GAMER with VS Code.
+
+## Developing with VS Code
+
+Before running the tasks below, **set the working directory** by selecting `Terminal` > `Run Task...` > `set-working-bin` and entering the name of the working directory under `bin/` where the input files are located.
+
+### Configuring GAMER
+
+Select `Terminal` > `Run Task...` > `config-GAMER` to configure GAMER using the `generate_make.sh` script in your working directory.
+
+### Building GAMER
+
+After configuring GAMER with [configure.py](https://github.com/gamer-project/gamer/wiki/Installation%3A-Configure.py) or `generate_make.sh`, select `Terminal` > `Run Task...` > `build-GAMER` to start the build process. This updates the macros and enables IntelliSense highlighting.
+
+> [!TIP]
+> To configure and build GAMER in one step, select `Terminal` > `Run Build Task...` or press `Ctrl + Shift + B` to run `config-GAMER` and `build-GAMER` sequentially.
+
+### Debugging GAMER
+
+To start debugging, select `Run` > `Start Debugging` or press `F5`. After entering the working directory, the debugger will launch. See the [official documentation](https://code.visualstudio.com/docs/editor/debugging) to learn more about debugging in VS Code.
+
+> [!IMPORTANT]
+> Ensure the compiler flags in `Makefile` are set to `-g -O0` for debugging. (TBD: Add an argument to `configure.py` to set these flags.)
+
+> [!NOTE]
+> If `gdb` is not supported on macOS, you can set up `lldb` as the debugger. Ensure [`lldb-mi`](https://github.com/lldb-tools/lldb-mi) is installed, then select `Terminal` > `Run Task...` > `updated_mac_launch`. This updates the debugger path in `launch.json` to your `lldb-mi` installation.
+> For manual setup or additional details, refer to the [official documentation](https://code.visualstudio.com/docs/cpp/launch-json-reference).
+
+### Cleaning the Working Directory
+
+Select `Terminal` > `Run Task...` > `clean-work-dir` to clean the working directory using the `clean.sh` script in your working directory.
+
+## Understanding Configuration Files
+
+The following configuration files are copied to the `.vscode` directory:
+- `c_cpp_properties.json`: Configures IntelliSense settings, including include paths and macros. See the [schema reference](https://code.visualstudio.com/docs/cpp/c-cpp-properties-schema-reference) and [IntelliSense documentation](https://code.visualstudio.com/docs/editor/intellisense) for details.
+- `launch.json`: Defines debugging configurations such as executable paths and arguments. See the [official documentation](https://code.visualstudio.com/docs/cpp/launch-json-reference) for more information.
+- `settings.json`: Specifies editor settings, such as indentation spaces and file types for extensions. See the [VS Code settings guide](https://code.visualstudio.com/docs/editor/settings) for details.
+- `tasks.json`: Defines build and auxiliary tasks. Learn more about [tasks in VS Code](https://code.visualstudio.com/docs/editor/tasks).
+- `gamercpp.natvis`: Customizes data structure visualizations in the debugger. Learn more about [customizing native object views](https://learn.microsoft.com/en-us/visualstudio/debugger/create-custom-views-of-native-objects?view=vs-2022).
diff --git a/doc/wiki/Home.md b/doc/wiki/Home.md
index 3e7cece293..6916e16db9 100644
--- a/doc/wiki/Home.md
+++ b/doc/wiki/Home.md
@@ -1,6 +1,6 @@
**GAMER** is a **G**PU-accelerated **A**daptive **ME**sh **R**efinement
-code for astrophysics. It features extremely high performance and
-parallel scalability and supports a rich set of physics modules.
+code for astrophysics. It features high computational performance and
+strong parallel scalability, and supports a rich set of physics modules.
***
@@ -10,19 +10,22 @@ parallel scalability and supports a rich set of physics modules.
* [Magnetohydrodynamics](https://iopscience.iop.org/article/10.3847/1538-4365/aac49e/meta)
* [Special relativistic hydrodynamics](https://academic.oup.com/mnras/article/504/3/3298/6224873)
* Self-gravity and external gravity
-* Particles
+* Particles and tracers
* Chemistry and radiative processes with [GRACKLE](http://grackle.readthedocs.io/en/latest/index.html)
* General equation of state
+* Passively advected scalars
+* Star formation
* [Cosmic rays with anisotropic diffusion](https://iopscience.iop.org/article/10.3847/1538-4357/ad50c5#apjad50c5app2)
-* Fuzzy (Wave) dark matter: [Nature Physics paper](http://www.nature.com/nphys/journal/v10/n7/covers/index.html), [code paper](https://arxiv.org/abs/2411.17288)
+* Fuzzy dark matter (FDM): [Nature Physics paper](http://www.nature.com/nphys/journal/v10/n7/covers/index.html), [code paper](https://arxiv.org/abs/2411.17288)
### Other Features
* Adaptive mesh refinement
* Adaptive timestep
-* Hybrid MPI/OpenMP/GPU parallelization (also support a CPU-only mode)
+* Hybrid MPI/OpenMP/GPU parallelization (also supports a CPU-only mode)
* Load balancing with Hilbert space-filling curve
* Bitwise reproducibility
-* [HDF5](https://support.hdfgroup.org/HDF5) output
+* Simple compilation using `configure.py` and machine files
+* [HDF5](https://www.hdfgroup.org/solutions/hdf5/) output
* Data analysis with [yt](http://yt-project.org)
* In situ Python analysis with [libyt](https://github.com/yt-project/libyt)
* Source-term interface
@@ -37,6 +40,7 @@ parallel scalability and supports a rich set of physics modules.
* Mailing list: [GAMER Google Group](https://groups.google.com/forum/#!forum/gamer-amr)
* Live chat: [GAMER Slack](https://join.slack.com/t/gamer-project/shared_invite/enQtNTUwMDA5ODAwMTMzLTc3ZWY2MWE2YTlmMDI0MTQ4M2JjOTg2NmU4OWVkOGY1ZTI3MmY5NjUxOTk1ZjM5ZjNjOGViMGY3ZGExMDdiYzU)
* Code papers:
+[GAMER-1](https://iopscience.iop.org/article/10.1088/0067-0049/186/2/457),
[GAMER-2](https://academic.oup.com/mnras/article/481/4/4815/5106358) ,
[GAMER-MHD](http://iopscience.iop.org/article/10.3847/1538-4365/aac49e/meta) ,
[GAMER-SR](https://academic.oup.com/mnras/article/504/3/3298/6224873) ,
diff --git a/doc/wiki/How-to-Release-Code.md b/doc/wiki/How-to-Release-Code.md
new file mode 100644
index 0000000000..911c5d746e
--- /dev/null
+++ b/doc/wiki/How-to-Release-Code.md
@@ -0,0 +1,35 @@
+## Major/Minor release
+In this example, we release version `gamer-2.3.0`, with the previous version being `gamer-2.2.1`.
+
+1. Create a new branch `v2.3.x`
+2. Bump `VERSION` from `gamer-2.3.dev` to `gamer-2.4.dev` in `include/Macro.h` at the `main` branch
+3. Bump `VERSION` from `gamer-2.3.dev` to `gamer-2.3.0` in `include/Macro.h` at the `v2.3.x` branch
+4. Go to `Release`
+[[/images/Release.png | alt=Release]]
+5. Click `Draft a new release`
+[[/images/DraftANewRelease.png | alt=DraftANewRelease]]
+6. Create a new tag `gamer-2.3.0` and select `v2.3.x` as the target branch
+[[/images/ChooseNewTag.png | alt=ChooseNewTag]]
+7. Click `Generate release notes`; polish the content if needed
+[[/images/GenerateReleaseNote_Major_Minor.png | alt=GenerateReleaseNote_Major_Minor]]
+8. Add the release title `GAMER 2.3.0`
+9. Click `Publish release`
+[[/images/PublishRelease.png | alt=PublishRelease]]
+
+
+## Maintenance release
+In this example, we release version `gamer-2.2.2`, with the previous version being `gamer-2.2.1`.
+
+1. Cherry-pick the necessary commits from the latest `main` branch into `v2.2.x`
+2. Bump `VERSION` from `gamer-2.2.1` to `gamer-2.2.2` in `include/Macro.h` at the `v2.2.x` branch
+3. Go to `Release`
+[[/images/Release.png | alt=Release]]
+4. Click `Draft a new release`
+[[/images/DraftANewRelease.png | alt=DraftANewRelease]]
+5. Create a new tag `gamer-2.2.2` and select `v2.2.x` as the target branch
+[[/images/ChooseNewTag.png | alt=ChooseNewTag]]
+6. Click `Generate release notes`; polish the content if needed
+[[/images/GenerateReleaseNote_Maintenance.png | alt=GenerateReleaseNote_Maintenance]]
+7. Add the release title `GAMER 2.2.2`
+8. Click `Publish release`
+[[/images/PublishRelease.png | alt=PublishRelease]]
diff --git a/doc/wiki/Hypre-related/[Hypre]-Develop-in-GAMER.md b/doc/wiki/Hypre-related/[Hypre]-Develop-in-GAMER.md
new file mode 100644
index 0000000000..e69de29bb2
diff --git a/doc/wiki/Hypre-related/[Hypre]-Instllation.md b/doc/wiki/Hypre-related/[Hypre]-Instllation.md
new file mode 100644
index 0000000000..e69de29bb2
diff --git a/doc/wiki/Hypre-related/[Hypre]-Introduction.md b/doc/wiki/Hypre-related/[Hypre]-Introduction.md
new file mode 100644
index 0000000000..df6d147e09
--- /dev/null
+++ b/doc/wiki/Hypre-related/[Hypre]-Introduction.md
@@ -0,0 +1,14 @@
+## Introduction
+
+## Installation example
+
+## How to setup in GAMER
+
+### Compilation options
+
+### Runtime parameters
+
+## Examples
+
+## Future works
+Please contact us if you would like to contribute any of the features.
\ No newline at end of file
diff --git a/doc/wiki/In-Situ-Python-Analysis-related/In-Situ-Python-Analysis:-Plummer.md b/doc/wiki/In-Situ-Python-Analysis-related/In-Situ-Python-Analysis:-Plummer.md
index b6da6416d4..5dec6bef21 100644
--- a/doc/wiki/In-Situ-Python-Analysis-related/In-Situ-Python-Analysis:-Plummer.md
+++ b/doc/wiki/In-Situ-Python-Analysis-related/In-Situ-Python-Analysis:-Plummer.md
@@ -3,7 +3,7 @@ This test problem Plummer demonstrates the in situ Python analysis feature in GA
***
1. Install the required packages for this demo.
- * [libyt and yt_libyt](https://yt-project.github.io/libyt/HowToInstall.html#how-to-install): the in situ analysis library and its yt frontend.
+ * [libyt and yt_libyt](https://libyt.readthedocs.io/en/latest/how-to-install/how-to-install.html): the in situ analysis library and its yt frontend.
> [!IMPORTANT]
> This example assumes it is using libyt interactive mode. Please compile libyt in interactive mode.
* [yt](https://yt-project.org/): the core analysis tool.
diff --git a/doc/wiki/Installation-related/Installation:-External-Libraries.md b/doc/wiki/Installation-related/Installation:-External-Libraries.md
index 2b4d1a84ea..fa726e517d 100644
--- a/doc/wiki/Installation-related/Installation:-External-Libraries.md
+++ b/doc/wiki/Installation-related/Installation:-External-Libraries.md
@@ -2,6 +2,7 @@
* [GRACKLE](#GRACKLE)
* [HDF5](#HDF5)
* [LIBYT](#LIBYT)
+* [HYPRE](#HYPRE)
## FFTW
GAMER supports both FFTW2 and FFTW3 for various calculations (e.g., the root-level Poisson solver).
@@ -91,6 +92,14 @@ Please refer to [libyt -- How to Install](https://libyt.readthedocs.io/en/latest
Set `LIBYT_PATH` to the folder that contains subfolders `include` and `lib`.
+## HYPRE (`HYPRE_PATH`)
+GAMER support to use [hypre](https://github.com/hypre-space/hypre) to Poisson equation.
+See [[ Hypre in GAMER | [Hypre]-Introduction ]] for details.
+
+Please refer to [hypre: General Information](https://hypre.readthedocs.io/en/latest/ch-misc.html).
+
+Set `HYPRE_PATH` to the folder that contains subfolders `include` and `lib`.
+
## Links
diff --git a/doc/wiki/Installation-related/Installation:-Option-List.md b/doc/wiki/Installation-related/Installation:-Option-List.md
index b832283fd6..f7725edc03 100644
--- a/doc/wiki/Installation-related/Installation:-Option-List.md
+++ b/doc/wiki/Installation-related/Installation:-Option-List.md
@@ -71,10 +71,10 @@ disabled). See the "Restriction" of each option carefully.
|
Option
|
Value
|
Default
|
Description
|
Restriction
| Corresponding symbolic constant |
|:---:|:---:|:---:|---|---|---|
-| `--pot_scheme` | `SOR`, `MG` | `SOR` | Poisson solver. `SOR`: successive-overrelaxation (recommended); `MG`: multigrid | Must be set when `--gravity` is enabled | `POT_SCHEME` |
-| `--store_pot_ghost` | `true`, `false` | `true` | Store the ghost-zone potential (recommended when `--particle` is enabled) | Must be enabled when both `--star_formation` and `--store_par_acc` are adopted | `STORE_POT_GHOST` |
-| `--unsplit_gravity` | `true`, `false` | Depend | Use operator-unsplit method to couple gravity to the adopted physical model (recommended) | Not supported for `--model=ELBDM` | `UNSPLIT_GRAVITY` |
-| `--comoving` | `true`, `false` | `false` | Cosmological simulation | - | `COMOVING` |
+| `--pot_scheme` | `SOR`, `MG`, `POI_HYPRE | `SOR` | Poisson solver. `SOR`: successive-overrelaxation (recommended); `MG`: multigrid | Must be set when `--gravity` is enabled; `POI_HYPRE`: use hypre library | `POT_SCHEME` |
+| `--store_pot_ghost` | `true`, `false` | `true` | Store the ghost-zone potential (recommended when `--particle` is enabled) | Must be enabled when both `--star_formation` and `--store_par_acc` are adopted | `STORE_POT_GHOST` |
+| `--unsplit_gravity` | `true`, `false` | Depend | Use operator-unsplit method to couple gravity to the adopted physical model (recommended) | Not supported for `--model=ELBDM` | `UNSPLIT_GRAVITY` |
+| `--comoving` | `true`, `false` | `false` | Cosmological simulation | - | `COMOVING` |
## Particle Options
-- see [[Particles]] for the related runtime parameters and other settings; must enable `--particle=true`
@@ -135,6 +135,7 @@ disabled). See the "Restriction" of each option carefully.
| `--gsl` | `true`, `false` | `false` | Enable GNU scientific library | May need to set `GSL_PATH` in [[configuration file \| Installation:-Machine-Configuration-File#1-Library-paths]] | `SUPPORT_GSL` |
| `--fftw` | `OFF`, `FFTW2`, `FFTW3` | `OFF` | Enable FFTW | May need to set `FFTW2/3_PATH` in [[configuration file \| Installation:-Machine-Configuration-File#1-Library-paths]] | `SUPPORT_FFTW` |
| `--spectral_interpolation` | `true`, `false` | `false` | Enable spectral interpolation | Must enable `--gsl` and set `--fftw` | `SUPPORT_SPECTRAL_INT` |
+| `--hypre` | `true`, `false` | `false` | Enable hypre library | May need to set `HYPRE_PATH` in [[configuration file \| Installation:-Machine-Configuration-File#1-Library-paths]] | `SUPPORT_HYPRE` |
| `--rng` | `RNG_GNU_EXT`, `RNG_CPP11` | `RNG_GNU_EXT` | Random number generators. `RNG_GNU_EXT`: GNU extension `drand48_r`; `RNG_CPP11`: c++11 `` | Use `RNG_GNU_EXT` for compilers supporting GNU extensions (they may not be supported on macOS); use `RNG_CPP11` for compilers supporting c++11 (one may need to add `-std=c++11` to `CXXFLAG` → see [[compilation flags \| Installation:-Machine-Configuration-File#3-Compilation-flags]]) | `RANDOM_NUMBER` |
> [!CAUTION]
> On macOS, we recommend using the GNU compiler and set `--rng=RNG_CPP11`.
diff --git a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md
index 0288725ab3..c51a0e899e 100644
--- a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md
+++ b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md
@@ -139,6 +139,15 @@ For variables with `Default/Min/Max` labeled as `Depend`, click the parameter na
| Name | Default | Min | Max | Short description |
| :--- | :--- | :--- | :--- | :--- |
| [[ HUBBLE0 \| Runtime-Parameters:-Cosmology#HUBBLE0 ]] | -1.0 | 2.22507386e-308 | 1.0 | dimensionless Hubble parameter (currently only for converting ELBDM_MASS to code units) |
+| [[ HYPRE_ABS_TOL \| [Runtime-Parameters]-Hypre#HYPRE_ABS_TOL ]] | -1.0 | None | None | |
+| [[ HYPRE_ENABLE_LOGGING \| [Runtime-Parameters]-Hypre#HYPRE_ENABLE_LOGGING ]] | 1 | 0 | 1 | |
+| [[ HYPRE_INIT_GUESS \| [Runtime-Parameters]-Hypre#HYPRE_INIT_GUESS ]] | 1 | None | None | |
+| [[ HYPRE_MAX_ITER \| [Runtime-Parameters]-Hypre#HYPRE_MAX_ITER ]] | -1 | None | None | |
+| HYPRE_NPOST_RELAX | 1 | None | None | |
+| HYPRE_NPRE_RELAX | 1 | None | None | |
+| [[ HYPRE_PRINT_LEVEL \| [Runtime-Parameters]-Hypre#HYPRE_PRINT_LEVEL ]] | 2 | 0 | None | |
+| [[ HYPRE_REL_TOL \| [Runtime-Parameters]-Hypre#HYPRE_REL_TOL ]] | -1.0 | None | None | |
+| [[ HYPRE_SOLVER \| [Runtime-Parameters]-Hypre#HYPRE_SOLVER ]] | 1 | 1 | 2 | |
# I
| Name | Default | Min | Max | Short description |
diff --git a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Particles.md b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Particles.md
index ed6a212ba0..87947a1e17 100644
--- a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Particles.md
+++ b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Particles.md
@@ -13,7 +13,8 @@ Parameters described on this page:
[PAR_IMPROVE_ACC](#PAR_IMPROVE_ACC),
[PAR_PREDICT_POS](#PAR_PREDICT_POS),
[PAR_REMOVE_CELL](#PAR_REMOVE_CELL),
-[OPT__FREEZE_PAR](#OPT__FREEZE_PAR)
+[OPT__FREEZE_PAR](#OPT__FREEZE_PAR),
+[OPT__PAR_INIT_CHECK](#OPT__PAR_INIT_CHECK)
Parameters below are shown in the format: **`Name` (Valid Values) [Default Value]**
@@ -150,6 +151,13 @@ Do not update particle position and velocity (except for tracer particles).
It can be useful for evolving fluid in a static gravitational potential of particles.
* **Restriction:**
+
+* #### `OPT__PAR_INIT_CHECK` (0=off, 1=on) [1]
+ * **Description:**
+Check particle initialization by calling `Par_Aux_InitCheck()` function after particle initialization (only works for [PAR_INIT](#PAR_INIT)!=`2`).
+Disable this check when particles are initialized _after_ setting grid fields, such as by `Init_User_Ptr()`.
+ * **Restriction:**
+
## Remarks
diff --git a/doc/wiki/Runtime-Parameters-related/[Runtime-Parameters]-Hypre.md b/doc/wiki/Runtime-Parameters-related/[Runtime-Parameters]-Hypre.md
new file mode 100644
index 0000000000..ecb767266b
--- /dev/null
+++ b/doc/wiki/Runtime-Parameters-related/[Runtime-Parameters]-Hypre.md
@@ -0,0 +1,84 @@
+Parameters described on this page:
+Parameters described on this page:
+[HYPRE_SOLVER](#HYPRE_SOLVER),
+[HYPRE_INIT_GUESS](#HYPRE_INIT_GUESS),
+[HYPRE_PRINT_LEVEL](#HYPRE_PRINT_LEVEL),
+[HYPRE_ENABLE_LOGGING](#HYPRE_ENABLE_LOGGING),
+[HYPRE_MAX_ITER](#HYPRE_MAX_ITER),
+[HYPRE_REL_TOL](#HYPRE_REL_TOL),
+[HYPRE_ABS_TOL](#HYPRE_ABS_TOL)
+
+
+Parameters below are shown in the format: **`Name` (Valid Values) [Default Value]**
+
+* #### `HYPRE_SOLVER` (1: SysPFMG) [1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_INIT_GUESS` (0=off, 1=on) [1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_PRINT_LEVEL` () [2]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_ENABLE_LOGGING` (0=off, 1=on) [1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_MAX_ITER` (<0 → set to default) [-1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_NPRE_RELAX` () [1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_NPOST_RELAX` () [1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_REL_TOL` (<0 → set to default) [-1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+* #### `HYPRE_ABS_TOL` (<0 → set to default) [-1]
+ * **Description:**
+TBF
+ * **Restriction:**
+TBF
+
+
+## Remarks
+
+
+
+
+## Links
+* [[ Main page of Runtime Parameters | Runtime Parameters ]]
+* [[ Main page of Hypre | [Hypre]-Intorduction ]]
diff --git a/doc/wiki/_Sidebar.md b/doc/wiki/_Sidebar.md
index 3b8aca343a..8a6e40e3fe 100644
--- a/doc/wiki/_Sidebar.md
+++ b/doc/wiki/_Sidebar.md
@@ -45,8 +45,11 @@
* [[Adding Parameters | Adding-Parameters]]
* [[AMR Structure | AMR-Structure]]
* [[configure.py | configure.py]]
+* [[ Hypre | [Hypre]-Introduction ]]
* [[ELBDM Hybrid Scheme | ELBDM-Hybrid-Scheme]]
* [[ELBDM Spectral Solver | ELBDM-Spectral-Solver]]
* [[ELBDM Spectral Interpolation | ELBDM-Spectral-Interpolation]]
* [[Style Guide | Style-Guide]]
* [[How to Contribute | How-to-Contribute]]
+* [[How to Release Code | How-to-Release-Code]]
+* [[Developing with VS Code | Developing-with-VS-Code]]
diff --git a/doc/wiki/images/ChooseNewTag.png b/doc/wiki/images/ChooseNewTag.png
new file mode 100644
index 0000000000..7b121bc1d2
Binary files /dev/null and b/doc/wiki/images/ChooseNewTag.png differ
diff --git a/doc/wiki/images/DraftANewRelease.png b/doc/wiki/images/DraftANewRelease.png
new file mode 100644
index 0000000000..988dc010f4
Binary files /dev/null and b/doc/wiki/images/DraftANewRelease.png differ
diff --git a/doc/wiki/images/GenerateReleaseNote_Maintenance.png b/doc/wiki/images/GenerateReleaseNote_Maintenance.png
new file mode 100644
index 0000000000..52af25ff1a
Binary files /dev/null and b/doc/wiki/images/GenerateReleaseNote_Maintenance.png differ
diff --git a/doc/wiki/images/GenerateReleaseNote_Major_Minor.png b/doc/wiki/images/GenerateReleaseNote_Major_Minor.png
new file mode 100644
index 0000000000..ad3e3c9e30
Binary files /dev/null and b/doc/wiki/images/GenerateReleaseNote_Major_Minor.png differ
diff --git a/doc/wiki/images/PublishRelease.png b/doc/wiki/images/PublishRelease.png
new file mode 100644
index 0000000000..c18a2482b0
Binary files /dev/null and b/doc/wiki/images/PublishRelease.png differ
diff --git a/doc/wiki/images/Release.png b/doc/wiki/images/Release.png
new file mode 100644
index 0000000000..25e3d34097
Binary files /dev/null and b/doc/wiki/images/Release.png differ
diff --git a/example/test_problem/Hydro/JeansInstability/Input__Flag_Rho b/example/test_problem/Hydro/JeansInstability/Input__Flag_Rho
new file mode 100644
index 0000000000..4747aa812d
--- /dev/null
+++ b/example/test_problem/Hydro/JeansInstability/Input__Flag_Rho
@@ -0,0 +1,13 @@
+# Level Density
+ 0 2.0000009
+ 1 16.0
+ 2 64.0
+ 3 256.0
+ 4 1024.0
+ 5 4096.0
+ 6 16384.0
+ 7 65536.0
+ 8 262144.0
+ 9 1048576.0
+ 10 4194304.0
+ 11 16777216.0
diff --git a/example/test_problem/Hydro/JeansInstability/Input__Parameter b/example/test_problem/Hydro/JeansInstability/Input__Parameter
index b370bc1fc7..f07e30a060 100644
--- a/example/test_problem/Hydro/JeansInstability/Input__Parameter
+++ b/example/test_problem/Hydro/JeansInstability/Input__Parameter
@@ -15,12 +15,12 @@
# simulation scale
BOX_SIZE 1.0 # box size along the longest side (in Mpc/h if COMOVING is adopted)
-NX0_TOT_X 64 # number of base-level cells along x
-NX0_TOT_Y 64 # number of base-level cells along y
-NX0_TOT_Z 64 # number of base-level cells along z
+NX0_TOT_X 32 # number of base-level cells along x
+NX0_TOT_Y 16 # number of base-level cells along y
+NX0_TOT_Z 16 # number of base-level cells along z
OMP_NTHREAD -1 # number of OpenMP threads (<=0=auto) [-1] ##OPENMP ONLY##
END_T -1.0 # end physical time (<0=auto -> must be set by test problems or restart) [-1.0]
-END_STEP -1 # end step (<0=auto -> must be set by test problems or restart) [-1]
+END_STEP 0 # end step (<0=auto -> must be set by test problems or restart) [-1]
# test problems
@@ -52,7 +52,7 @@ DT__SYNC_PARENT_LV 0.1 # dt criterion: allow dt to adjust by
DT__SYNC_CHILDREN_LV 0.1 # dt criterion: allow dt to adjust by (1.0-DT__SYNC_CHILDREN) in order to synchronize
# with the children level (for OPT__DT_LEVEL==3 only; 0=off) [0.1]
OPT__DT_USER 0 # dt criterion: user-defined -> edit "Mis_GetTimeStep_UserCriteria.cpp" [0]
-OPT__DT_LEVEL 3 # dt at different AMR levels (1=shared, 2=differ by two, 3=flexible) [3]
+OPT__DT_LEVEL 1 # dt at different AMR levels (1=shared, 2=differ by two, 3=flexible) [3]
OPT__RECORD_DT 1 # record info of the dt determination [1]
AUTO_REDUCE_DT 1 # reduce dt automatically when the program fails (for OPT__DT_LEVEL==3 only) [1]
AUTO_REDUCE_DT_FACTOR 1.0 # reduce dt by a factor of AUTO_REDUCE_DT_FACTOR when the program fails [1.0]
@@ -64,12 +64,12 @@ AUTO_REDUCE_INT_MONO_MIN 1.0e-2 # minimum allowed INT_MONO_COEFF(_B) a
# grid refinement (examples of Input__Flag_XXX tables are put at "example/input/")
-REGRID_COUNT 4 # refine every REGRID_COUNT sub-step [4]
+REGRID_COUNT 2 # refine every REGRID_COUNT sub-step [4]
FLAG_BUFFER_SIZE -1 # number of buffer cells for the flag operation (0~PATCH_SIZE; <0=auto -> PATCH_SIZE) [-1]
FLAG_BUFFER_SIZE_MAXM1_LV -1 # FLAG_BUFFER_SIZE at the level MAX_LEVEL-1 (<0=auto -> REGRID_COUNT) [-1]
FLAG_BUFFER_SIZE_MAXM2_LV -1 # FLAG_BUFFER_SIZE at the level MAX_LEVEL-2 (<0=auto) [-1]
-MAX_LEVEL 0 # maximum refinement level (0~NLEVEL-1) [NLEVEL-1]
-OPT__FLAG_RHO 0 # flag: density (Input__Flag_Rho) [0]
+MAX_LEVEL 1 # maximum refinement level (0~NLEVEL-1) [NLEVEL-1]
+OPT__FLAG_RHO 1 # flag: density (Input__Flag_Rho) [0]
OPT__FLAG_RHO_GRADIENT 0 # flag: density gradient (Input__Flag_RhoGradient) [0]
OPT__FLAG_PRES_GRADIENT 0 # flag: pressure gradient (Input__Flag_PresGradient) [0] ##HYDRO ONLY##
OPT__FLAG_VORTICITY 0 # flag: vorticity (Input__Flag_Vorticity) [0] ##HYDRO ONLY##
@@ -115,7 +115,7 @@ GPU_NSTREAM -1 # number of CUDA streams for the async
OPT__FIXUP_FLUX 1 # correct coarse grids by the fine-grid boundary fluxes [1] ##HYDRO and ELBDM ONLY##
OPT__FIXUP_ELECTRIC 1 # correct coarse grids by the fine-grid boundary electric field [1] ##MHD ONLY##
OPT__FIXUP_RESTRICT 1 # correct coarse grids by averaging the fine-grid data [1]
-OPT__CORR_AFTER_ALL_SYNC -1 # apply various corrections after all levels are synchronized (see "Flu_CorrAfterAllSync"):
+OPT__CORR_AFTER_ALL_SYNC 2 # apply various corrections after all levels are synchronized (see "Flu_CorrAfterAllSync"):
# (-1=auto, 0=off, 1=every step, 2=before dump) [-1]
OPT__NORMALIZE_PASSIVE 1 # ensure "sum(passive_scalar_density) == gas_density" [1]
OPT__OVERLAP_MPI 0 # overlap MPI communication with CPU/GPU computations [0] ##NOT SUPPORTED YET##
@@ -167,7 +167,7 @@ MONO_MAX_ITER 10 # maximum number of iterations to redu
# data dump
-OPT__OUTPUT_TOTAL 0 # output the simulation snapshot: (0=off, 1=HDF5, 2=C-binary) [1]
+OPT__OUTPUT_TOTAL 1 # output the simulation snapshot: (0=off, 1=HDF5, 2=C-binary) [1]
OPT__OUTPUT_PART 0 # output a single line or slice: (0=off, 1=xy, 2=yz, 3=xz, 4=x, 5=y, 6=z, 7=diag, 8=entire box) [0]
OPT__OUTPUT_USER 1 # output the user-specified data -> edit "Output_User.cpp" [0]
OPT__OUTPUT_BASEPS 0 # output the base-level power spectrum [0]
@@ -185,7 +185,7 @@ OPT__OUTPUT_LORENTZ 0 # output Lorentz factor [0] ##SRHD ONL
OPT__OUTPUT_3VELOCITY 0 # output 3-velocities [0] ##SRHD ONLY##
OPT__OUTPUT_USER_FIELD 0 # output user-defined derived fields [0] -> edit "Flu_DerivedField_User.cpp"
OPT__OUTPUT_MODE 1 # (1=const step, 2=const dt, 3=dump table) -> edit "Input__DumpTable" for 3
-OUTPUT_STEP 5 # output data every OUTPUT_STEP step ##OPT__OUTPUT_MODE==1 ONLY##
+OUTPUT_STEP 1 # output data every OUTPUT_STEP step ##OPT__OUTPUT_MODE==1 ONLY##
OUTPUT_DT 1.0 # output data every OUTPUT_DT time interval ##OPT__OUTPUT_MODE==2 ONLY##
OUTPUT_PART_X 0.0 # x coordinate for OPT__OUTPUT_PART [-1.0]
OUTPUT_PART_Y 0.0 # y coordinate for OPT__OUTPUT_PART [-1.0]
@@ -195,7 +195,7 @@ OUTPUT_DIR . # set the output directory [.]
# miscellaneous
-OPT__VERBOSE 0 # output the simulation progress in detail [0]
+OPT__VERBOSE 1 # output the simulation progress in detail [0]
OPT__TIMING_BARRIER -1 # synchronize before timing -> more accurate, but may slow down the run (<0=auto) [-1]
OPT__TIMING_BALANCE 0 # record the max/min elapsed time in various code sections for checking load balance [0]
OPT__TIMING_MPI 0 # record the MPI bandwidth achieved in various code sections [0] ##LOAD_BALANCE ONLY##
@@ -222,3 +222,11 @@ OPT__CK_NEGATIVE 0 # check the negative values: (0=off, 1
OPT__CK_MEMFREE 1.0 # check the free memory in GB (0=off, >0=threshold) [1.0]
OPT__CK_INTERFACE_B 0 # check the consistency of patch interface B field [0] ##MHD ONLY##
OPT__CK_DIVERGENCE_B 0 # check the divergence-free constraint on B field (0=off, 1=on, 2=on+verbose) [0] ##MHD ONLY##
+
+# Hypre
+HYPRE_SOLVER 1
+HYPRE_MAX_ITER 100
+HYPRE_REL_TOL 1.e-10
+
+DT__MAX 1.0
+PT__FREEZE_FLUID 1
diff --git a/example/test_problem/Hydro/JeansInstability/Input__TestProb b/example/test_problem/Hydro/JeansInstability/Input__TestProb
index 085bf572c4..5c5cdb5c2f 100644
--- a/example/test_problem/Hydro/JeansInstability/Input__TestProb
+++ b/example/test_problem/Hydro/JeansInstability/Input__TestProb
@@ -4,6 +4,6 @@ Jeans_Rho1 1.0e-6 # density perturbation amplitude
Jeans_P0 1.0e-2 # background pressure
Jeans_v0 0.0 # background velocity [0.0]
Jeans_B0 1.0e-2 # background magnetic field amplitude
-Jeans_Dir 3 # wave direction: (0/1/2/3) --> (x/y/z/diagonal)
+Jeans_Dir 0 # wave direction: (0/1/2/3) --> (x/y/z/diagonal)
Jeans_Sign +1.0 # (>0/<0) --> (stable:right/left-moving wave; unstable:growing/decaying mode) [+1]
-Jeans_Phase0 0.2 # initial phase shift [0.0]
+Jeans_Phase0 1.5 #0.2 # initial phase shift [0.0]
diff --git a/example/test_problem/Template/Input__Parameter b/example/test_problem/Template/Input__Parameter
index ad73dfe255..79a0253114 100644
--- a/example/test_problem/Template/Input__Parameter
+++ b/example/test_problem/Template/Input__Parameter
@@ -421,6 +421,18 @@ YT_FIG_BASENAME Fig # figure basename [Fig]
YT_JUPYTER_USE_CONNECTION_FILE 0 # use user-provided connection file when using libyt Jupyter UI [0]
+# Hypre
+HYPRE_SOLVER 1 #
+HYPRE_INIT_GUESS 1 #
+HYPRE_PRINT_LEVEL 2 #
+HYPRE_ENABLE_LOGGING 1 #
+HYPRE_MAX_ITER 100 #
+HYPRE_NPRE_RELAX 1 #
+HYPRE_NPOST_RELAX 1 #
+HYPRE_ABS_TOL -1 #
+HYPRE_REL_TOL -1 #
+
+
# miscellaneous
OPT__VERBOSE 0 # output the simulation progress in detail [0]
OPT__TIMING_BARRIER -1 # synchronize before timing -> more accurate, but may slow down the run (<0=auto) [-1]
diff --git a/include/GAMER.h b/include/GAMER.h
index eb1aae3622..26cce41747 100644
--- a/include/GAMER.h
+++ b/include/GAMER.h
@@ -61,6 +61,10 @@ extern "C" {
# include
#endif
+#ifdef SUPPORT_HYPRE
+# include "Hypre.h"
+#endif
+
#include "Macro.h"
#include "Typedef.h"
#include "AMR.h"
diff --git a/include/GatherTree.h b/include/GatherTree.h
index cf10c165d5..3d35faae67 100644
--- a/include/GatherTree.h
+++ b/include/GatherTree.h
@@ -19,6 +19,8 @@ class NonCopyable
struct LB_GlobalPatch
{
int corner[3];
+ int cornerL[3];
+ int cornerR[3];
int sibling[26];
int father;
int son;
@@ -61,6 +63,8 @@ struct LB_LocalPatchExchangeList : private NonCopyable
long *LBIdxList_Local [NLEVEL]; // load balance ids
int (*CrList_Local [NLEVEL])[3]; // patch corners
+ int (*CrLList_Local [NLEVEL])[3]; // patch left corners
+ int (*CrRList_Local [NLEVEL])[3]; // patch right corners
int *FaList_Local [NLEVEL]; // father GIDs
int *SonList_Local [NLEVEL]; // son GIDs
int (*SibList_Local [NLEVEL])[26]; // sibling GIDs
@@ -88,6 +92,8 @@ struct LB_GlobalPatchExchangeList : private NonCopyable
long *LBIdxList_AllLv; // load balance ids
int (*CrList_AllLv )[3]; // patch corners
+ int (*CrLList_AllLv )[3]; // patch corners
+ int (*CrRList_AllLv )[3]; // patch corners
int *FaList_AllLv; // father GIDs
int *SonList_AllLv; // son GIDs
int (*SibList_AllLv )[26]; // sibling GIDs
diff --git a/include/Global.h b/include/Global.h
index e61e8da28e..6a19eb32e7 100644
--- a/include/Global.h
+++ b/include/Global.h
@@ -8,6 +8,9 @@
#ifdef SUPPORT_LIBYT
#include
#endif
+#ifdef SUPPORT_HYPRE
+#include "Hypre.h"
+#endif
#include "GatherTree.h"
@@ -181,7 +184,9 @@ extern double AveDensity_Init; // initial average mass density (in al
extern int Pot_ParaBuf; // number of parallel buffers to exchange potential for the
// Poisson/Gravity solvers and the potential refinement
extern int Rho_ParaBuf; // number of parallel buffers to exchange density for the Poisson solver
-extern real *GreenFuncK;
+extern bool FFTW_Inited[NLEVEL];
+extern bool GreenFuncK_Inited[NLEVEL];
+extern real *GreenFuncK[NLEVEL];
extern double GFUNC_COEFF0;
extern double DT__GRAVITY;
extern double NEWTON_G;
@@ -365,7 +370,9 @@ extern bool FB_Any;
extern int FB_ParaBuf;
#endif
+
// (2-13) spectral interpolation
+// =======================================================================================================
#ifdef SUPPORT_SPECTRAL_INT
extern char SPEC_INT_TABLE_PATH[MAX_STRING];
extern int SPEC_INT_GHOST_BOUNDARY;
@@ -378,7 +385,7 @@ extern InterpolationHandler Int_InterpolationHandler;
#endif // #ifdef SUPPORT_SPECTRAL_INT
-// (2-13) cosmic ray
+// (2-14) cosmic ray
// =======================================================================================================
#ifdef COSMIC_RAY
extern double GAMMA_CR;
@@ -387,7 +394,7 @@ extern double FlagTable_CRay[NLEVEL-1];
#endif
-// (2-14) microphysics
+// (2-15) microphysics
// =======================================================================================================
extern MicroPhy_t MicroPhy;
#ifdef CR_DIFFUSION
@@ -398,6 +405,24 @@ extern double CR_DIFF_MIN_B;
#endif
+// (2-16) Hypre
+// =======================================================================================================
+#ifdef SUPPORT_HYPRE
+extern Hypre_Solver_t HYPRE_SOLVER;
+extern bool HYPRE_INIT_GUESS;
+extern int HYPRE_PRINT_LEVEL, HYPRE_ENABLE_LOGGING;
+extern int HYPRE_MAX_ITER, HYPRE_NPRE_RELAX, HYPRE_NPOST_RELAX;
+extern double HYPRE_REL_TOL, HYPRE_ABS_TOL;
+extern HYPRE_SStructGrid Hypre_grid;
+extern HYPRE_SStructGraph Hypre_graph;
+extern HYPRE_SStructStencil Hypre_stencil;
+extern HYPRE_SStructMatrix Hypre_A;
+extern HYPRE_SStructVector Hypre_b;
+extern HYPRE_SStructVector Hypre_x;
+extern HYPRE_SStructSolver Hypre_solver;
+#endif
+
+
// 3. CPU (host) arrays for transferring data between CPU and GPU
// ============================================================================================================
diff --git a/include/HDF5_Typedef.h b/include/HDF5_Typedef.h
index f6169399b3..45bc48d501 100644
--- a/include/HDF5_Typedef.h
+++ b/include/HDF5_Typedef.h
@@ -153,6 +153,7 @@ struct Makefile_t
int LibYTJupyter;
# endif
int SupportGrackle;
+ int SupportHypre;
int RandomNumber;
# ifdef GRAVITY
@@ -847,6 +848,17 @@ struct InputPara_t
int Yt_JupyterUseConnectionFile;
# endif
+// Hypre
+# ifdef SUPPORT_HYPRE
+ int Hypre_Solver;
+ int Hypre_InitGuess;
+ int Hypre_PrintLevel;
+ int Hypre_EnableLogging;
+ int Hypre_MaxIter;
+ double Hypre_RelTol;
+ double Hypre_AbsTol;
+# endif
+
// miscellaneous
int Opt__Verbose;
int Opt__TimingBarrier;
diff --git a/include/Hypre.h b/include/Hypre.h
new file mode 100644
index 0000000000..b21c07757d
--- /dev/null
+++ b/include/Hypre.h
@@ -0,0 +1,47 @@
+#ifndef __HYPRE_H__
+#define __HYPRE_H__
+
+
+
+#ifdef GAMER_DEBUG
+# define DEBUG_HYPRE
+#endif
+
+#ifndef SERIAL
+#define HYPRE_MPI_COMM MPI_COMM_WORLD
+#else
+#define HYPRE_MPI_COMM NULL_INT
+#endif
+
+#ifdef GPU
+#include
+#endif
+
+
+#include "HYPRE_config.h"
+#include "HYPRE_sstruct_ls.h"
+#include "HYPRE_sstruct_mv.h"
+#include "HYPRE_struct_ls.h"
+#include "HYPRE_struct_mv.h"
+
+
+// macro for checking Hypre function return
+#ifdef DEBUG_HYPRE
+
+# define HYPRE_CHECK_FUNC( call ) \
+ { \
+ int hypre_ierr = call; \
+ \
+ if ( hypre_ierr ) \
+ Aux_Error( ERROR_INFO, "hypre error code: %d !!\n", hypre_ierr ); \
+ }
+
+#else
+
+# define HYPRE_CHECK_FUNC( call ) call
+
+#endif
+
+
+
+#endif // #ifndef __HYPRE_H__
diff --git a/include/Macro.h b/include/Macro.h
index a0aa3716ed..fd6d8a698a 100644
--- a/include/Macro.h
+++ b/include/Macro.h
@@ -131,6 +131,7 @@
// Poisson solvers
#define SOR 1
#define MG 2
+#define HYPRE_POI 3
// load-balance parallelization
diff --git a/include/Patch.h b/include/Patch.h
index 5a4c7badcf..f43debf981 100644
--- a/include/Patch.h
+++ b/include/Patch.h
@@ -78,6 +78,8 @@ long LB_Corner2Index( const int lv, const int Corner[], const Check_t Check );
// convert corner to that of the corresponding real patch
// --> For example, external patches can have corner < 0
// --> Different from EdgeL/R, which always assume periodicity
+// cornerL[3] : Grid indices of the cell at the left patch corner on the current level
+// cornerR[3] : Grid indices of the cell at the right patch corner on the current level
// sibling[26] : Patch IDs of the 26 sibling patches (-1->no sibling; -1XX->external)
//
// NOTE FOR NON-PERIODIC BOUNDARY CONDITIONS:
@@ -225,6 +227,8 @@ struct patch_t
# endif
int corner[3];
+ int cornerL[3];
+ int cornerR[3];
int sibling[26];
int father;
int son;
@@ -378,6 +382,9 @@ struct patch_t
EdgeL[d] = BoxEdgeL[d] + (double)( ( corner[d] + BoxScale[d] ) % BoxScale[d] )*dh_min;
EdgeR[d] = BoxEdgeL[d] + (double)( ( corner[d] + BoxScale[d] ) % BoxScale[d] + PScale )*dh_min;
+ cornerL[d] = corner[d] / (1<<(NLEVEL-1-lv));
+ cornerR[d] = cornerL[d] + PS1 - 1;
+
// do no use the following non-periodic version anymore --> it does not work with the current particle implementation
/*
EdgeL[d] = BoxEdgeL[d] + (double)corner[d]*dh_min;
diff --git a/include/Prototype.h b/include/Prototype.h
index af1ba95c9b..ad14208c51 100644
--- a/include/Prototype.h
+++ b/include/Prototype.h
@@ -259,14 +259,15 @@ void Init_ByRestart_HDF5( const char *FileName );
#endif
#ifdef SUPPORT_FFTW
void End_FFTW();
-void Init_FFTW();
+void Init_FFTW( const int lv );
void Patch2Slab( real *VarS, real *SendBuf_Var, real *RecvBuf_Var, long *SendBuf_SIdx, long *RecvBuf_SIdx,
int **List_PID, int **List_k, long *List_NSend_Var, long *List_NRecv_Var,
const int *List_z_start, const int local_nz, const int FFT_Size[], const int NRecvSlice,
- const double PrepTime, const long TVar, const bool InPlacePad, const bool ForPoisson, const bool AddExtraMass );
+ const double PrepTime, const long TVar, const bool InPlacePad, const bool ForPoisson,
+ const bool AddExtraMass, const int lv );
void Slab2Patch( const real *VarS, real *SendBuf, real *RecvBuf, const int SaveSg, const long *List_SIdx,
int **List_PID, int **List_k, long *List_NSend, long *List_NRecv, const int local_nz, const int FFT_Size[],
- const int NSendSlice, const long TVar, const bool InPlacePad );
+ const int NSendSlice, const long TVar, const bool InPlacePad, const int lv );
#endif // #ifdef SUPPORT_FFTW
void Microphysics_Init();
void Microphysics_End();
@@ -404,12 +405,12 @@ void CPU_PoissonGravitySolver( const real h_Rho_Array [][RHO_NXT][RHO_NXT][RH
const bool SelfGravity, const OptExtPot_t ExtPot, const OptExtAcc_t ExtAcc,
const double TimeNew, const double TimeOld, const real MinEint,
const bool UseWaveFlag );
-void CPU_ExtPotSolver_BaseLevel( const ExtPot_t Func, const double AuxArray_Flt[], const int AuxArray_Int[],
- const real Table[], void **GenePtr,
- const double Time, const bool PotIsInit, const int SaveSg );
+void CPU_ExtPotSolver_FullyRefinedLevel( const ExtPot_t Func, const double AuxArray_Flt[], const int AuxArray_Int[],
+ const real Table[], void **GenePtr,
+ const double Time, const bool PotIsInit, const int SaveSg, const int lv );
#ifdef SUPPORT_FFTW
-void CPU_PoissonSolver_FFT( const real Poi_Coeff, const int SaveSg, const double PrepTime );
-void Init_GreenFuncK();
+void CPU_PoissonSolver_FFT( const real Poi_Coeff, const int SaveSg, const double PrepTime, const int lv );
+void Init_GreenFuncK( const int lv );
#endif
void End_MemFree_PoissonGravity();
void Gra_AdvanceDt( const int lv, const double TimeNew, const double TimeOld, const double dt,
@@ -844,6 +845,24 @@ int FB_Aux_CellPatchRelPos( const int ijk[] );
#endif
+// Hypre
+#ifdef SUPPORT_HYPRE
+void Hypre_Init();
+void Hypre_Solver( const Hypre_SolveType_t SolveType, const int lv, const double TimeNew, const double TimeOld,
+ const double dt_in, const double Poi_Coeff, const int SaveSg_Flu, const int SaveSg_Mag,
+ const int SaveSg_Pot );
+void Hypre_PrepareSingleLevel( const Hypre_SolveType_t, const int lv );
+void Hypre_FillArrays( const Hypre_SolveType_t SolveType, const int lv, const double TimeNew, const real Poi_Coeff,
+ const int SaveSg_Pot );
+void Hypre_Solve( const Hypre_Solver_t Solver, int *N_iter, real_hypre *final_res_norm );
+void Hypre_UpdateArrays( const Hypre_SolveType_t SolveType, const int lv, const int SaveSg_Flu,
+ const int SaveSg_Mag, const int SaveSg_Pot );
+void Hypre_Free();
+void Hypre_Aux_Record( const Hypre_SolveType_t SolveType, const int lv, const int iteration, const real residual );
+void Hypre_End();
+#endif
+
+
// EoS in hydrodynamics
#if ( MODEL == HYDRO )
void EoS_Init();
diff --git a/include/Typedef.h b/include/Typedef.h
index 6fd92c29ad..41426f424f 100644
--- a/include/Typedef.h
+++ b/include/Typedef.h
@@ -45,6 +45,15 @@ typedef float real_che;
#endif
#endif // #ifdef SUPPORT_GRACKLE
+#ifdef SUPPORT_HYPRE
+#include "HYPRE_config.h"
+#ifdef HYPRE_SINGLE
+typedef float real_hypre;
+#else
+typedef double real_hypre;
+#endif
+#endif // #ifdef SUPPORT_HYPRE
+
#if ( GRAMFE_SCHEME == GRAMFE_FFT )
#ifdef GRAMFE_FFT_FLOAT8
typedef double gramfe_fft_float;
@@ -562,6 +571,20 @@ const LoadParaMode_t
LOAD_HDF5_OUTPUT = 2;
+// Hypre
+#ifdef SUPPORT_HYPRE
+// solve type
+typedef int Hypre_SolveType_t;
+const Hypre_SolveType_t
+ HYPRE_SOLVE_TYPE_POISSON = 1;
+
+// solver
+typedef int Hypre_Solver_t;
+const Hypre_Solver_t
+ HYPRE_SOLVER_SSTRUCT_SYS_PFMG = 1,
+ HYPRE_SOLVER_SSTRUCT_SPLIT = 2;
+#endif
+
// function pointers
typedef real (*EoS_GUESS_t) ( const real Con[], real* const Constant, const double AuxArray_Flt[],
const int AuxArray_Int[], const real *const Table[EOS_NTABLE_MAX] );
diff --git a/src/Auxiliary/Aux_Check_Parameter.cpp b/src/Auxiliary/Aux_Check_Parameter.cpp
index a018831b4b..6259bdccd2 100644
--- a/src/Auxiliary/Aux_Check_Parameter.cpp
+++ b/src/Auxiliary/Aux_Check_Parameter.cpp
@@ -1524,8 +1524,16 @@ void Aux_Check_Parameter()
# error : ERROR : SUPPORT_FFTW != FFTW2/FFTW3 !!
# endif
-# if ( POT_SCHEME != SOR && POT_SCHEME != MG )
-# error : ERROR : unsupported Poisson solver in the makefile (SOR/MG) !!
+# if ( POT_SCHEME != SOR && POT_SCHEME != MG && POT_SCHEME != HYPRE_POI )
+# error : ERROR : unsupported Poisson solver in the makefile (SOR/MG/HYPRE_POI) !!
+# endif
+
+# if ( POT_SCHEME == HYPRE_POI && !defined SUPPORT_HYPRE )
+# error : ERROR : must enable SUPPORT_HYPRE for POT_SCHEME == HYPRE_POI !!
+# endif
+
+# if ( POT_SCHEME == HYPRE_POI && defined STORE_POT_GHOST )
+# error : STORE_POT_GHOST is useless for POT_SCHEME == HYPRE_POI (if you really need this feature please comment out this line) !!
# endif
# if ( POT_GHOST_SIZE <= GRA_GHOST_SIZE )
diff --git a/src/Auxiliary/Aux_TakeNote.cpp b/src/Auxiliary/Aux_TakeNote.cpp
index 1b3a165c3a..c14c1abda5 100644
--- a/src/Auxiliary/Aux_TakeNote.cpp
+++ b/src/Auxiliary/Aux_TakeNote.cpp
@@ -73,6 +73,8 @@ void Aux_TakeNote()
fprintf( Note, "POT_SCHEME SOR\n" );
# elif ( POT_SCHEME == MG )
fprintf( Note, "POT_SCHEME MG\n" );
+# elif ( POT_SCHEME == HYPRE_POI )
+ fprintf( Note, "POT_SCHEME HYPRE_POI\n" );
# elif ( POT_SCHEME == NONE )
fprintf( Note, "POT_SCHEME NONE\n" );
# else
@@ -480,6 +482,12 @@ void Aux_TakeNote()
fprintf( Note, "SUPPORT_LIBYT OFF\n" );
# endif // #ifdef SUPPORT_LIBYT ... else ...
+# ifdef SUPPORT_HYPRE
+ fprintf( Note, "SUPPORT_HYPRE ON\n" );
+# else
+ fprintf( Note, "SUPPORT_HYPRE OFF\n" );
+# endif
+
# if ( RANDOM_NUMBER == RNG_GNU_EXT )
fprintf( Note, "RANDOM_NUMBER RNG_GNU_EXT\n" );
# elif ( RANDOM_NUMBER == RNG_CPP11 )
@@ -1663,6 +1671,26 @@ void Aux_TakeNote()
# endif
+// record the parameters of Hypre
+# ifdef SUPPORT_HYPRE
+ fprintf( Note, "Parameters of Hypre\n" );
+ fprintf( Note, "***********************************************************************************\n" );
+ fprintf( Note, "HYPRE_SOLVER %s\n", ( HYPRE_SOLVER == HYPRE_SOLVER_SSTRUCT_SYS_PFMG ) ? "SSTRUCT_SYS_PFMG" :
+ ( HYPRE_SOLVER == HYPRE_SOLVER_SSTRUCT_SPLIT ) ? "SSTRUCT_SPLIT" :
+ "UNKNOWN" );
+ fprintf( Note, "HYPRE_INIT_GUESS % d\n", HYPRE_INIT_GUESS );
+ fprintf( Note, "HYPRE_PRINT_LEVEL % d\n", HYPRE_PRINT_LEVEL );
+ fprintf( Note, "HYPRE_ENABLE_LOGGING % d\n", HYPRE_ENABLE_LOGGING );
+ fprintf( Note, "HYPRE_MAX_ITER % d\n", HYPRE_MAX_ITER );
+ fprintf( Note, "HYPRE_NPRE_RELAX % d\n", HYPRE_NPRE_RELAX );
+ fprintf( Note, "HYPRE_NPOST_RELAX % d\n", HYPRE_NPOST_RELAX );
+ fprintf( Note, "HYPRE_REL_TOL % 21.14e\n", HYPRE_REL_TOL );
+ fprintf( Note, "HYPRE_ABS_TOL % 21.14e\n", HYPRE_ABS_TOL );
+ fprintf( Note, "***********************************************************************************\n" );
+ fprintf( Note, "\n\n");
+# endif
+
+
// record the parameters of miscellaneous purposes
fprintf( Note, "Parameters of Miscellaneous Purposes\n" );
fprintf( Note, "***********************************************************************************\n" );
diff --git a/src/GPU_API/CUAPI_Asyn_PoissonGravitySolver.cu b/src/GPU_API/CUAPI_Asyn_PoissonGravitySolver.cu
index 138718bd04..90a82a70b5 100644
--- a/src/GPU_API/CUAPI_Asyn_PoissonGravitySolver.cu
+++ b/src/GPU_API/CUAPI_Asyn_PoissonGravitySolver.cu
@@ -414,6 +414,11 @@ void CUAPI_Asyn_PoissonGravitySolver( const real h_Rho_Array [][RHO_NXT][RHO_
dh, MG_Max_Iter, MG_NPre_Smooth, MG_NPost_Smooth, MG_Tolerated_Error,
Poi_Coeff, IntScheme );
+
+# elif ( POT_SCHEME == HYPRE_POI )
+
+ Aux_Error( ERROR_INFO, "POT_SCHEME == HYPRE_POI is not supported yet !!\n" );
+
# else
# error : unsupported GPU Poisson solver
diff --git a/src/Hypre/Hypre_Aux_Record.cpp b/src/Hypre/Hypre_Aux_Record.cpp
new file mode 100644
index 0000000000..e49b277043
--- /dev/null
+++ b/src/Hypre/Hypre_Aux_Record.cpp
@@ -0,0 +1,112 @@
+#include "GAMER.h"
+
+
+
+
+#ifdef SUPPORT_HYPRE
+
+
+
+//-------------------------------------------------------------------------------------------------------
+// Function : Hypre_Aux_Record
+// Description : TBF
+//
+// Note : TBF
+//
+// Parameter : SolveType :
+// lv :
+// iteration :
+// residual :
+//-------------------------------------------------------------------------------------------------------
+void Hypre_Aux_Record( const Hypre_SolveType_t SolveType, const int lv, const int iteration, const real residual )
+{
+
+ static bool FirstTime = true;
+ char FileName[2*MAX_STRING];
+ sprintf( FileName, "%s/Record__Hypre", OUTPUT_DIR );
+
+ char SolveName[MAX_STRING];
+ switch ( SolveType )
+ {
+# ifdef GRAVITY
+ case HYPRE_SOLVE_TYPE_POISSON: sprintf( SolveName, "Poisson" ); break;
+# endif
+ default :
+ Aux_Error( ERROR_INFO, "incorrect parameter %s = %d !!\n", "SolveType", SolveType );
+ } // switch ( SolveType )
+
+ if ( FirstTime )
+ {
+ if ( MPI_Rank == 0 )
+ {
+ if ( Aux_CheckFileExist( FileName ) )
+ Aux_Message( stderr, "WARNING : file \"%s\" already exists !!\n", FileName );
+
+ FILE *File = fopen( FileName, "a" );
+ fprintf( File, "# ====================================================================================================\n");
+ fprintf( File, "# Hypre package information\n" );
+ fprintf( File, "# ====================================================================================================\n");
+# ifdef HYPRE_RELEASE_DATE
+ fprintf( File, "# %-30s : %s\n", "HYPRE_RELEASE_DATE", HYPRE_RELEASE_DATE );
+# endif
+# ifdef HYPRE_RELEASE_NUMBER
+ fprintf( File, "# %-30s : %d\n", "HYPRE_RELEASE_NUMBER", HYPRE_RELEASE_NUMBER );
+# endif
+# ifdef HYPRE_RELEASE_VERSION
+ fprintf( File, "# %-30s : %s\n", "HYPRE_RELEASE_VERSION", HYPRE_RELEASE_VERSION );
+# endif
+
+# ifdef HYPRE_SINGLE
+ fprintf( File, "# %-30s : %d\n", "HYPRE_SINGLE", HYPRE_SINGLE );
+# else
+ fprintf( File, "# %-30s : %d\n", "HYPRE_SINGLE", 0 );
+# endif
+# ifdef HYPRE_DEBUG
+ fprintf( File, "# %-30s : %d\n", "HYPRE_DEBUG", HYPRE_DEBUG );
+# else
+ fprintf( File, "# %-30s : %d\n", "HYPRE_DEBUG", 0 );
+# endif
+# ifdef HYPRE_HAVE_MPI
+ fprintf( File, "# %-30s : %d\n", "HYPRE_HAVE_MPI", HYPRE_HAVE_MPI );
+# else
+ fprintf( File, "# %-30s : %d\n", "HYPRE_HAVE_MPI", 0 );
+# endif
+# ifdef HYPRE_USING_OPENMP
+ fprintf( File, "# %-30s : %d\n", "HYPRE_USING_OPENMP", HYPRE_USING_OPENMP );
+# else
+ fprintf( File, "# %-30s : %d\n", "HYPRE_USING_OPENMP", 0 );
+# endif
+# ifdef HYPRE_USING_CUDA
+ fprintf( File, "# %-30s : %d\n", "HYPRE_USING_CUDA", HYPRE_USING_CUDA );
+# else
+ fprintf( File, "# %-30s : %d\n", "HYPRE_USING_CUDA", 0 );
+# endif
+
+ fprintf( File, "# ====================================================================================================\n");
+ fprintf( File, "# Hypre runtime parameters\n" );
+ fprintf( File, "# ====================================================================================================\n");
+ fprintf( File, "# Maximum iteration: %d\n", HYPRE_MAX_ITER );
+ fprintf( File, "# Relative tolerance: %22.16e\n", HYPRE_REL_TOL );
+ fprintf( File, "# Absolute tolerance: %22.16e\n", HYPRE_ABS_TOL );
+ fprintf( File, "# ====================================================================================================\n");
+ fprintf( File, "#%19s %5s %10s %14s", "Solver", "Lv", "Step", "Counter" );
+ fprintf( File, " %20s %24s", "Iteration", "Residual" );
+ fprintf( File, "\n" );
+ fclose( File );
+ }
+
+ FirstTime = false;
+ } // if ( FirstTime )
+
+ if ( MPI_Rank == 0 )
+ {
+ FILE *File = fopen( FileName, "a" );
+ fprintf( File, "%20s %5d %10ld %14ld %20d %24.16e\n", SolveName, lv, Step, AdvanceCounter[lv], iteration, residual );
+ fclose( File );
+ }
+
+} // FUNCTION : Hypre_Aux_Record
+
+
+
+#endif // #ifdef SUPPORT_HYPRE
diff --git a/src/Hypre/Hypre_End.cpp b/src/Hypre/Hypre_End.cpp
new file mode 100644
index 0000000000..0f4094682c
--- /dev/null
+++ b/src/Hypre/Hypre_End.cpp
@@ -0,0 +1,29 @@
+#include "GAMER.h"
+
+
+
+
+#ifdef SUPPORT_HYPRE
+//-------------------------------------------------------------------------------------------------------
+// Function : Hypre_End
+// Description : Free the resources used by the Hypre
+//
+// Note : 1. Invoked by End_GAMER()
+//
+// Parameter : None
+//
+// Return : None
+//-------------------------------------------------------------------------------------------------------
+void Hypre_End()
+{
+
+ if ( MPI_Rank == 0 ) Aux_Message( stdout, "%s ...\n", __FUNCTION__ );
+
+// Hypre_Free();
+// finalize Hypre
+ HYPRE_CHECK_FUNC( HYPRE_Finalize() );
+
+ if ( MPI_Rank == 0 ) Aux_Message( stdout, "%s ... done\n", __FUNCTION__ );
+
+} // FUNCTION : Hypre_End
+#endif // #ifdef SUPPORT_HYPRE
diff --git a/src/Hypre/Hypre_FillArrays.cpp b/src/Hypre/Hypre_FillArrays.cpp
new file mode 100644
index 0000000000..e3adb5e752
--- /dev/null
+++ b/src/Hypre/Hypre_FillArrays.cpp
@@ -0,0 +1,450 @@
+#include "GAMER.h"
+
+
+
+#ifdef SUPPORT_HYPRE
+static void Hypre_FillArrays_Poisson( const int lv, const int NExtend, const double TimeNew, const real Poi_Coeff,
+ const int SaveSg_Pot );
+#define POT_NXT_INT ( (POT_NXT-2)*2 ) // size of the array "Pot_Int_Array"
+
+
+//-------------------------------------------------------------------------------------------------------
+// Function : Hypre_FillArrays
+// Description : TBF
+//
+// Note : TBF
+//
+// Parameter : TBF
+//-------------------------------------------------------------------------------------------------------
+void Hypre_FillArrays( const Hypre_SolveType_t SolveType, const int lv, const double TimeNew, const real Poi_Coeff,
+ const int SaveSg_Pot )
+{
+
+ switch ( SolveType )
+ {
+# ifdef GRAVITY
+ case HYPRE_SOLVE_TYPE_POISSON:
+ const int NExtend = 1;
+ Hypre_FillArrays_Poisson( lv, NExtend, TimeNew, Poi_Coeff, SaveSg_Pot );
+ break;
+# endif
+ default :
+ Aux_Error( ERROR_INFO, "incorrect parameter %s = %d !!\n", "SolveType", SolveType );
+ } // switch ( SolveType )
+
+} // FUNCTION : Hypre_FillArrays
+
+
+
+#ifdef GRAVITY
+//-------------------------------------------------------------------------------------------------------
+// Function : Hypre_FillArrays_Poisson
+// Description : TBF
+//
+// Note : TBF
+//
+// Parameter : TBF
+//-------------------------------------------------------------------------------------------------------
+void Hypre_FillArrays_Poisson( const int lv, const int NExtend, const double TimeNew, const real Poi_Coeff,
+ const int SaveSg_Pot )
+{
+
+ const bool In_time = ( amr->PotSgTime[lv][SaveSg_Pot] == TimeNew || amr->PotSgTime[lv][1-SaveSg_Pot] == TimeNew );
+ const double Time_tmp = amr->PotSgTime[lv][SaveSg_Pot];
+ if ( !In_time ) amr->PotSgTime[lv][SaveSg_Pot] = TimeNew;
+
+ const bool IntPhase_No = false;
+ const bool DE_Consistency_No = false;
+ const real MinDens_No = -1.0;
+ const real MinPres_No = -1.0;
+ const real MinTemp_No = -1.0;
+ const real MinEntr_No = -1.0;
+ const int NDim = 3;
+ const int NEntries = 2*NDim + 1;
+ const int var = 0;
+ const int part = 0;
+ const double dh = amr->dh[lv];
+ const double dh2 = SQR( dh );
+ const double coeff = - Poi_Coeff * dh2;
+
+ int stencil_indices[NEntries];
+
+ for (int i=0; iNPatchComma[lv][1] / 8 ];
+ int *PID0_BC_List = new int [ amr->NPatchComma[lv][1] / 8 ];
+ int NPG_BC = 0;
+ int NPG = 0;
+
+ for (int PID0=0; PID0NPatchComma[lv][1]; PID0+=8)
+ {
+ bool atBC = false;
+
+ for (int PID=PID0; PIDpatch[0][lv][PID]->sibling[s];
+
+ if ( SibPID >= 0 ) continue;
+ if ( SibPID < SIB_OFFSET_NONPERIODIC && OPT__BC_POT == BC_POT_PERIODIC ) continue;
+
+ atBC = true;
+ break;
+ }
+ if ( atBC ) break;
+ } // for (int PID=PID0; PIDNPatchComma[lv][1] / 8; PID0+=8)
+
+ real_hypre *Matrix_Laplace, *pote, *dens, *bcBox, *bcVal;
+# ifdef GPU
+ cudaMallocManaged( &Matrix_Laplace, sizeof(real_hypre) * CUBE(PS1) * NEntries, cudaMemAttachGlobal );
+ cudaMallocManaged( &pote, sizeof(real_hypre) * CUBE(PS1), cudaMemAttachGlobal );
+ cudaMallocManaged( &dens, sizeof(real_hypre) * CUBE(PS1), cudaMemAttachGlobal );
+ cudaMallocManaged( &bcBox, sizeof(real_hypre) * SQR(PS1) * 1, cudaMemAttachGlobal );
+ cudaMallocManaged( &bcVal, sizeof(real_hypre) * SQR(PS1), cudaMemAttachGlobal );
+# else
+ Matrix_Laplace = new real_hypre [NEntries*CUBE(PS1)];
+ pote = new real_hypre [CUBE(PS1)];
+ dens = new real_hypre [CUBE(PS1)];
+ bcBox = new real_hypre [1*SQR(PS1)];
+ bcVal = new real_hypre [SQR(PS1)];
+# endif
+ real (*Pot_Array)[CUBE(PS1+2)] = NULL;
+ real (*Dens_Array)[CUBE(PS1)] = NULL;
+
+// fill common Laplace matrix
+ for (int k=0; k Pot_Int_Array
+# pragma omp for schedule( runtime )
+ for (int P=0; P<8*NPG; P++)
+ {
+ switch ( OPT__POT_INT_SCHEME )
+ {
+ /*
+ case INT_CENTRAL :
+ {
+ for (int k=1; kpatch[0][lv][PID]->EdgeL[2] + (0.5+k)*dh;
+ for (int j=0; jpatch[0][lv][PID]->EdgeL[1] + (0.5+j)*dh;
+ for (int i=0; ipatch[0][lv][PID]->EdgeL[0] + (0.5+i)*dh;
+ const int idx = IDX321( i, j, k, PS1, PS1 );
+ const int idx_i = i + (POT_NXT_INT-PS1)/2;
+ const int idx_j = j + (POT_NXT_INT-PS1)/2;
+ const int idx_k = k + (POT_NXT_INT-PS1)/2;
+
+ real_hypre Dens = (real_hypre)Dens_Array[8*PG+LocalID][idx] - RhoSubtract;
+
+// add extra mass source for gravity if required
+ if ( OPT__GRAVITY_EXTRA_MASS )
+ Dens += Poi_AddExtraMassForGravity_Ptr( x, y, z, Time[lv], lv, NULL );
+
+ dens[idx] = (real_hypre)coeff * Dens;
+ if ( HYPRE_INIT_GUESS ) pote[idx] = (real_hypre)Pot_Int_Array[8*PG+LocalID][idx_k][idx_j][idx_i];
+ else pote[idx] = (real_hypre)0.0;
+ }}} // i, j, k
+
+ HYPRE_CHECK_FUNC( HYPRE_SStructVectorSetBoxValues( Hypre_b, part, amr->patch[0][lv][PID]->cornerL, amr->patch[0][lv][PID]->cornerR, var, dens ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructVectorSetBoxValues( Hypre_x, part, amr->patch[0][lv][PID]->cornerL, amr->patch[0][lv][PID]->cornerR, var, pote ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructMatrixSetBoxValues( Hypre_A, part, amr->patch[0][lv][PID]->cornerL, amr->patch[0][lv][PID]->cornerR, var, NEntries, stencil_indices, Matrix_Laplace ) );
+ } // for (int PID=0; PIDNPatchComma[lv][1]; PID++)
+
+// prepare boundary potential
+ if ( NPG_BC > 0 )
+ {
+ Pot_Array = new real [8*NPG_BC][ CUBE(PS1+2) ];
+
+ Prepare_PatchData( lv, TimeNew, Pot_Array[0], NULL,
+ NExtend, NPG_BC, PID0_BC_List, _POTE, _NONE, OPT__GRA_INT_SCHEME, INT_NONE, UNIT_PATCH, NSIDE_06,
+ IntPhase_No, OPT__BC_FLU, OPT__BC_POT, MinDens_No, MinPres_No, MinTemp_No, MinEntr_No, DE_Consistency_No );
+ }
+
+// fill boundary condition
+ for (int PG=0; PGpatch[0][lv][PID]->cornerL[0], amr->patch[0][lv][PID]->cornerL[1], amr->patch[0][lv][PID]->cornerL[2] };
+ int cornerR[3] = { amr->patch[0][lv][PID]->cornerR[0], amr->patch[0][lv][PID]->cornerR[1], amr->patch[0][lv][PID]->cornerR[2] };
+
+// remove the connection and update A due to Dirichlet boundary condition
+ for (int s=0; s<6; s++)
+ {
+ const int SibPID = amr->patch[0][lv][PID]->sibling[s];
+
+ if ( SibPID >= 0 ) continue;
+ if ( SibPID < SIB_OFFSET_NONPERIODIC && OPT__BC_POT == BC_POT_PERIODIC ) continue;
+
+ const int d = s/2;
+ const int TDir1 = (d+1) % 3;
+ const int TDir2 = (d+2) % 3;
+ const int lr = s%2;
+
+ int stencil_indices_bc[1] = { s+1 };
+ int cornerL_bc [3] = { cornerL[0], cornerL[1], cornerL[2] };
+ int cornerR_bc [3] = { cornerR[0], cornerR[1], cornerR[2] };
+
+ if ( lr ) cornerL_bc[d] = cornerR_bc[d];
+ else cornerR_bc[d] = cornerL_bc[d];
+
+ const int Nk = cornerR_bc[2] - cornerL_bc[2] + 1;
+ const int Nj = cornerR_bc[1] - cornerL_bc[1] + 1;
+ const int Ni = cornerR_bc[0] - cornerL_bc[0] + 1;
+
+ HYPRE_CHECK_FUNC( HYPRE_SStructVectorGetBoxValues( Hypre_b, part, cornerL_bc, cornerR_bc, var, bcVal ) );
+
+ for (int k=0; kPotSgTime[lv][SaveSg_Pot] = Time_tmp;
+
+} // FUNCITON : Hypre_FillArrays_Poisson
+#endif // #ifdef GRAVITY
+
+
+
+#endif // #ifdef SUPPORT_HYPRE
diff --git a/src/Hypre/Hypre_Free.cpp b/src/Hypre/Hypre_Free.cpp
new file mode 100644
index 0000000000..ffe33f173d
--- /dev/null
+++ b/src/Hypre/Hypre_Free.cpp
@@ -0,0 +1,23 @@
+#include "GAMER.h"
+
+
+
+#ifdef SUPPORT_HYPRE
+//-------------------------------------------------------------------------------------------------------
+// Function : Hypre_Free
+// Description : Destory the Hypre arrays
+//
+// Parameter : None
+//-------------------------------------------------------------------------------------------------------
+void Hypre_Free()
+{
+
+ HYPRE_CHECK_FUNC( HYPRE_SStructGridDestroy( Hypre_grid ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructGraphDestroy( Hypre_graph ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructStencilDestroy( Hypre_stencil ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructMatrixDestroy( Hypre_A ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructVectorDestroy( Hypre_x ) );
+ HYPRE_CHECK_FUNC( HYPRE_SStructVectorDestroy( Hypre_b ) );
+
+} // FUNCITON : Hypre_Free
+#endif // #ifdef SUPPORT_HYPRE
diff --git a/src/Hypre/Hypre_Init.cpp b/src/Hypre/Hypre_Init.cpp
new file mode 100644
index 0000000000..a29970d7c1
--- /dev/null
+++ b/src/Hypre/Hypre_Init.cpp
@@ -0,0 +1,73 @@
+#include "GAMER.h"
+
+
+
+
+#ifdef SUPPORT_HYPRE
+//-------------------------------------------------------------------------------------------------------
+// Function : Hypre_Init
+// Description : Initialize Hypre
+//
+// Note : 1. Invoked by Init_GAMER()
+//
+// Parameter : None
+//
+// Return : None
+//-------------------------------------------------------------------------------------------------------
+void Hypre_Init()
+{
+
+ if ( MPI_Rank == 0 ) Aux_Message( stdout, "%s ...\n", __FUNCTION__ );
+
+// checks
+# if ( !defined SERIAL && HYPRE_HAVE_MPI != 1 )
+ Aux_Error( ERROR_INFO, "GAMER MPI (enabled) != HYPRE MPI (disable) !!\n" );
+# endif
+
+# if ( defined SERIAL && HYPRE_HAVE_MPI == 1 )
+ Aux_Error( ERROR_INFO, "GAMER MPI (disabled) != HYPRE MPI (enable) !!\n" );
+# endif
+
+# if ( defined GPU && HYPRE_USING_CUDA != 1 )
+ Aux_Error( ERROR_INFO, "GAMER GPU (enabled) != HYPRE GPU (disable) !!\n" );
+# endif
+
+# if ( !defined GPU && HYPRE_USING_CUDA == 1 )
+ Aux_Error( ERROR_INFO, "GAMER GPU (disabled) != HYPRE GPU (enable) !!\n" );
+# endif
+
+
+// warning
+# if ( defined OPENMP && HYPRE_USING_OPENMP != 1 )
+ Aux_Message( stdout, "WARNING : GAMER OPENMP (enabled) != HYPRE OPENMP (disable) !!\n" );
+# endif
+
+# if ( !defined OPENMP && HYPRE_USING_OPENMP == 1 )
+ Aux_Message( stdout, "WARNING : GAMER OPENMP (disabled) != HYPRE OPENMP (enable) !!\n" );
+# endif
+
+
+// initialize Hypre
+ HYPRE_CHECK_FUNC( HYPRE_Initialize() );
+
+# ifdef GPU
+ for (int r=0; rNPatchComma[lv][1]; PID++)
+ {
+# ifdef DEBUG_HYPRE
+ for (int d=0; d<3; d++)
+ if ( amr->patch[0][lv][PID]->cornerR[d] - amr->patch[0][lv][PID]->cornerL[d] != PS1-1 ) Aux_Error( ERROR_INFO, "cornerR[%d]-cornerL[%d] != PS1-1 !!\n", d, d );
+# endif
+
+ // TODO: Do not include the patch has son (not this level) ? If yes, need to update how to fill the array boundary
+ HYPRE_CHECK_FUNC( HYPRE_SStructGridSetExtents( Hypre_grid, part, amr->patch[0][lv][PID]->cornerL, amr->patch[0][lv][PID]->cornerR ) );
+ } // for (int PID=0; PIDNPatchComma[lv][1]; PID++)
+
+// solve one cell-centered variables
+ HYPRE_SStructVariable vartypes[] = { HYPRE_SSTRUCT_VARIABLE_CELL };
+
+ HYPRE_CHECK_FUNC( HYPRE_SStructGridSetVariables( Hypre_grid, part, NVars, vartypes ) );
+
+// set periodic
+ bool SetPeriodic = false;
+ int periodicity[3];
+ for (int d=0; d<3; d++)
+ {
+ periodicity[d] = ( Periodic[d] ) ? NX0_TOT[d]*(int)(1L<NPatchComma[lv][1]; PID++)
+ {
+ HYPRE_CHECK_FUNC( HYPRE_SStructVectorGetBoxValues( Hypre_x, part, amr->patch[0][lv][PID]->cornerL, amr->patch[0][lv][PID]->cornerR, var, pote ) );
+
+ for (int k=0; kpatch[ SaveSg_Pot ][lv][PID]->pot[k][j][i] = (real)pote[idx];
+ } // i, j, k
+ } // for (int PID=0; PIDNPatchComma[lv][1]; PID++)
+
+# ifdef GPU
+ cudaFree( pote );
+# else
+ delete [] pote;
+# endif
+
+} // FUNCTION : Hypre_UpdateArrays_Poisson
+#endif // #ifdef GRAVITY
+
+
+
+#endif // #ifdef SUPPORT_HYPRE
diff --git a/src/Init/End_GAMER.cpp b/src/Init/End_GAMER.cpp
index 58aac3d1ea..6639ebfdeb 100644
--- a/src/Init/End_GAMER.cpp
+++ b/src/Init/End_GAMER.cpp
@@ -37,6 +37,10 @@ void End_GAMER()
Grackle_End();
# endif
+# ifdef SUPPORT_HYPRE
+ Hypre_End();
+# endif
+
# if ( MODEL == HYDRO )
EoS_End();
# endif
diff --git a/src/Init/Init_ByRestart_HDF5.cpp b/src/Init/Init_ByRestart_HDF5.cpp
index 88c3af4c93..a2476658a3 100644
--- a/src/Init/Init_ByRestart_HDF5.cpp
+++ b/src/Init/Init_ByRestart_HDF5.cpp
@@ -1607,6 +1607,7 @@ void Check_Makefile( const char *FileName, const int FormatVersion )
LoadField( "LibYTJupyter", &RS.LibYTJupyter, SID, TID, NonFatal, &RT.LibYTJupyter, 1, NonFatal );
# endif
LoadField( "SupportGrackle", &RS.SupportGrackle, SID, TID, NonFatal, &RT.SupportGrackle, 1, NonFatal );
+ LoadField( "SupportHypre", &RS.SupportHypre, SID, TID, NonFatal, &RT.SupportHypre, 1, NonFatal );
LoadField( "RandomNumber", &RS.RandomNumber, SID, TID, NonFatal, &RT.RandomNumber, 1, NonFatal );
LoadField( "NLevel", &RS.NLevel, SID, TID, NonFatal, &RT.NLevel, 1, NonFatal );
@@ -2379,6 +2380,16 @@ void Check_InputPara( const char *FileName, const int FormatVersion )
LoadField( "Yt_JupyterUseConnectionFile", &RS.Yt_JupyterUseConnectionFile, SID, TID, NonFatal, &RT.Yt_JupyterUseConnectionFile, 1, NonFatal );
# endif
+// Hypre
+# ifdef SUPPORT_HYPRE
+ LoadField( "Hypre_Solver", &RS.Hypre_Solver, SID, TID, NonFatal, &RT.Hypre_Solver, 1, NonFatal );
+ LoadField( "Hypre_PrintLevel", &RS.Hypre_PrintLevel, SID, TID, NonFatal, &RT.Hypre_PrintLevel, 1, NonFatal );
+ LoadField( "Hypre_EnableLogging", &RS.Hypre_EnableLogging, SID, TID, NonFatal, &RT.Hypre_EnableLogging, 1, NonFatal );
+ LoadField( "Hypre_MaxIter", &RS.Hypre_MaxIter, SID, TID, NonFatal, &RT.Hypre_MaxIter, 1, NonFatal );
+ LoadField( "Hypre_RelTol", &RS.Hypre_RelTol, SID, TID, NonFatal, &RT.Hypre_RelTol, 1, NonFatal );
+ LoadField( "Hypre_AbsTol", &RS.Hypre_AbsTol, SID, TID, NonFatal, &RT.Hypre_AbsTol, 1, NonFatal );
+# endif
+
// miscellaneous
LoadField( "Opt__Verbose", &RS.Opt__Verbose, SID, TID, NonFatal, &RT.Opt__Verbose, 1, NonFatal );
LoadField( "Opt__TimingBarrier", &RS.Opt__TimingBarrier, SID, TID, NonFatal, &RT.Opt__TimingBarrier, 1, NonFatal );
diff --git a/src/Init/Init_FFTW.cpp b/src/Init/Init_FFTW.cpp
index acc1b78a1d..9eb15fe703 100644
--- a/src/Init/Init_FFTW.cpp
+++ b/src/Init/Init_FFTW.cpp
@@ -4,15 +4,30 @@
static int ZIndex2Rank( const int IndexZ, const int *List_z_start, const int TRank_Guess );
-root_fftw::real_plan_nd FFTW_Plan_PS; // PS : plan for calculating the power spectrum
+// PS : plan for calculating the power spectrum
+root_fftw::real_plan_nd FFTW_Plan_PS;
+static void Init_FFTW_PowerSpectrum( const int StartupFlag );
+static void End_FFTW_PowerSpectrum();
+
+// Poi : plan for the self-gravity Poisson solver
#ifdef GRAVITY
-root_fftw::real_plan_nd FFTW_Plan_Poi, FFTW_Plan_Poi_Inv; // Poi : plan for the self-gravity Poisson solver
+root_fftw::real_plan_nd FFTW_Plan_Poi[NLEVEL], FFTW_Plan_Poi_Inv[NLEVEL];
+static void Init_FFTW_Poisson( const int StartupFlag, const int lv );
+static void End_FFTW_Poisson();
#endif // #ifdef GRAVITY
+
+// Psi : plan for the ELBDM spectral solver
#if ( MODEL == ELBDM )
-root_fftw::complex_plan_nd FFTW_Plan_Psi, FFTW_Plan_Psi_Inv; // Psi : plan for the ELBDM spectral solver
+root_fftw::complex_plan_nd FFTW_Plan_Psi, FFTW_Plan_Psi_Inv;
+static void Init_FFTW_ELBDMSpectral( const int StartupFlag );
+static void End_FFTW_ELBDMSpectral();
+
+// ExtPsi : plan for the Gram Fourier extension solver
#if ( WAVE_SCHEME == WAVE_GRAMFE )
-gramfe_fftw::complex_plan_1d FFTW_Plan_ExtPsi, FFTW_Plan_ExtPsi_Inv; // ExtPsi : plan for the Gram Fourier extension solver
-#endif // #if (WAVE_SCHEME == WAVE_GRAMFE)
+gramfe_fftw::complex_plan_1d FFTW_Plan_ExtPsi, FFTW_Plan_ExtPsi_Inv;
+static void Init_FFTW_GramFE( const int StartupFlag );
+static void End_FFTW_GramFE();
+#endif // #if ( WAVE_SCHEME == WAVE_GRAMFE )
#endif // #if ( MODEL == ELBDM )
@@ -26,8 +41,11 @@ gramfe_fftw::complex_plan_1d FFTW_Plan_ExtPsi, FFTW_Plan_ExtPsi_Inv; // ExtPsi
//
// Return : length of array that is large enough to store FFT input and output
//-------------------------------------------------------------------------------------------------------
-int ComputePaddedTotalSize(int* size) {
+int ComputePaddedTotalSize( int* size )
+{
+
return 2*(size[0]/2+1)*size[1]*size[2];
+
} // FUNCTION : ComputePaddedTotalSize
@@ -40,8 +58,11 @@ int ComputePaddedTotalSize(int* size) {
//
// Return : length of array that is large enough to store FFT input and output
//-------------------------------------------------------------------------------------------------------
-int ComputeTotalSize(int* size) {
+int ComputeTotalSize( int* size )
+{
+
return size[0]*size[1]*size[2];
+
} // FUNCTION : ComputeTotalSize
@@ -49,181 +70,101 @@ int ComputeTotalSize(int* size) {
//-------------------------------------------------------------------------------------------------------
// Function : Init_FFTW
// Description : Create the FFTW plans
+//
+// Parameter : lv : Target level
//-------------------------------------------------------------------------------------------------------
-void Init_FFTW()
+void Init_FFTW( const int lv )
{
if ( MPI_Rank == 0 ) Aux_Message( stdout, "%s ... ", __FUNCTION__ );
+ if ( FFTW_Inited[lv] ) return;
-# if ( SUPPORT_FFTW == FFTW3 )
- FFTW3_Double_OMP_Enabled = false;
- FFTW3_Single_OMP_Enabled = false;
-
-
-# ifdef OPENMP
-
-# ifndef SERIAL
-// check the level of MPI thread support
- int MPI_Thread_Status;
- MPI_Query_thread( &MPI_Thread_Status );
-
-// enable multithreading if possible
- FFTW3_Double_OMP_Enabled = MPI_Thread_Status >= MPI_THREAD_FUNNELED;
- FFTW3_Single_OMP_Enabled = FFTW3_Double_OMP_Enabled;
-# else // # ifndef SERIAL
-
-// always enable multithreading in serial mode with openmp
- FFTW3_Double_OMP_Enabled = true;
- FFTW3_Single_OMP_Enabled = true;
-# endif // # ifndef SERIAL ... # else
-
-// initialise fftw multithreading
- if (FFTW3_Double_OMP_Enabled) {
- FFTW3_Double_OMP_Enabled = fftw_init_threads();
- if ( !FFTW3_Double_OMP_Enabled ) Aux_Error( ERROR_INFO, "fftw_init_threads() failed !!\n" );
- }
- if (FFTW3_Single_OMP_Enabled) {
- FFTW3_Single_OMP_Enabled = fftwf_init_threads();
- if ( !FFTW3_Single_OMP_Enabled ) Aux_Error( ERROR_INFO, "fftwf_init_threads() failed !!\n" );
- }
-# endif // # ifdef OPENMP
-
-// initialise fftw mpi support
-# ifndef SERIAL
- fftw_mpi_init();
- fftwf_mpi_init();
-# endif // # ifndef SERIAL
-
-// tell all subsequent fftw3 planners to use OMP_NTHREAD threads
-# ifdef OPENMP
- if (FFTW3_Double_OMP_Enabled) fftw_plan_with_nthreads (OMP_NTHREAD);
- if (FFTW3_Single_OMP_Enabled) fftwf_plan_with_nthreads(OMP_NTHREAD);
-# endif // # ifdef OPENMP
-# endif // # if ( SUPPORT_FFTW == FFTW3 )
-
-
-
-// determine the FFT size for the power spectrum
- int PS_FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
-
-// determine the FFT size for the self-gravity solver
-# ifdef GRAVITY
- int Gravity_FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
-
-// the zero-padding method is adopted for the isolated BC.
- if ( OPT__BC_POT == BC_POT_ISOLATED )
- for (int d=0; d<3; d++) Gravity_FFT_Size[d] *= 2;
-
-// check
- if ( MPI_Rank == 0 )
- for (int d=0; d<3; d++)
- {
- if ( Gravity_FFT_Size[d] <= 0 ) Aux_Error( ERROR_INFO, "Gravity_FFT_Size[%d] = %d < 0 !!\n", d, Gravity_FFT_Size[d] );
- }
-# endif // # ifdef GRAVITY
-
-// determine the FFT size for the base-level FFT wave solver
-# if ( MODEL == ELBDM )
- int Psi_FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
-# if ( defined(SERIAL) || SUPPORT_FFTW == FFTW3 )
- int InvPsi_FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
-# else // # ifdef SERIAL || FFTW3
-// Note that the dimensions of the inverse transform in FFTW2,
-// which are given by the dimensions of the output of the forward transform,
-// are Ny*Nz*Nx because we are using "FFTW_TRANSPOSED_ORDER" in fftwnd_mpi().
- int InvPsi_FFT_Size[3] = { NX0_TOT[0], NX0_TOT[2], NX0_TOT[1] };
-# endif // # ifdef SERIAL || FFTW3 ... # else
-
-# if ( WAVE_SCHEME == WAVE_GRAMFE )
- int ExtPsi_FFT_Size = GRAMFE_FLU_NXT;
-# endif // # if ( WAVE_SCHEME == WAVE_GRAMFE )
-# endif // # if ( MODEL == ELBDM )
- real* PS = NULL;
- real* RhoK = NULL;
- real* PsiK = NULL;
-# if ( WAVE_SCHEME == WAVE_GRAMFE )
- gramfe_fftw::fft_complex* ExtPsiK = NULL;
-# endif // # if ( WAVE_SCHEME == WAVE_GRAMFE )
+ static bool FirstTime = true;
// determine how to initialise fftw plans
int StartupFlag;
switch ( OPT__FFTW_STARTUP )
{
- case FFTW_STARTUP_ESTIMATE: StartupFlag = FFTW_ESTIMATE; break;
- case FFTW_STARTUP_MEASURE: StartupFlag = FFTW_MEASURE; break;
+ case FFTW_STARTUP_ESTIMATE: StartupFlag = FFTW_ESTIMATE; break;
+ case FFTW_STARTUP_MEASURE: StartupFlag = FFTW_MEASURE; break;
# if ( SUPPORT_FFTW == FFTW3 )
- case FFTW_STARTUP_PATIENT: StartupFlag = FFTW_PATIENT; break;
+ case FFTW_STARTUP_PATIENT: StartupFlag = FFTW_PATIENT; break;
# endif // # if ( SUPPORT_FFTW == FFTW3 )
- default: Aux_Error( ERROR_INFO, "unrecognised FFTW startup option %d !!\n", OPT__FFTW_STARTUP );
+ default: Aux_Error( ERROR_INFO, "unrecognised FFTW startup option %d !!\n", OPT__FFTW_STARTUP );
} // switch ( OPT__FFTW_STARTUP )
-// allocate memory for arrays in fftw3
-# if ( SUPPORT_FFTW == FFTW3 )
- PS = (real*) root_fftw::fft_malloc(ComputePaddedTotalSize(PS_FFT_Size ) * sizeof(real));
-# ifdef GRAVITY
- RhoK = (real*) root_fftw::fft_malloc(ComputePaddedTotalSize(Gravity_FFT_Size) * sizeof(real));
-# endif // # ifdef GRAVITY
-# if ( MODEL == ELBDM )
- PsiK = (real*) root_fftw::fft_malloc( ComputeTotalSize ( Psi_FFT_Size ) * sizeof(real) * 2 ); // 2 * real for size of complex number
-# endif // # if ( MODEL == ELBDM )
-# if ( WAVE_SCHEME == WAVE_GRAMFE )
- ExtPsiK = (gramfe_fftw::fft_complex*) gramfe_fftw::fft_malloc( ExtPsi_FFT_Size * sizeof(gramfe_fftw::fft_complex) );
-# endif // # if ( WAVE_SCHEME == WAVE_GRAMFE )
-# endif // # if ( SUPPORT_FFTW == FFTW3 )
+ if ( FirstTime )
+ {
+# if ( SUPPORT_FFTW == FFTW3 )
+ FFTW3_Double_OMP_Enabled = false;
+ FFTW3_Single_OMP_Enabled = false;
-// create plans for power spectrum and the self-gravity solver
- FFTW_Plan_PS = root_fftw_create_3d_r2c_plan(PS_FFT_Size, PS, StartupFlag);
-# ifdef GRAVITY
- FFTW_Plan_Poi = root_fftw_create_3d_r2c_plan(Gravity_FFT_Size, RhoK, StartupFlag);
- FFTW_Plan_Poi_Inv = root_fftw_create_3d_c2r_plan(Gravity_FFT_Size, RhoK, StartupFlag);
-# endif // # ifdef GRAVITY
-# if ( MODEL == ELBDM )
- FFTW_Plan_Psi = root_fftw_create_3d_forward_c2c_plan ( Psi_FFT_Size, PsiK, StartupFlag );
- FFTW_Plan_Psi_Inv = root_fftw_create_3d_backward_c2c_plan( InvPsi_FFT_Size, PsiK, StartupFlag );
+# ifdef OPENMP
-# if ( WAVE_SCHEME == WAVE_GRAMFE )
+# ifndef SERIAL
+// check the level of MPI thread support
+ int MPI_Thread_Status;
+ MPI_Query_thread( &MPI_Thread_Status );
-// the Gram-Fourier extension planners only use one thread because OMP parallelisation evolves different patches parallely
-// From the FFTW3 documentation: https://www.fftw.org/fftw3_doc/Usage-of-Multi_002dthreaded-FFTW.html
-// "You can call fftw_plan_with_nthreads, create some plans,
-// call fftw_plan_with_nthreads again with a different argument, and create some more plans for a new number of threads."
+// enable multithreading if possible
+ FFTW3_Double_OMP_Enabled = MPI_Thread_Status >= MPI_THREAD_FUNNELED;
+ FFTW3_Single_OMP_Enabled = FFTW3_Double_OMP_Enabled;
+# else // # ifndef SERIAL
-# if ( defined(SUPPORT_FFTW3) && defined(OPENMP) )
- if (FFTW3_Double_OMP_Enabled) fftw_plan_with_nthreads(1);
- if (FFTW3_Single_OMP_Enabled) fftwf_plan_with_nthreads(1);
-# endif // # if ( defined(SUPPORT_FFTW3) && defined(OPENMP) )
+// always enable multithreading in serial mode with openmp
+ FFTW3_Double_OMP_Enabled = true;
+ FFTW3_Single_OMP_Enabled = true;
+# endif // # ifndef SERIAL ... # else
- FFTW_Plan_ExtPsi = gramfe_fftw_create_1d_forward_c2c_plan ( ExtPsi_FFT_Size, ExtPsiK, StartupFlag );
- FFTW_Plan_ExtPsi_Inv = gramfe_fftw_create_1d_backward_c2c_plan( ExtPsi_FFT_Size, ExtPsiK, StartupFlag );
+// initialise fftw multithreading
+ if ( FFTW3_Double_OMP_Enabled )
+ {
+ FFTW3_Double_OMP_Enabled = fftw_init_threads();
+ if ( !FFTW3_Double_OMP_Enabled ) Aux_Error( ERROR_INFO, "fftw_init_threads() failed !!\n" );
+ }
+ if ( FFTW3_Single_OMP_Enabled )
+ {
+ FFTW3_Single_OMP_Enabled = fftwf_init_threads();
+ if ( !FFTW3_Single_OMP_Enabled ) Aux_Error( ERROR_INFO, "fftwf_init_threads() failed !!\n" );
+ }
+# endif // # ifdef OPENMP
+
+// initialise fftw mpi support
+# ifndef SERIAL
+ fftw_mpi_init();
+ fftwf_mpi_init();
+# endif // # ifndef SERIAL
+
+// tell all subsequent fftw3 planners to use OMP_NTHREAD threads
+# ifdef OPENMP
+ if ( FFTW3_Double_OMP_Enabled ) fftw_plan_with_nthreads ( OMP_NTHREAD );
+ if ( FFTW3_Single_OMP_Enabled ) fftwf_plan_with_nthreads( OMP_NTHREAD );
+# endif // # ifdef OPENMP
+# endif // # if ( SUPPORT_FFTW == FFTW3 )
-// restore regular settings
-# if ( defined(SUPPORT_FFTW3) && defined(OPENMP) )
- if (FFTW3_Double_OMP_Enabled) fftw_plan_with_nthreads(OMP_NTHREAD);
- if (FFTW3_Single_OMP_Enabled) fftwf_plan_with_nthreads(OMP_NTHREAD);
-# endif // # if ( defined(SUPPORT_FFTW3) && defined(OPENMP) )
-# endif // # if ( WAVE_SCHEME == WAVE_GRAMFE )
+ FirstTime = false;
+ } // if ( FirstTime )
-# endif // # if ( MODEL == ELBDM )
+ if ( lv == 0 ) Init_FFTW_PowerSpectrum( StartupFlag );
-// free memory for arrays in fftw3
-# if ( SUPPORT_FFTW == FFTW3 )
- root_fftw::fft_free(PS);
# ifdef GRAVITY
- root_fftw::fft_free(RhoK);
-# endif // # ifdef GRAVITY
+ Init_FFTW_Poisson( StartupFlag, lv );
+# endif
+
# if ( MODEL == ELBDM )
- root_fftw::fft_free( PsiK );
+ if ( lv == 0 ) Init_FFTW_ELBDMSpectral( StartupFlag );
+
# if ( WAVE_SCHEME == WAVE_GRAMFE )
- gramfe_fftw::fft_free( ExtPsiK );
-# endif // # if ( WAVE_SCHEME == WAVE_GRAMFE )
-# endif // # if ( MODEL == ELBDM )
-# endif // # if ( SUPPORT_FFTW == FFTW3 )
+ if ( lv == 0 ) Init_FFTW_GramFE( StartupFlag );
+# endif // #if ( WAVE_SCHEME == WAVE_GRAMFE )
+# endif // #if ( MODEL == ELBDM )
+ FFTW_Inited[lv] = true;
if ( MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
@@ -240,24 +181,22 @@ void End_FFTW()
if ( MPI_Rank == 0 ) Aux_Message( stdout, "%s ... ", __FUNCTION__ );
- root_fftw::destroy_real_plan_nd ( FFTW_Plan_PS );
-# ifdef GRAVITY
- root_fftw::destroy_real_plan_nd ( FFTW_Plan_Poi );
- root_fftw::destroy_real_plan_nd ( FFTW_Plan_Poi_Inv );
-# endif // # ifdef GRAVITY
+ End_FFTW_PowerSpectrum();
+# ifdef GRAVITY
+ End_FFTW_Poisson();
+# endif
# if ( MODEL == ELBDM )
- root_fftw::destroy_complex_plan_nd ( FFTW_Plan_Psi );
- root_fftw::destroy_complex_plan_nd ( FFTW_Plan_Psi_Inv );
+ End_FFTW_ELBDMSpectral();
# if ( WAVE_SCHEME == WAVE_GRAMFE )
- gramfe_fftw::destroy_complex_plan_1d ( FFTW_Plan_ExtPsi );
- gramfe_fftw::destroy_complex_plan_1d ( FFTW_Plan_ExtPsi_Inv );
+ End_FFTW_GramFE();
# endif // # if ( WAVE_SCHEME == WAVE_GRAMFE )
# endif // #if ( MODEL == ELBDM )
+
# if ( SUPPORT_FFTW == FFTW3 )
# ifdef OPENMP
if ( FFTW3_Double_OMP_Enabled ) fftw_cleanup_threads();
@@ -281,6 +220,253 @@ void End_FFTW()
+//-------------------------------------------------------------------------------------------------------
+// Function : Init_FFTW_PowerSpectrum
+// Description : Create the FFTW plan for the power spectrum
+//
+// Parameter : StartupFlag : Initialize FFTW plan method
+//
+// Return : none
+//-------------------------------------------------------------------------------------------------------
+void Init_FFTW_PowerSpectrum( const int StartupFlag )
+{
+
+// determine the FFT size
+ int PS_FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
+
+ real* PS = NULL;
+
+// allocate memory for arrays in fftw3
+# if ( SUPPORT_FFTW == FFTW3 )
+ PS = (real*)root_fftw::fft_malloc( ComputePaddedTotalSize( PS_FFT_Size ) * sizeof(real) );
+# endif
+
+// create plans
+ FFTW_Plan_PS = root_fftw_create_3d_r2c_plan( PS_FFT_Size, PS, StartupFlag );
+
+// free memory for arrays in fftw3
+# if ( SUPPORT_FFTW == FFTW3 )
+ root_fftw::fft_free( PS );
+# endif
+
+} // FUNCITON : Init_FFTW_PowerSpectrum
+
+
+
+//-------------------------------------------------------------------------------------------------------
+// Function : End_FFTW_PowerSpectrum
+// Description : Delete the FFTW plan for the power spectrum
+//
+// Return : none
+//-------------------------------------------------------------------------------------------------------
+void End_FFTW_PowerSpectrum()
+{
+
+ root_fftw::destroy_real_plan_nd( FFTW_Plan_PS );
+
+} // FUNCITON : End_FFTW_PowerSpectrum
+
+
+
+#ifdef GRAVITY
+//-------------------------------------------------------------------------------------------------------
+// Function : Init_FFTW_Poisson
+// Description : Create the FFTW plan for the self-gravity Poisson solver
+//
+// Parameter : StartupFlag : Initialize FFTW plan method
+// lv : Target level
+//
+// Return : none
+//-------------------------------------------------------------------------------------------------------
+void Init_FFTW_Poisson( const int StartupFlag, const int lv )
+{
+
+// determine the FFT size
+ int Gravity_FFT_Size[3] = { NX0_TOT[0]*(1L< "You can call fftw_plan_with_nthreads, create some plans, call fftw_plan_with_nthreads
+// again with a different argument, and create some more plans for a new number of threads."
+//
+// Parameter : StartupFlag : Initialize FFTW plan method
+//
+// Return : none
+//-------------------------------------------------------------------------------------------------------
+void Init_FFTW_GramFE( const int StartupFlag )
+{
+
+// determine the FFT size
+ int ExtPsi_FFT_Size = GRAMFE_FLU_NXT;
+
+ gramfe_fftw::fft_complex* ExtPsiK = NULL;
+
+// allocate memory for arrays in fftw3
+# if ( SUPPORT_FFTW == FFTW3 )
+ ExtPsiK = (gramfe_fftw::fft_complex*)gramfe_fftw::fft_malloc( ExtPsi_FFT_Size * sizeof(gramfe_fftw::fft_complex) );
+# endif
+
+// See note 1.
+# if ( defined( SUPPORT_FFTW3 ) && defined( OPENMP ) )
+ if ( FFTW3_Double_OMP_Enabled ) fftw_plan_with_nthreads( 1 );
+ if ( FFTW3_Single_OMP_Enabled ) fftwf_plan_with_nthreads( 1 );
+# endif
+
+// create plans
+ FFTW_Plan_ExtPsi = gramfe_fftw_create_1d_forward_c2c_plan ( ExtPsi_FFT_Size, ExtPsiK, StartupFlag );
+ FFTW_Plan_ExtPsi_Inv = gramfe_fftw_create_1d_backward_c2c_plan( ExtPsi_FFT_Size, ExtPsiK, StartupFlag );
+
+// restore regular settings
+# if ( defined( SUPPORT_FFTW3 ) && defined( OPENMP ) )
+ if ( FFTW3_Double_OMP_Enabled ) fftw_plan_with_nthreads( OMP_NTHREAD );
+ if ( FFTW3_Single_OMP_Enabled ) fftwf_plan_with_nthreads( OMP_NTHREAD );
+# endif
+
+// free memory for arrays in fftw3
+# if ( SUPPORT_FFTW == FFTW3 )
+ gramfe_fftw::fft_free( ExtPsiK );
+# endif
+
+} // FUNCTION : Init_FFTW_Poisson
+
+
+
+//-------------------------------------------------------------------------------------------------------
+// Function : End_FFTW_GramFE
+// Description : Delete the FFTW plan for the Gram Fourier extension solver
+//
+// Return : none
+//-------------------------------------------------------------------------------------------------------
+void End_FFTW_GramFE()
+{
+
+ gramfe_fftw::destroy_complex_plan_1d( FFTW_Plan_ExtPsi );
+ gramfe_fftw::destroy_complex_plan_1d( FFTW_Plan_ExtPsi_Inv );
+
+} // FUNCITON : End_FFTW_Poisson
+#endif // #if ( WAVE_SCHEME == WAVE_GRAMFE )
+#endif // #if ( MODEL == ELBDM )
+
+
+
//-------------------------------------------------------------------------------------------------------
// Function : Patch2Slab
// Description : Patch-based data --> slab domain decomposition
@@ -306,13 +492,16 @@ void End_FFTW()
// InPlacePad : Whether or not to pad the array size for in-place real-to-complex FFT
// ForPoisson : Preparing the density field for the Poisson solver
// AddExtraMass : Adding an extra density field for computing gravitational potential (only works with ForPoisson)
+// lv : Target level
//-------------------------------------------------------------------------------------------------------
void Patch2Slab( real *VarS, real *SendBuf_Var, real *RecvBuf_Var, long *SendBuf_SIdx, long *RecvBuf_SIdx,
int **List_PID, int **List_k, long *List_NSend_Var, long *List_NRecv_Var,
const int *List_z_start, const int local_nz, const int FFT_Size[], const int NRecvSlice,
- const double PrepTime, const long TVar, const bool InPlacePad, const bool ForPoisson, const bool AddExtraMass )
+ const double PrepTime, const long TVar, const bool InPlacePad, const bool ForPoisson,
+ const bool AddExtraMass, const int lv )
{
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
// check
// check only single field
if ( TVar == 0 || TVar & (TVar-1) )
@@ -333,13 +522,13 @@ void Patch2Slab( real *VarS, real *SendBuf_Var, real *RecvBuf_Var, long *SendBuf
local_nz, List_z_start[MPI_Rank+1] - List_z_start[MPI_Rank] );
# endif // GAMER_DEBUG
-
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
const int SSize[2] = { ( InPlacePad ? 2*(FFT_Size[0]/2+1) : FFT_Size[0] ), FFT_Size[1] }; // padded slab size in the x and y directions
const int PSSize = PS1*PS1; // patch slice size
-// const int MemUnit = amr->NPatchComma[0][1]*PS1/MPI_NRank; // set arbitrarily
- const int MemUnit = amr->NPatchComma[0][1]*PS1; // set arbitrarily
+// const int MemUnit = amr->NPatchComma[lv][1]*PS1/MPI_NRank; // set arbitrarily
+ const int MemUnit = amr->NPatchComma[lv][1]*PS1; // set arbitrarily
const int AveNz = FFT_Size[2]/MPI_NRank + ( ( FFT_Size[2]%MPI_NRank == 0 ) ? 0 : 1 ); // average slab thickness
- const int Scale0 = amr->scale[0];
+ const int Scale = amr->scale[lv];
int Cr[3]; // corner coordinates of each patch normalized to the base-level grid size
int BPos_z; // z coordinate of each patch slice in the simulation box
@@ -353,6 +542,7 @@ void Patch2Slab( real *VarS, real *SendBuf_Var, real *RecvBuf_Var, long *SendBuf
int TRank, TRank_Guess, MemSize[MPI_NRank], idx;
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
// 1. set memory allocation unit
for (int r=0; rNPatchComma[0][1]; PID0+=8)
+ for (int PID0=0; PID0NPatchComma[lv][1]; PID0+=8)
{
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
// even with NSIDE_00 and GhostSize=0, we still need OPT__BC_FLU to determine whether periodic BC is adopted
// for depositing particle mass onto grids.
// also note that we do not check minimum density here since no ghost zones are required
- Prepare_PatchData( 0, PrepTime, VarPatch[0][0][0], NULL, GhostSize, NPG, &PID0, TVar, _NONE,
+ Prepare_PatchData( lv, PrepTime, VarPatch[0][0][0], NULL, GhostSize, NPG, &PID0, TVar, _NONE,
IntScheme, INT_NONE, UNIT_PATCH, NSide_None, IntPhase_No, OPT__BC_FLU, PotBC_None,
MinDens_No, MinPres_No, MinTemp_No, MinEntr_No, DE_Consistency_No );
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
# ifdef GRAVITY
// add extra mass source for gravity if required
if ( ForPoisson && AddExtraMass )
{
- const double dh = amr->dh[0];
+ const double dh = amr->dh[lv];
for (int PID=PID0, LocalID=0; PIDpatch[0][0][PID]->EdgeL[0] + 0.5*dh;
- const double y0 = amr->patch[0][0][PID]->EdgeL[1] + 0.5*dh;
- const double z0 = amr->patch[0][0][PID]->EdgeL[2] + 0.5*dh;
+ const double x0 = amr->patch[0][lv][PID]->EdgeL[0] + 0.5*dh;
+ const double y0 = amr->patch[0][lv][PID]->EdgeL[1] + 0.5*dh;
+ const double z0 = amr->patch[0][lv][PID]->EdgeL[2] + 0.5*dh;
double x, y, z;
for (int k=0; kpatch[0][0][PID]->corner[d] / Scale0;
+ for (int d=0; d<3; d++) Cr[d] = amr->patch[0][lv][PID]->corner[d] / Scale;
for (int k=0; kNPatchComma[0][1]; PID0+=8)
delete [] VarPatch;
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
// 3. prepare the send buffer
int Send_Disp_SIdx[MPI_NRank], Recv_Disp_SIdx[MPI_NRank];
@@ -503,12 +699,13 @@ void Patch2Slab( real *VarS, real *SendBuf_Var, real *RecvBuf_Var, long *SendBuf
Recv_Disp_Var [r] = Recv_Disp_Var [r-1] + List_NRecv_Var [r-1];
}
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
// check
# ifdef GAMER_DEBUG
const long NSend_Total = Send_Disp_Var[MPI_NRank-1] + List_NSend_Var[MPI_NRank-1];
const long NRecv_Total = Recv_Disp_Var[MPI_NRank-1] + List_NRecv_Var[MPI_NRank-1];
- const long NSend_Expect = (long)amr->NPatchComma[0][1]*(long)CUBE(PS1);
- const long NRecv_Expect = (long)NX0_TOT[0]*(long)NX0_TOT[1]*(long)NRecvSlice;
+ const long NSend_Expect = (long)amr->NPatchComma[lv][1]*(long)CUBE(PS1);
+ const long NRecv_Expect = (long)NX0_TOT[0]*(long)(1<patch[SaveSg][0][PID]->fluid[TVarIdx][k], RecvPtr, PSSize*sizeof(real) );
+ memcpy( amr->patch[SaveSg][lv][PID]->fluid[TVarIdx][k], RecvPtr, PSSize*sizeof(real) );
# ifdef GRAVITY
else if ( TVarIdx == NCOMP_TOTAL+NDERIVE ) // TVar == _POTE
- memcpy( amr->patch[SaveSg][0][PID]->pot[k], RecvPtr, PSSize*sizeof(real) );
+ memcpy( amr->patch[SaveSg][lv][PID]->pot[k], RecvPtr, PSSize*sizeof(real) );
# endif
else
Aux_Error( ERROR_INFO, "incorrect target variable index %s = %d !!\n", "TVarIdx", TVarIdx );
diff --git a/src/Init/Init_GAMER.cpp b/src/Init/Init_GAMER.cpp
index fc0bfa840a..5b0308852a 100644
--- a/src/Init/Init_GAMER.cpp
+++ b/src/Init/Init_GAMER.cpp
@@ -88,9 +88,15 @@ void Init_GAMER( int *argc, char ***argv )
# endif
+// initialize Hypre
+# ifdef SUPPORT_HYPRE
+ Hypre_Init();
+# endif
+
+
# ifdef SUPPORT_FFTW
// initialize FFTW
- Init_FFTW();
+ Init_FFTW( 0 );
# endif
@@ -266,7 +272,7 @@ void Init_GAMER( int *argc, char ***argv )
{
# ifdef SUPPORT_FFTW
// initialize the k-space Green's function for the isolated BC.
- if ( OPT__SELF_GRAVITY && OPT__BC_POT == BC_POT_ISOLATED ) Init_GreenFuncK();
+ if ( OPT__SELF_GRAVITY && OPT__BC_POT == BC_POT_ISOLATED ) Init_GreenFuncK( 0 );
# endif
diff --git a/src/Init/Init_Load_Parameter.cpp b/src/Init/Init_Load_Parameter.cpp
index b19f86e3d4..838f79a7b0 100644
--- a/src/Init/Init_Load_Parameter.cpp
+++ b/src/Init/Init_Load_Parameter.cpp
@@ -566,6 +566,21 @@ void Init_Load_Parameter()
# endif
# endif
+
+// Hypre
+# ifdef SUPPORT_HYPRE
+ ReadPara->Add( "HYPRE_SOLVER", &HYPRE_SOLVER, 1, 1, 2 );
+ ReadPara->Add( "HYPRE_INIT_GUESS", &HYPRE_INIT_GUESS, true, Useless_bool, Useless_bool );
+ ReadPara->Add( "HYPRE_PRINT_LEVEL", &HYPRE_PRINT_LEVEL, 2, 0, NoMax_int );
+ ReadPara->Add( "HYPRE_ENABLE_LOGGING", &HYPRE_ENABLE_LOGGING, 1, 0, 1 );
+ ReadPara->Add( "HYPRE_MAX_ITER", &HYPRE_MAX_ITER, -1, NoMin_int, NoMax_int );
+ ReadPara->Add( "HYPRE_NPRE_RELAX", &HYPRE_NPRE_RELAX, 1, NoMin_int, NoMax_int );
+ ReadPara->Add( "HYPRE_NPOST_RELAX", &HYPRE_NPOST_RELAX, 1, NoMin_int, NoMax_int );
+ ReadPara->Add( "HYPRE_REL_TOL", &HYPRE_REL_TOL, -1.0, NoMin_double, NoMax_double );
+ ReadPara->Add( "HYPRE_ABS_TOL", &HYPRE_ABS_TOL, -1.0, NoMin_double, NoMax_double );
+# endif
+
+
// miscellaneous
ReadPara->Add( "OPT__VERBOSE", &OPT__VERBOSE, false, Useless_bool, Useless_bool );
// do not check OPT__TIMING_BARRIER since it depends on other options
diff --git a/src/Init/Init_ResetParameter.cpp b/src/Init/Init_ResetParameter.cpp
index 4ecd59889f..331dcba80a 100644
--- a/src/Init/Init_ResetParameter.cpp
+++ b/src/Init/Init_ResetParameter.cpp
@@ -1292,6 +1292,39 @@ void Init_ResetParameter()
# endif
+// Hypre
+# ifdef SUPPORT_HYPRE
+ if ( HYPRE_MAX_ITER < 0 )
+ {
+ HYPRE_MAX_ITER = 50;
+
+ PRINT_RESET_PARA( HYPRE_MAX_ITER, FORMAT_INT, "" );
+ }
+
+ if ( HYPRE_REL_TOL < 0.0 )
+ {
+# ifdef FLOAT8
+ HYPRE_REL_TOL = 1.e-14;
+# else
+ HYPRE_REL_TOL = 1.e-6;
+# endif
+
+ PRINT_RESET_PARA( HYPRE_REL_TOL, FORMAT_REAL, "" );
+ }
+
+ if ( HYPRE_ABS_TOL < 0.0 )
+ {
+# ifdef FLOAT8
+ HYPRE_ABS_TOL = 1.e-14;
+# else
+ HYPRE_ABS_TOL = 1.e-6;
+# endif
+
+ PRINT_RESET_PARA( HYPRE_ABS_TOL, FORMAT_REAL, "" );
+ }
+# endif // #ifdef SUPPORT_HYPRE
+
+
// must set OPT__FFTW_STARTUP = FFTW_STARTUP_ESTIMATE for BITWISE_REPRODUCIBILITY
// --> even when disabling BITWISE_REPRODUCIBILITY, we still use FFTW_STARTUP_ESTIMATE
// by default since otherwise the FFT results can vary in each run on the level
diff --git a/src/Interpolation/Interpolate.cpp b/src/Interpolation/Interpolate.cpp
index a14bc3b1bd..df8ba76738 100644
--- a/src/Interpolation/Interpolate.cpp
+++ b/src/Interpolation/Interpolate.cpp
@@ -339,6 +339,19 @@ void Interpolate_Iterate( real CData[], const int CSize[3], const int CStart[3],
# endif
+// use dual-energy fix before the general check
+# ifdef DUAL_ENERGY
+ const bool CheckMinPres_No = false;
+ const real UseDual2FixEngy = HUGE_NUMBER;
+ char dummy; // we do not record the dual-energy status here
+
+ if ( !FData_is_Prim )
+ Hydro_DualEnergyFix( Temp[DENS], Temp[MOMX], Temp[MOMY], Temp[MOMZ], Temp[ENGY], Temp[DUAL],
+ dummy, EoS_AuxArray_Flt[1], EoS_AuxArray_Flt[2],
+ CheckMinPres_No, NULL_REAL, UseDual2FixEngy, Emag );
+# endif
+
+
// 5-2. general check
bool Fail_ThisCell
= Hydro_IsUnphysical( (FData_is_Prim)?UNPHY_MODE_PRIM:UNPHY_MODE_CONS, Temp, NULL,
diff --git a/src/LoadBalance/LB_GatherTree.cpp b/src/LoadBalance/LB_GatherTree.cpp
index f74d8c4d46..e2579add1a 100644
--- a/src/LoadBalance/LB_GatherTree.cpp
+++ b/src/LoadBalance/LB_GatherTree.cpp
@@ -8,7 +8,7 @@ Instructions for adding new patch_t members to GatherTree:
- add lists to LB_LocalPatchExchangeList and LB_GlobalPatchExchangeList
-> modify "LB_GatherTree.cpp"
- read new member in LB_FillLocalExchangeList
-- transfer member in LB_FillGlobalExchangeList
+- transfer member in LB_FillGlobalPatchExchangeList
- write member to LB_GlobalPatch in LB_ConstructGlobalTree
*/
@@ -44,6 +44,8 @@ LB_LocalPatchExchangeList::LB_LocalPatchExchangeList() : isInitialised(false), L
{
LBIdxList_Local [lv] = new long [ amr->NPatchComma[lv][1] ];
CrList_Local [lv] = new int [ amr->NPatchComma[lv][1] ][3];
+ CrLList_Local [lv] = new int [ amr->NPatchComma[lv][1] ][3];
+ CrRList_Local [lv] = new int [ amr->NPatchComma[lv][1] ][3];
FaList_Local [lv] = new int [ amr->NPatchComma[lv][1] ];
SonList_Local [lv] = new int [ amr->NPatchComma[lv][1] ];
SibList_Local [lv] = new int [ amr->NPatchComma[lv][1] ][26];
@@ -69,6 +71,8 @@ LB_LocalPatchExchangeList::~LB_LocalPatchExchangeList() {
{
delete [] LBIdxList_Local[lv];
delete [] CrList_Local[lv];
+ delete [] CrLList_Local[lv];
+ delete [] CrRList_Local[lv];
delete [] FaList_Local[lv];
delete [] SonList_Local[lv];
delete [] SibList_Local[lv];
@@ -99,6 +103,8 @@ LB_GlobalPatchExchangeList::LB_GlobalPatchExchangeList( LB_PatchCount& pc, int r
if ( root < 0 || root == MPI_Rank) {
LBIdxList_AllLv = new long [ pc.NPatchAllLv ];
CrList_AllLv = new int [ pc.NPatchAllLv ][3];
+ CrLList_AllLv = new int [ pc.NPatchAllLv ][3];
+ CrRList_AllLv = new int [ pc.NPatchAllLv ][3];
FaList_AllLv = new int [ pc.NPatchAllLv ];
SonList_AllLv = new int [ pc.NPatchAllLv ];
SibList_AllLv = new int [ pc.NPatchAllLv ][26];
@@ -115,6 +121,8 @@ LB_GlobalPatchExchangeList::LB_GlobalPatchExchangeList( LB_PatchCount& pc, int r
} else {
LBIdxList_AllLv = NULL;
CrList_AllLv = NULL;
+ CrLList_AllLv = NULL;
+ CrRList_AllLv = NULL;
FaList_AllLv = NULL;
SonList_AllLv = NULL;
SibList_AllLv = NULL;
@@ -136,6 +144,8 @@ LB_GlobalPatchExchangeList::~LB_GlobalPatchExchangeList() {
if ( isAllocated ) {
delete [] LBIdxList_AllLv;
delete [] CrList_AllLv;
+ delete [] CrLList_AllLv;
+ delete [] CrRList_AllLv;
delete [] FaList_AllLv;
delete [] SonList_AllLv;
delete [] SibList_AllLv;
@@ -317,7 +327,11 @@ void LB_FillLocalPatchExchangeList( LB_PatchCount& pc, LB_LocalPatchExchangeList
// 2. corner
for (int d=0; d<3; d++)
- lel.CrList_Local[lv][PID][d] = amr->patch[0][lv][PID]->corner[d];
+ {
+ lel.CrList_Local[lv][PID][d] = amr->patch[0][lv][PID]->corner[d];
+ lel.CrLList_Local[lv][PID][d] = amr->patch[0][lv][PID]->cornerL[d];
+ lel.CrRList_Local[lv][PID][d] = amr->patch[0][lv][PID]->cornerR[d];
+ }
// 3. father GID
@@ -519,6 +533,8 @@ void LB_FillGlobalPatchExchangeList( LB_PatchCount& pc, LB_LocalPatchExchangeLis
// sending and receiving lists for MPI communication
int RecvCount_Cr [MPI_NRank], RecvDisp_Cr [MPI_NRank];
+ int RecvCount_CrL [MPI_NRank], RecvDisp_CrL [MPI_NRank];
+ int RecvCount_CrR [MPI_NRank], RecvDisp_CrR [MPI_NRank];
int RecvCount_Fa [MPI_NRank], RecvDisp_Fa [MPI_NRank];
int RecvCount_Son [MPI_NRank], RecvDisp_Son [MPI_NRank];
int RecvCount_Sib [MPI_NRank], RecvDisp_Sib [MPI_NRank];
@@ -540,6 +556,8 @@ void LB_FillGlobalPatchExchangeList( LB_PatchCount& pc, LB_LocalPatchExchangeLis
RecvCount_Son [r] = RecvCount_Fa[r];
RecvCount_Sib [r] = RecvCount_Fa[r]*26;
RecvCount_Cr [r] = RecvCount_Fa[r]*3;
+ RecvCount_CrL [r] = RecvCount_Fa[r]*3;
+ RecvCount_CrR [r] = RecvCount_Fa[r]*3;
RecvCount_EdgeL [r] = RecvCount_Fa[r]*3;
RecvCount_EdgeR [r] = RecvCount_Fa[r]*3;
RecvCount_PaddedCr1D[r] = RecvCount_Fa[r];
@@ -552,6 +570,8 @@ void LB_FillGlobalPatchExchangeList( LB_PatchCount& pc, LB_LocalPatchExchangeLis
RecvDisp_Son [r] = RecvDisp_Fa[r];
RecvDisp_Sib [r] = RecvDisp_Fa[r]*26;
RecvDisp_Cr [r] = RecvDisp_Fa[r]*3;
+ RecvDisp_CrL [r] = RecvDisp_Fa[r]*3;
+ RecvDisp_CrR [r] = RecvDisp_Fa[r]*3;
RecvDisp_EdgeL [r] = RecvDisp_Fa[r]*3;
RecvDisp_EdgeR [r] = RecvDisp_Fa[r]*3;
RecvDisp_PaddedCr1D [r] = RecvDisp_Fa[r];
@@ -575,6 +595,12 @@ void LB_FillGlobalPatchExchangeList( LB_PatchCount& pc, LB_LocalPatchExchangeLis
MPI_Allgatherv( lel.CrList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_INT,
(gel.CrList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_Cr, RecvDisp_Cr, MPI_INT, MPI_COMM_WORLD );
+ MPI_Allgatherv( lel.CrLList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_INT,
+ (gel.CrLList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_CrL, RecvDisp_CrL, MPI_INT, MPI_COMM_WORLD );
+
+ MPI_Allgatherv( lel.CrRList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_INT,
+ (gel.CrRList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_CrR, RecvDisp_CrR, MPI_INT, MPI_COMM_WORLD );
+
MPI_Allgatherv( lel.EdgeLList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_DOUBLE,
(gel.EdgeLList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_EdgeL, RecvDisp_EdgeL, MPI_DOUBLE, MPI_COMM_WORLD );
@@ -604,6 +630,12 @@ void LB_FillGlobalPatchExchangeList( LB_PatchCount& pc, LB_LocalPatchExchangeLis
MPI_Gatherv( lel.CrList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_INT,
(gel.CrList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_Cr, RecvDisp_Cr, MPI_INT, root, MPI_COMM_WORLD );
+ MPI_Gatherv( lel.CrLList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_INT,
+ (gel.CrLList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_CrL, RecvDisp_CrL, MPI_INT, root, MPI_COMM_WORLD );
+
+ MPI_Gatherv( lel.CrLList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_INT,
+ (gel.CrRList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_CrR, RecvDisp_CrR, MPI_INT, root, MPI_COMM_WORLD );
+
MPI_Gatherv( lel.EdgeLList_Local[lv][0], amr->NPatchComma[lv][1]*3, MPI_DOUBLE,
(gel.EdgeLList_AllLv+pc.GID_LvStart[lv])[0], RecvCount_EdgeL, RecvDisp_EdgeL, MPI_DOUBLE, root, MPI_COMM_WORLD );
@@ -675,6 +707,10 @@ LB_GlobalPatch* LB_ConstructGlobalTree( LB_PatchCount& pc, LB_GlobalPatchExchang
for (int c=0; c<3 ; c++)
global_tree[MyGID].corner[c] = gel.CrList_AllLv [MyGID][c];
for (int c=0; c<3 ; c++)
+ global_tree[MyGID].cornerL[c] = gel.CrList_AllLv [MyGID][c];
+ for (int c=0; c<3 ; c++)
+ global_tree[MyGID].cornerR[c] = gel.CrList_AllLv [MyGID][c];
+ for (int c=0; c<3 ; c++)
global_tree[MyGID].EdgeL[c] = gel.EdgeLList_AllLv [MyGID][c];
for (int c=0; c<3 ; c++)
global_tree[MyGID].EdgeR[c] = gel.EdgeRList_AllLv [MyGID][c];
diff --git a/src/Main/EvolveLevel.cpp b/src/Main/EvolveLevel.cpp
index 3305379464..bb3202529e 100644
--- a/src/Main/EvolveLevel.cpp
+++ b/src/Main/EvolveLevel.cpp
@@ -19,6 +19,7 @@ extern Timer_t *Timer_Par_2Son [NLEVEL];
#endif
bool AutoReduceDt_Continue;
+int Evolve_stage;
@@ -140,6 +141,19 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
# endif
// ===============================================================================================
+ Evolve_stage = 1;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+ if ( OPT__VERBOSE && MPI_Rank == 0 )
+ Aux_Message( stdout, " Lv %2d: Mis_UserWorkBeforeNextSubstep %4s... ", lv, "" );
+
+// use the same timer as the fluid solver for now
+ TIMING_FUNC( Mis_UserWorkBeforeNextSubstep_Ptr( lv, TimeNew, TimeOld, dt_SubStep ),
+ Timer_Flu_Advance[lv], TIMER_ON );
+
+ if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
+ }
+
// 2. fluid solver
// ===============================================================================================
@@ -312,6 +326,18 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
// ===============================================================================================
+ Evolve_stage = 2;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+ if ( OPT__VERBOSE && MPI_Rank == 0 )
+ Aux_Message( stdout, " Lv %2d: Mis_UserWorkBeforeNextSubstep %4s... ", lv, "" );
+
+// use the same timer as the fluid solver for now
+ TIMING_FUNC( Mis_UserWorkBeforeNextSubstep_Ptr( lv, TimeNew, TimeOld, dt_SubStep ),
+ Timer_Flu_Advance[lv], TIMER_ON );
+
+ if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
+ }
// 3. update particles (prediction for KDK) and exchange particles
@@ -352,7 +378,11 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
if ( OPT__VERBOSE && MPI_Rank == 0 )
Aux_Message( stdout, " Lv %2d: Gra_AdvanceDt, counter = %8ld ... ", lv, AdvanceCounter[lv] );
- if ( lv == 0 )
+ bool FullRefinedLv = false;
+ if ( (long)NPatchTotal[lv] == NX0_TOT[0]*NX0_TOT[1]*NX0_TOT[2]/512*(1L< 0
@@ -418,23 +448,45 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
// --> we will do this after all other operations (e.g., star formation) if OPT__MINIMIZE_MPI_BARRIER is adopted
// --> assuming that all remaining operations do not need to access the potential in the buffer patches
// --> one must enable both STORE_POT_GHOST and PAR_IMPROVE_ACC for this purpose
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
if ( UsePot && !OPT__MINIMIZE_MPI_BARRIER )
TIMING_FUNC( Buf_GetBufferData( lv, NULL_INT, NULL_INT, SaveSg_Pot, POT_FOR_POISSON,
_POTE, _NONE, Pot_ParaBuf, USELB_YES ),
Timer_GetBuf[lv][1], TIMER_ON );
+// Aux_Message( stdout, "DEBUG: %s %d\n", __FILE__, __LINE__ );
} // if ( OPT__OVERLAP_MPI ) ... else ...
+ // Aux_Message( stdout, "\n" );
+ // Aux_Message( stdout, "Track flu prepare %d %24.16e %d %24.16e %s %d\n", amr->FluSg[lv], amr->FluSgTime[lv][amr->FluSg[lv]], 1-amr->FluSg[lv], amr->FluSgTime[lv][1-amr->FluSg[lv]], __FUNCTION__, __LINE__ );
+ // Aux_Message( stdout, "Track pot prepare %d %24.16e %d %24.16e %s %d\n", amr->PotSg[lv], amr->PotSgTime[lv][amr->PotSg[lv]], 1-amr->PotSg[lv], amr->PotSgTime[lv][1-amr->PotSg[lv]], __FUNCTION__, __LINE__ );
if ( UsePot )
{
amr->PotSg [lv] = SaveSg_Pot;
amr->PotSgTime[lv][SaveSg_Pot] = TimeNew;
}
+ // Aux_Message( stdout, "\n" );
+ // Aux_Message( stdout, "Track flu prepare %d %24.16e %d %24.16e %s %d\n", amr->FluSg[lv], amr->FluSgTime[lv][amr->FluSg[lv]], 1-amr->FluSg[lv], amr->FluSgTime[lv][1-amr->FluSg[lv]], __FUNCTION__, __LINE__ );
+ // Aux_Message( stdout, "Track pot prepare %d %24.16e %d %24.16e %s %d\n", amr->PotSg[lv], amr->PotSgTime[lv][amr->PotSg[lv]], 1-amr->PotSg[lv], amr->PotSgTime[lv][1-amr->PotSg[lv]], __FUNCTION__, __LINE__ );
} // if ( lv == 0 ) ... else ...
if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
# endif // #ifdef GRAVITY
// ===============================================================================================
+ // Aux_Message( stdout, "Evolve after gra flu prepare %d %24.16e %d %24.16e\n", amr->FluSg[lv], amr->FluSgTime[lv][amr->FluSg[lv]], 1-amr->FluSg[lv], amr->FluSgTime[lv][1-amr->FluSg[lv]] );
+ // Aux_Message( stdout, "Evolve after gra pot prepare %d %24.16e %d %24.16e\n", amr->PotSg[lv], amr->PotSgTime[lv][amr->PotSg[lv]], 1-amr->PotSg[lv], amr->PotSgTime[lv][1-amr->PotSg[lv]] );
+ Evolve_stage = 4;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+ if ( OPT__VERBOSE && MPI_Rank == 0 )
+ Aux_Message( stdout, " Lv %2d: Mis_UserWorkBeforeNextSubstep %4s... ", lv, "" );
+
+// use the same timer as the fluid solver for now
+ TIMING_FUNC( Mis_UserWorkBeforeNextSubstep_Ptr( lv, TimeNew, TimeOld, dt_SubStep ),
+ Timer_Flu_Advance[lv], TIMER_ON );
+
+ if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
+ }
// 5. correct particles velocity and send particles to lv+1
@@ -697,6 +749,20 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
amr->NUpdateLv[lv] ++;
if ( AdvanceCounter[lv] >= __LONG_MAX__ ) Aux_Message( stderr, "WARNING : AdvanceCounter overflow !!\n" );
+// ===============================================================================================
+ Evolve_stage = 10;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+ if ( OPT__VERBOSE && MPI_Rank == 0 )
+ Aux_Message( stdout, " Lv %2d: Mis_UserWorkBeforeNextSubstep %4s... ", lv, "" );
+
+// use the same timer as the fluid solver for now
+ TIMING_FUNC( Mis_UserWorkBeforeNextSubstep_Ptr( lv, TimeNew, TimeOld, dt_SubStep ),
+ Timer_Flu_Advance[lv], TIMER_ON );
+
+ if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
+ }
+// ===============================================================================================
if ( lv != TOP_LEVEL && NPatchTotal[lv+1] != 0 )
@@ -716,6 +782,20 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
Timer_Lv[lv]->Start();
# endif
// ===============================================================================================
+// ===============================================================================================
+ Evolve_stage = 11;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+ if ( OPT__VERBOSE && MPI_Rank == 0 )
+ Aux_Message( stdout, " Lv %2d: Mis_UserWorkBeforeNextSubstep %4s... ", lv, "" );
+
+// use the same timer as the fluid solver for now
+ TIMING_FUNC( Mis_UserWorkBeforeNextSubstep_Ptr( lv, TimeNew, TimeOld, dt_SubStep ),
+ Timer_Flu_Advance[lv], TIMER_ON );
+
+ if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
+ }
+// ===============================================================================================
// 12. correct the data at the current level with the data at the next finer level
@@ -818,8 +898,24 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
// ===============================================================================================
-
} // if ( lv != TOP_LEVEL && NPatchTotal[lv+1] != 0 )
+// ===============================================================================================
+ Evolve_stage = 12;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+ if ( OPT__VERBOSE && MPI_Rank == 0 )
+ Aux_Message( stdout, " Lv %2d: Mis_UserWorkBeforeNextSubstep %4s... ", lv, "" );
+
+// use the same timer as the fluid solver for now
+ TIMING_FUNC( Mis_UserWorkBeforeNextSubstep_Ptr( lv, TimeNew, TimeOld, dt_SubStep ),
+ Timer_Flu_Advance[lv], TIMER_ON );
+
+ if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, "done\n" );
+ }
+// ===============================================================================================
+
+
+
// 13. refine to higher level(s)
@@ -937,9 +1033,9 @@ void EvolveLevel( const int lv, const double dTime_FaLv )
} // if ( lv != TOP_LEVEL && AdvanceCounter[lv] % REGRID_COUNT == 0 )
// ===============================================================================================
-
// 14. user-specified operations before proceeding to the next sub-step
// ===============================================================================================
+ Evolve_stage = 13;
if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
{
if ( OPT__VERBOSE && MPI_Rank == 0 )
diff --git a/src/Main/InterpolateGhostZone.cpp b/src/Main/InterpolateGhostZone.cpp
index 7506672dac..afe83bd48b 100644
--- a/src/Main/InterpolateGhostZone.cpp
+++ b/src/Main/InterpolateGhostZone.cpp
@@ -1407,7 +1407,10 @@ void InterpolateGhostZone( const int lv, const int PID, real IntData_CC[], real
// determine which variables require **monotonic** interpolation
- bool Monotonicity_CC[NVarCC_Flu];
+ // bool Monotonicity_CC[NVarCC_Flu];
+ bool *Monotonicity_CC;
+ if ( NVarCC_Flu == 0 ) Monotonicity_CC = new bool [ 1 ];
+ else Monotonicity_CC = new bool [ NVarCC_Flu ];
for (int v=0; vdh[JEANS_MIN_PRES_LEVEL])/(GAMMA*M_PI) : NULL_REAL;
+# else
const real JeansMinPres_Coeff = ( JEANS_MIN_PRES ) ?
- NEWTON_G*SQR(JEANS_MIN_PRES_NCELL*amr->dh[JEANS_MIN_PRES_LEVEL])/(GAMMA*M_PI) : NULL_REAL;
+ NEWTON_G*SQR(JEANS_MIN_PRES_NCELL*amr->dh[JEANS_MIN_PRES_LEVEL])/(GAMMA*M_PI) : NULL_REAL;
+# endif
# else
const real JEANS_MIN_PRES = false;
const real JeansMinPres_Coeff = NULL_REAL;
diff --git a/src/Main/Main.cpp b/src/Main/Main.cpp
index 703305e807..5896201462 100644
--- a/src/Main/Main.cpp
+++ b/src/Main/Main.cpp
@@ -170,7 +170,9 @@ bool ELBDM_BASE_SPECTRAL;
double AveDensity_Init = -1.0; // initialize it as <= 0 to check if it is properly set later
int Pot_ParaBuf, Rho_ParaBuf;
-real *GreenFuncK = NULL;
+bool FFTW_Inited[NLEVEL] = { false };
+bool GreenFuncK_Inited[NLEVEL] = { false };
+real *GreenFuncK[NLEVEL] = { NULL };
double GFUNC_COEFF0;
double DT__GRAVITY;
double NEWTON_G;
@@ -368,6 +370,23 @@ double DT__CR_DIFFUSION;
double CR_DIFF_MIN_B;
#endif
+// (2-16) Hypre
+#ifdef SUPPORT_HYPRE
+Hypre_Solver_t HYPRE_SOLVER;
+bool HYPRE_INIT_GUESS;
+int HYPRE_PRINT_LEVEL, HYPRE_ENABLE_LOGGING;
+int HYPRE_MAX_ITER, HYPRE_NPRE_RELAX, HYPRE_NPOST_RELAX;
+double HYPRE_REL_TOL, HYPRE_ABS_TOL;
+// TODO : each level has one set of these?
+HYPRE_SStructGrid Hypre_grid;
+HYPRE_SStructGraph Hypre_graph;
+HYPRE_SStructStencil Hypre_stencil;
+HYPRE_SStructMatrix Hypre_A;
+HYPRE_SStructVector Hypre_b;
+HYPRE_SStructVector Hypre_x;
+HYPRE_SStructSolver Hypre_solver;
+#endif
+
// 3. CPU (host) arrays for transferring data between CPU and GPU
// =======================================================================================================
@@ -594,6 +613,7 @@ Timer_t Timer_OutputWalltime;
+extern int Evolve_stage;
//-------------------------------------------------------------------------------------------------------
// Function : main
@@ -697,6 +717,14 @@ int main( int argc, char *argv[] )
if ( OPT__CORR_AFTER_ALL_SYNC == CORR_AFTER_SYNC_EVERY_STEP )
TIMING_FUNC( Flu_CorrAfterAllSync(), Timer_Main[6], TIMER_ON );
+ Evolve_stage = 14;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+// use the same timer as the fluid solver for now
+ Mis_UserWorkBeforeNextSubstep_Ptr( 0, NULL_REAL, NULL_REAL, NULL_REAL );
+ Mis_UserWorkBeforeNextSubstep_Ptr( 1, NULL_REAL, NULL_REAL, NULL_REAL );
+ }
+
# if ( MODEL == HYDRO && defined MHD )
if ( OPT__SAME_INTERFACE_B )
{
@@ -716,6 +744,13 @@ int main( int argc, char *argv[] )
// 3. output data and execute auxiliary functions
// ---------------------------------------------------------------------------------------------------
TIMING_FUNC( Output_DumpData( 1 ), Timer_Main[3], TIMER_ON );
+ Evolve_stage = 15;
+ if ( Mis_UserWorkBeforeNextSubstep_Ptr != NULL )
+ {
+// use the same timer as the fluid solver for now
+ Mis_UserWorkBeforeNextSubstep_Ptr( 0, NULL_REAL, NULL_REAL, NULL_REAL );
+ Mis_UserWorkBeforeNextSubstep_Ptr( 1, NULL_REAL, NULL_REAL, NULL_REAL );
+ }
if ( OPT__PATCH_COUNT == 1 )
TIMING_FUNC( Aux_Record_PatchCount(), Timer_Main[4], TIMER_ON );
@@ -750,6 +785,7 @@ int main( int argc, char *argv[] )
// ---------------------------------------------------------------------------------------------------
+
// 4. perform yt inline analysis
// ---------------------------------------------------------------------------------------------------
# ifdef SUPPORT_LIBYT
diff --git a/src/Makefile_base b/src/Makefile_base
index 05a5ff17b0..73ff15b9d3 100644
--- a/src/Makefile_base
+++ b/src/Makefile_base
@@ -53,6 +53,7 @@ GRACKLE_PATH := @@@GRACKLE_PATH@@@
GSL_PATH := @@@GSL_PATH@@@
LIBYT_PATH := @@@LIBYT_PATH@@@
CUFFTDX_PATH := @@@CUFFTDX_PATH@@@
+HYPRE_PATH := @@@HYPRE_PATH@@@
@@ -203,7 +204,7 @@ GPU_FILE += CUPOT_PoissonSolver_SOR.cu \
CUPOT_ExtPotSolver.cu CUPOT_ExtPot_Tabular.cu
CPU_FILE += CPU_PoissonGravitySolver.cpp CPU_PoissonSolver_SOR.cpp CPU_PoissonSolver_FFT.cpp \
- CPU_PoissonSolver_MG.cpp CPU_ExtPotSolver.cpp CPU_ExtPotSolver_BaseLevel.cpp
+ CPU_PoissonSolver_MG.cpp CPU_ExtPotSolver.cpp CPU_ExtPotSolver_FullyRefinedLevel.cpp
CPU_FILE += Gra_Close.cpp Gra_Prepare_Flu.cpp Gra_Prepare_Pot.cpp Gra_Prepare_Corner.cpp \
Gra_AdvanceDt.cpp Poi_Close.cpp Poi_Prepare_Pot.cpp Poi_Prepare_Rho.cpp \
@@ -303,6 +304,16 @@ vpath %.cpp YT
endif # SUPPORT_LIBYT
+# HYPRE source files (included only if "SUPPORT_HYPRE" is turned on)
+# ------------------------------------------------------------------------------------
+ifeq "$(filter -DSUPPORT_HYPRE, $(SIMU_OPTION))" "-DSUPPORT_HYPRE"
+CPU_FILE += Hypre_Init.cpp Hypre_End.cpp Hypre_Main.cpp Hypre_Solve.cpp Hypre_Aux_Record.cpp \
+ Hypre_PrepareSingleLevel.cpp Hypre_FillArrays.cpp Hypre_UpdateArrays.cpp Hypre_Free.cpp
+
+vpath %.cpp Hypre
+endif # SUPPORT_HYPRE
+
+
# local source terms source files
# ------------------------------------------------------------------------------------
GPU_FILE += CUAPI_Asyn_SrcSolver.cu CUSRC_SrcSolver_IterateAllCells.cu CUSRC_Src_Deleptonization.cu \
@@ -467,6 +478,14 @@ LIB += -L$(LIBYT_PATH)/lib -lyt
LIB += -Wl,-rpath,$(LIBYT_PATH)/lib
endif
+ifeq "$(filter -DSUPPORT_HYPRE, $(SIMU_OPTION))" "-DSUPPORT_HYPRE"
+LIB += -L$(HYPRE_PATH)/lib -lHYPRE
+LIB += -Wl,-rpath=$(HYPRE_PATH)/lib
+ifeq "$(filter -DGPU, $(SIMU_OPTION))" "-DGPU"
+LIB += -lcusparse -lcublas -lcurand
+endif
+endif
+
# headers
# -------------------------------------------------------------------------------
@@ -504,6 +523,13 @@ ifeq "$(filter -DSUPPORT_LIBYT, $(SIMU_OPTION))" "-DSUPPORT_LIBYT"
INCLUDE += -I$(LIBYT_PATH)/include
endif
+ifeq "$(filter -DSUPPORT_HYPRE, $(SIMU_OPTION))" "-DSUPPORT_HYPRE"
+INCLUDE += -I$(HYPRE_PATH)/include
+ifeq "$(filter -DGPU, $(SIMU_OPTION))" "-DGPU"
+INCLUDE += -I$(CUDA_PATH)/include
+endif
+endif
+
ifeq "$(filter -DGPU, $(SIMU_OPTION))" "-DGPU"
ifeq "$(filter -DWAVE_SCHEME=WAVE_GRAMFE, $(SIMU_OPTION))" "-DWAVE_SCHEME=WAVE_GRAMFE"
ifeq "$(filter -DGRAMFE_SCHEME=GRAMFE_FFT, $(SIMU_OPTION))" "-DGRAMFE_SCHEME=GRAMFE_FFT"
@@ -575,6 +601,7 @@ GRACKLE_PATH := $(strip $(GRACKLE_PATH))
GSL_PATH := $(strip $(GSL_PATH))
LIBYT_PATH := $(strip $(LIBYT_PATH))
CUFFTDX_PATH := $(strip $(CUFFTDX_PATH))
+HYPRE_PATH := $(strip $(HYPRE_PATH))
# implicit rules (do NOT modify the order of the following rules)
diff --git a/src/Model_ELBDM/CPU_ELBDM/CPU_ELBDMSolver_FFT.cpp b/src/Model_ELBDM/CPU_ELBDM/CPU_ELBDMSolver_FFT.cpp
index 63c014bce6..80319a1a29 100644
--- a/src/Model_ELBDM/CPU_ELBDM/CPU_ELBDMSolver_FFT.cpp
+++ b/src/Model_ELBDM/CPU_ELBDM/CPU_ELBDMSolver_FFT.cpp
@@ -136,6 +136,8 @@ void Psi_Advance_FFT( real *PsiR, real *PsiI, const int j_start, const int dj, c
void CPU_ELBDMSolver_FFT( const real dt, const double PrepTime, const int SaveSg )
{
+ const int lv = 0;
+
// determine the FFT size
const int FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
@@ -199,9 +201,9 @@ void CPU_ELBDMSolver_FFT( const real dt, const double PrepTime, const int SaveSg
real *PsiR = (real*)root_fftw::fft_malloc( sizeof(real)*total_local_size ); // array storing real and imaginary parts of wave function
real *PsiI = (real*)root_fftw::fft_malloc( sizeof(real)*total_local_size );
- real *SendBuf = new real [ (long)amr->NPatchComma[0][1]*CUBE(PS1) ]; // MPI send buffer
+ real *SendBuf = new real [ (long)amr->NPatchComma[lv][1]*CUBE(PS1) ]; // MPI send buffer
real *RecvBuf = new real [ (long)NX0_TOT[0]*NX0_TOT[1]*NRecvSlice ]; // MPI recv buffer
- long *SendBuf_SIdx = new long [ (long)amr->NPatchComma[0][1]*PS1 ]; // MPI send buffer for 1D coordinate in slab
+ long *SendBuf_SIdx = new long [ (long)amr->NPatchComma[lv][1]*PS1 ]; // MPI send buffer for 1D coordinate in slab
long *RecvBuf_SIdx = new long [ (long)NX0_TOT[0]*NX0_TOT[1]*NRecvSlice/SQR(PS1) ]; // MPI recv buffer for 1D coordinate in slab
int *List_PID_R [MPI_NRank]; // PID of each patch slice sent to each rank for the real part
@@ -214,9 +216,9 @@ void CPU_ELBDMSolver_FFT( const real dt, const double PrepTime, const int SaveSg
// rearrange data from patch to slab
Patch2Slab( PsiR, SendBuf, RecvBuf, SendBuf_SIdx, RecvBuf_SIdx, List_PID_R, List_k_R, List_NSend, List_NRecv, List_z_start,
- local_nz, FFT_Size, NRecvSlice, PrepTime, _REAL, InPlacePad_No, ForPoisson_No, false );
+ local_nz, FFT_Size, NRecvSlice, PrepTime, _REAL, InPlacePad_No, ForPoisson_No, false, lv );
Patch2Slab( PsiI, SendBuf, RecvBuf, SendBuf_SIdx, RecvBuf_SIdx, List_PID_I, List_k_I, List_NSend, List_NRecv, List_z_start,
- local_nz, FFT_Size, NRecvSlice, PrepTime, _IMAG, InPlacePad_No, ForPoisson_No, false );
+ local_nz, FFT_Size, NRecvSlice, PrepTime, _IMAG, InPlacePad_No, ForPoisson_No, false, lv );
// advance wave function by exp( -i*dt*k^2/(2*ELBDM_ETA) ) in the k-space using FFT
@@ -225,23 +227,23 @@ void CPU_ELBDMSolver_FFT( const real dt, const double PrepTime, const int SaveSg
// rearrange data from slab back to patch
Slab2Patch( PsiR, RecvBuf, SendBuf, SaveSg, RecvBuf_SIdx, List_PID_R, List_k_R, List_NRecv, List_NSend,
- local_nz, FFT_Size, NRecvSlice, _REAL, InPlacePad_No );
+ local_nz, FFT_Size, NRecvSlice, _REAL, InPlacePad_No, lv );
Slab2Patch( PsiI, RecvBuf, SendBuf, SaveSg, RecvBuf_SIdx, List_PID_I, List_k_I, List_NRecv, List_NSend,
- local_nz, FFT_Size, NRecvSlice, _IMAG, InPlacePad_No );
+ local_nz, FFT_Size, NRecvSlice, _IMAG, InPlacePad_No, lv );
// update density according to the updated wave function
# pragma omp parallel for schedule( runtime )
- for (int PID=0; PIDNPatchComma[0][1]; PID++)
+ for (int PID=0; PIDNPatchComma[lv][1]; PID++)
for (int k=0; kpatch[SaveSg][0][PID]->fluid[REAL][k][j][i];
- const real NewImag = amr->patch[SaveSg][0][PID]->fluid[IMAG][k][j][i];
+ const real NewReal = amr->patch[SaveSg][lv][PID]->fluid[REAL][k][j][i];
+ const real NewImag = amr->patch[SaveSg][lv][PID]->fluid[IMAG][k][j][i];
const real NewDens = SQR( NewReal ) + SQR( NewImag );
- amr->patch[SaveSg][0][PID]->fluid[DENS][k][j][i] = NewDens;
+ amr->patch[SaveSg][lv][PID]->fluid[DENS][k][j][i] = NewDens;
} // PID,i,j,k
diff --git a/src/Output/Output_BasePowerSpectrum.cpp b/src/Output/Output_BasePowerSpectrum.cpp
index e9565aa517..d94f6e4a3c 100644
--- a/src/Output/Output_BasePowerSpectrum.cpp
+++ b/src/Output/Output_BasePowerSpectrum.cpp
@@ -32,6 +32,7 @@ void Output_BasePowerSpectrum( const char *FileName, const long TVar )
if ( NX0_TOT[0] != NX0_TOT[1] || NX0_TOT[0] != NX0_TOT[2] )
Aux_Error( ERROR_INFO, "%s only works with CUBIC domain !!\n", __FUNCTION__ );
+ const int lv = 0;
// 1. determine the FFT size
const int Nx_Padded = NX0_TOT[0]/2+1;
@@ -93,9 +94,9 @@ void Output_BasePowerSpectrum( const char *FileName, const long TVar )
double *PS_total = NULL;
real *VarK = (real*)root_fftw::fft_malloc( sizeof(real)*total_local_size ); // array storing data
- real *SendBuf = new real [ (long)amr->NPatchComma[0][1]*CUBE(PS1) ]; // MPI send buffer for data
+ real *SendBuf = new real [ (long)amr->NPatchComma[lv][1]*CUBE(PS1) ]; // MPI send buffer for data
real *RecvBuf = new real [ (long)NX0_TOT[0]*NX0_TOT[1]*NRecvSlice ]; // MPI recv buffer for data
- long *SendBuf_SIdx = new long [ (long)amr->NPatchComma[0][1]*PS1 ]; // MPI send buffer for 1D coordinate in slab
+ long *SendBuf_SIdx = new long [ (long)amr->NPatchComma[lv][1]*PS1 ]; // MPI send buffer for 1D coordinate in slab
long *RecvBuf_SIdx = new long [ (long)NX0_TOT[0]*NX0_TOT[1]*NRecvSlice/SQR(PS1) ]; // MPI recv buffer for 1D coordinate in slab
int *List_PID [MPI_NRank]; // PID of each patch slice sent to each rank
@@ -123,17 +124,17 @@ void Output_BasePowerSpectrum( const char *FileName, const long TVar )
# endif
if ( TVar == _TOTAL_DENS ) {
- Par_CollectParticle2OneLevel( 0, _PAR_MASS|_PAR_POSX|_PAR_POSY|_PAR_POSZ, _PAR_TYPE, PredictPos, Time[0],
+ Par_CollectParticle2OneLevel( lv, _PAR_MASS|_PAR_POSX|_PAR_POSY|_PAR_POSZ, _PAR_TYPE, PredictPos, Time[lv],
SibBufPatch, FaSibBufPatch, JustCountNPar_No, TimingSendPar_No );
- Prepare_PatchData_InitParticleDensityArray( 0, Time[0] );
+ Prepare_PatchData_InitParticleDensityArray( lv, Time[lv] );
} // if ( TVar == _TOTAL_DENS )
# endif // #ifdef MASSIVE_PARTICLES
// 4. rearrange data from patch to slab
Patch2Slab( VarK, SendBuf, RecvBuf, SendBuf_SIdx, RecvBuf_SIdx, List_PID, List_k, List_NSend, List_NRecv, List_z_start,
- local_nz, FFT_Size, NRecvSlice, Time[0], TVar, InPlacePad, ForPoisson, false );
+ local_nz, FFT_Size, NRecvSlice, Time[lv], TVar, InPlacePad, ForPoisson, false, lv );
// 5. evaluate the base-level power spectrum by FFT
@@ -185,9 +186,9 @@ void Output_BasePowerSpectrum( const char *FileName, const long TVar )
// free memory for collecting particles from other ranks and levels, and free density arrays with ghost zones (rho_ext)
# ifdef MASSIVE_PARTICLES
if ( TVar == _TOTAL_DENS ) {
- Par_CollectParticle2OneLevel_FreeMemory( 0, SibBufPatch, FaSibBufPatch );
+ Par_CollectParticle2OneLevel_FreeMemory( lv, SibBufPatch, FaSibBufPatch );
- Prepare_PatchData_FreeParticleDensityArray( 0 );
+ Prepare_PatchData_FreeParticleDensityArray( lv );
}
# endif
diff --git a/src/Output/Output_DumpData_Total_HDF5.cpp b/src/Output/Output_DumpData_Total_HDF5.cpp
index 046ebd810d..8b3bfc5908 100644
--- a/src/Output/Output_DumpData_Total_HDF5.cpp
+++ b/src/Output/Output_DumpData_Total_HDF5.cpp
@@ -79,7 +79,7 @@ Procedure for outputting new variables:
//-------------------------------------------------------------------------------------------------------
-// Function : Output_DumpData_Total_HDF5 (FormatVersion = 2503)
+// Function : Output_DumpData_Total_HDF5 (FormatVersion = 2510)
// Description : Output all simulation data in the HDF5 format, which can be used as a restart file
// or loaded by YT
//
@@ -277,6 +277,7 @@ Procedure for outputting new variables:
// 2502 : 2025/01/16 --> output ConRef[]
// 2503 : 2025/01/17 --> output user-defined parameters in "User/UserPara" and
// Input__TestProb parameters in "Info/InputTest"
+// 2510 : 2025/**/** --> output Hypre options
//-------------------------------------------------------------------------------------------------------
void Output_DumpData_Total_HDF5( const char *FileName )
{
@@ -1662,7 +1663,7 @@ void FillIn_KeyInfo( KeyInfo_t &KeyInfo, const int NFieldStored )
const time_t CalTime = time( NULL ); // calendar time
- KeyInfo.FormatVersion = 2503;
+ KeyInfo.FormatVersion = 2510;
KeyInfo.Model = MODEL;
KeyInfo.NLevel = NLEVEL;
KeyInfo.NCompFluid = NCOMP_FLUID;
@@ -1930,6 +1931,12 @@ void FillIn_Makefile( Makefile_t &Makefile )
Makefile.SupportGrackle = 0;
# endif
+# ifdef SUPPORT_HYPRE
+ Makefile.SupportHypre = 1;
+# else
+ Makefile.SupportHypre = 0;
+# endif
+
Makefile.RandomNumber = RANDOM_NUMBER;
Makefile.NLevel = NLEVEL;
@@ -2848,6 +2855,16 @@ void FillIn_InputPara( InputPara_t &InputPara, const int NFieldStored, char Fiel
InputPara.Yt_JupyterUseConnectionFile = YT_JUPYTER_USE_CONNECTION_FILE;
# endif
+# ifdef SUPPORT_HYPRE
+ InputPara.Hypre_Solver = HYPRE_SOLVER;
+ InputPara.Hypre_InitGuess = HYPRE_INIT_GUESS;
+ InputPara.Hypre_PrintLevel = HYPRE_PRINT_LEVEL;
+ InputPara.Hypre_EnableLogging = HYPRE_ENABLE_LOGGING;
+ InputPara.Hypre_MaxIter = HYPRE_MAX_ITER;
+ InputPara.Hypre_RelTol = HYPRE_REL_TOL;
+ InputPara.Hypre_AbsTol = HYPRE_ABS_TOL;
+# endif
+
// miscellaneous
InputPara.Opt__Verbose = OPT__VERBOSE;
InputPara.Opt__TimingBarrier = OPT__TIMING_BARRIER;
@@ -3109,6 +3126,7 @@ void GetCompound_Makefile( hid_t &H5_TypeID )
H5Tinsert( H5_TypeID, "LibYTJupyter", HOFFSET(Makefile_t,LibYTJupyter ), H5T_NATIVE_INT );
# endif
H5Tinsert( H5_TypeID, "SupportGrackle", HOFFSET(Makefile_t,SupportGrackle ), H5T_NATIVE_INT );
+ H5Tinsert( H5_TypeID, "SupportHypre", HOFFSET(Makefile_t,SupportHypre ), H5T_NATIVE_INT );
H5Tinsert( H5_TypeID, "RandomNumber", HOFFSET(Makefile_t,RandomNumber ), H5T_NATIVE_INT );
H5Tinsert( H5_TypeID, "NLevel", HOFFSET(Makefile_t,NLevel ), H5T_NATIVE_INT );
@@ -3896,6 +3914,17 @@ void GetCompound_InputPara( hid_t &H5_TypeID, const int NFieldStored )
H5Tinsert( H5_TypeID, "Yt_JupyterUseConnectionFile", HOFFSET(InputPara_t,Yt_JupyterUseConnectionFile), H5T_NATIVE_INT );
# endif
+// Hypre
+# ifdef SUPPORT_HYPRE
+ H5Tinsert( H5_TypeID, "Hypre_Solver", HOFFSET(InputPara_t,Hypre_Solver ), H5T_NATIVE_INT );
+ H5Tinsert( H5_TypeID, "Hypre_InitGuess", HOFFSET(InputPara_t,Hypre_InitGuess ), H5T_NATIVE_INT );
+ H5Tinsert( H5_TypeID, "Hypre_PrintLevel", HOFFSET(InputPara_t,Hypre_PrintLevel ), H5T_NATIVE_INT );
+ H5Tinsert( H5_TypeID, "Hypre_EnableLogging", HOFFSET(InputPara_t,Hypre_EnableLogging ), H5T_NATIVE_INT );
+ H5Tinsert( H5_TypeID, "Hypre_MaxIter", HOFFSET(InputPara_t,Hypre_MaxIter ), H5T_NATIVE_INT );
+ H5Tinsert( H5_TypeID, "Hypre_RelTol", HOFFSET(InputPara_t,Hypre_RelTol ), H5T_NATIVE_DOUBLE );
+ H5Tinsert( H5_TypeID, "Hypre_AbsTol", HOFFSET(InputPara_t,Hypre_AbsTol ), H5T_NATIVE_DOUBLE );
+# endif
+
// miscellaneous
H5Tinsert( H5_TypeID, "Opt__Verbose", HOFFSET(InputPara_t,Opt__Verbose ), H5T_NATIVE_INT );
H5Tinsert( H5_TypeID, "Opt__TimingBarrier", HOFFSET(InputPara_t,Opt__TimingBarrier ), H5T_NATIVE_INT );
diff --git a/src/Refine/Flag_Check.cpp b/src/Refine/Flag_Check.cpp
index bcd9de5419..ebd6801851 100644
--- a/src/Refine/Flag_Check.cpp
+++ b/src/Refine/Flag_Check.cpp
@@ -42,6 +42,7 @@ static bool Check_Radial( const int i, const int j, const int k, const int lv, c
// ParDens : Input array storing the particle mass density on each cell
// JeansCoeff : Pi*GAMMA/(SafetyFactor^2*G), where SafetyFactor = FlagTable_Jeans[lv]
// --> Flag if dh^2 > JeansCoeff*Pres/Dens^2
+// --> When COMOVING is on, G has been replaced by a*G, where a is the scale factor
// Interf_Var : Input array storing the density and phase for the interference condition
// Spectral_Cond : Input variable storing the spectral refinement condition
//
diff --git a/src/Refine/Flag_Real.cpp b/src/Refine/Flag_Real.cpp
index ec721d8492..5dd1dabef0 100644
--- a/src/Refine/Flag_Real.cpp
+++ b/src/Refine/Flag_Real.cpp
@@ -76,7 +76,11 @@ void Flag_Real( const int lv, const UseLBFunc_t UseLBFunc )
# endif // # if ( MODEL == ELBDM )
# if ( MODEL == HYDRO && defined GRAVITY )
- const real JeansCoeff_Factor = M_PI/( SQR(FlagTable_Jeans[lv])*NEWTON_G ); // flag if dh^2 > JeansCoeff_Factor*Gamma*Pres/Dens^2
+# ifdef COMOVING
+ const real JeansCoeff_Factor = M_PI/( SQR(FlagTable_Jeans[lv])*NEWTON_G*Time[lv] ); // flag if dh^2 > JeansCoeff_Factor*Gamma*Pres/Dens^2
+# else
+ const real JeansCoeff_Factor = M_PI/( SQR(FlagTable_Jeans[lv])*NEWTON_G );
+# endif
# endif
# ifndef GRAVITY
const OptPotBC_t OPT__BC_POT = BC_POT_NONE;
diff --git a/src/SelfGravity/CPU_Poisson/CPU_ExtPotSolver_BaseLevel.cpp b/src/SelfGravity/CPU_Poisson/CPU_ExtPotSolver_FullyRefinedLevel.cpp
similarity index 84%
rename from src/SelfGravity/CPU_Poisson/CPU_ExtPotSolver_BaseLevel.cpp
rename to src/SelfGravity/CPU_Poisson/CPU_ExtPotSolver_FullyRefinedLevel.cpp
index 6272823cbe..9990cbe267 100644
--- a/src/SelfGravity/CPU_Poisson/CPU_ExtPotSolver_BaseLevel.cpp
+++ b/src/SelfGravity/CPU_Poisson/CPU_ExtPotSolver_FullyRefinedLevel.cpp
@@ -6,8 +6,8 @@
//-----------------------------------------------------------------------------------------
-// Function : CPU_ExtPotSolver_BaseLevel
-// Description : Add external potential on the base level
+// Function : CPU_ExtPotSolver_FullyRefinedLevel
+// Description : Add external potential on the fully refined level
//
// Note : 1. External potential is specified by the input function Func()
// 2. Set PotIsInit to false if the base-level potential has not been initialized
@@ -23,12 +23,13 @@
// --> true : **add** to the original data
// false: **overwrite** the original data
// SaveSg : Sandglass to store the updated potential
+// lv : Target level
//
// Return : amr->patch->pot[]
//-----------------------------------------------------------------------------------------
-void CPU_ExtPotSolver_BaseLevel( const ExtPot_t Func, const double AuxArray_Flt[], const int AuxArray_Int[],
- const real Table[], void **GenePtr,
- const double Time, const bool PotIsInit, const int SaveSg )
+void CPU_ExtPotSolver_FullyRefinedLevel( const ExtPot_t Func, const double AuxArray_Flt[], const int AuxArray_Int[],
+ const real Table[], void **GenePtr,
+ const double Time, const bool PotIsInit, const int SaveSg, const int lv )
{
// check
@@ -50,7 +51,6 @@ void CPU_ExtPotSolver_BaseLevel( const ExtPot_t Func, const double AuxArray_Flt[
# endif
- const int lv = 0;
const double dh = amr->dh[lv];
const double dh_2 = 0.5*dh;
@@ -79,7 +79,7 @@ void CPU_ExtPotSolver_BaseLevel( const ExtPot_t Func, const double AuxArray_Flt[
} // for (int PID=0; PIDNPatchComma[lv][1]; PID++)
} // end of OpenMP parallel region
-} // FUNCTION : CPU_ExtPotSolver_BaseLevel
+} // FUNCTION : CPU_ExtPotSolver_FullyRefinedLevel
diff --git a/src/SelfGravity/CPU_Poisson/CPU_PoissonGravitySolver.cpp b/src/SelfGravity/CPU_Poisson/CPU_PoissonGravitySolver.cpp
index 1a6500929a..eb1231a2d2 100644
--- a/src/SelfGravity/CPU_Poisson/CPU_PoissonGravitySolver.cpp
+++ b/src/SelfGravity/CPU_Poisson/CPU_PoissonGravitySolver.cpp
@@ -189,6 +189,10 @@ void CPU_PoissonGravitySolver( const real h_Rho_Array [][RHO_NXT][RHO_NXT][RH
MG_Max_Iter, MG_NPre_Smooth, MG_NPost_Smooth, MG_Tolerated_Error,
Poi_Coeff, IntScheme );
+# elif ( POT_SCHEME == HYPRE_POI )
+
+ Aux_Error( ERROR_INFO, "POT_SCHEME == HYPRE_POISSON is not supported yet !!\n" );
+
# else
# error : ERROR : unsupported CPU Poisson solver !!
diff --git a/src/SelfGravity/CPU_Poisson/CPU_PoissonSolver_FFT.cpp b/src/SelfGravity/CPU_Poisson/CPU_PoissonSolver_FFT.cpp
index 76be56e434..346670f59e 100644
--- a/src/SelfGravity/CPU_Poisson/CPU_PoissonSolver_FFT.cpp
+++ b/src/SelfGravity/CPU_Poisson/CPU_PoissonSolver_FFT.cpp
@@ -2,10 +2,10 @@
#if ( defined GRAVITY && defined SUPPORT_FFTW )
-static void FFT_Periodic( real *RhoK, const real Poi_Coeff, const int j_start, const int dj, const long RhoK_Size );
-static void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const long RhoK_Size );
+static void FFT_Periodic( real *RhoK, const real Poi_Coeff, const int j_start, const int dj, const long RhoK_Size, const int lv );
+static void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const long RhoK_Size, const int lv );
-extern root_fftw::real_plan_nd FFTW_Plan_Poi, FFTW_Plan_Poi_Inv;
+extern root_fftw::real_plan_nd FFTW_Plan_Poi[NLEVEL], FFTW_Plan_Poi_Inv[NLEVEL];
@@ -22,21 +22,22 @@ extern root_fftw::real_plan_nd FFTW_Plan_Poi, FFTW_Plan_Poi_Inv;
// j_start : Starting j index
// dj : Size of array in the j (y) direction after the forward FFT
// RhoK_Size : Size of the array "RhoK"
+// lv : Target level
//-------------------------------------------------------------------------------------------------------
-void FFT_Periodic( real *RhoK, const real Poi_Coeff, const int j_start, const int dj, const long RhoK_Size )
+void FFT_Periodic( real *RhoK, const real Poi_Coeff, const int j_start, const int dj, const long RhoK_Size, const int lv )
{
- const int Nx = NX0_TOT[0];
- const int Ny = NX0_TOT[1];
- const int Nz = NX0_TOT[2];
- const int Nx_Padded = Nx/2 + 1;
- const real dh = amr->dh[0];
+ const int Nx = NX0_TOT[0]*(1L<dh[lv];
real Deno;
gamer_fftw::fft_complex *cdata;
// forward FFT
- root_fftw_r2c( FFTW_Plan_Poi, RhoK );
+ root_fftw_r2c( FFTW_Plan_Poi[lv], RhoK );
// the data are now complex, so typecast a pointer
cdata = (gamer_fftw::fft_complex*) RhoK;
@@ -102,7 +103,7 @@ void FFT_Periodic( real *RhoK, const real Poi_Coeff, const int j_start, const in
// backward FFT
- root_fftw_c2r( FFTW_Plan_Poi_Inv, RhoK );
+ root_fftw_c2r( FFTW_Plan_Poi_Inv[lv], RhoK );
// normalization
const real norm = dh*dh / ( (real)Nx*Ny*Nz );
@@ -124,8 +125,9 @@ void FFT_Periodic( real *RhoK, const real Poi_Coeff, const int j_start, const in
// Parameter : RhoK : Array storing the input density and output potential
// Poi_Coeff : Coefficient in front of density in the Poisson equation (4*Pi*Newton_G*a)
// RhoK_Size : Size of the array "RhoK"
+// lv : Target level
//-------------------------------------------------------------------------------------------------------
-void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const long RhoK_Size )
+void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const long RhoK_Size, const int lv )
{
gamer_fftw::fft_complex *RhoK_cplx = (gamer_fftw::fft_complex *)RhoK;
@@ -134,7 +136,7 @@ void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const l
// forward FFT
- root_fftw_r2c( FFTW_Plan_Poi, RhoK );
+ root_fftw_r2c( FFTW_Plan_Poi[lv], RhoK );
// multiply density and Green's function in the k space
@@ -151,7 +153,7 @@ void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const l
// backward FFT
- root_fftw_c2r( FFTW_Plan_Poi_Inv, RhoK );
+ root_fftw_c2r( FFTW_Plan_Poi_Inv[lv], RhoK );
// effect of "4*PI*NEWTON_G" has been included in gFuncK, but the scale factor in the comoving frame hasn't
# ifdef COMOVING
@@ -163,6 +165,7 @@ void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const l
} // FUNCTION : FFT_Isolated
+
//-------------------------------------------------------------------------------------------------------
// Function : CPU_PoissonSolver_FFT
// Description : Evaluate the base-level potential by FFT
@@ -172,12 +175,13 @@ void FFT_Isolated( real *RhoK, const real *gFuncK, const real Poi_Coeff, const l
// Parameter : Poi_Coeff : Coefficient in front of the RHS in the Poisson eq.
// SaveSg : Sandglass to store the updated data
// PrepTime : Physical time for preparing the density field
+// lv : Target level
//-------------------------------------------------------------------------------------------------------
-void CPU_PoissonSolver_FFT( const real Poi_Coeff, const int SaveSg, const double PrepTime )
+void CPU_PoissonSolver_FFT( const real Poi_Coeff, const int SaveSg, const double PrepTime, const int lv )
{
// determine the FFT size (the zero-padding method is adopted for the isolated BC)
- int FFT_Size[3] = { NX0_TOT[0], NX0_TOT[1], NX0_TOT[2] };
+ int FFT_Size[3] = { NX0_TOT[0]*(1L<NPatchComma[0][1]*CUBE(PS1) ]; // MPI send buffer for density and potential
- real *RecvBuf = new real [ (long)NX0_TOT[0]*NX0_TOT[1]*NRecvSlice ]; // MPI recv buffer for density and potentia
- long *SendBuf_SIdx = new long [ (long)amr->NPatchComma[0][1]*PS1 ]; // MPI send buffer for 1D coordinate in slab
- long *RecvBuf_SIdx = new long [ (long)NX0_TOT[0]*NX0_TOT[1]*NRecvSlice/SQR(PS1) ]; // MPI recv buffer for 1D coordinate in slab
+ real *RhoK = (real*)root_fftw::fft_malloc( sizeof(real)*total_local_size ); // array storing both density and potential
+ real *SendBuf = new real [ (long)amr->NPatchComma[lv][1]*CUBE(PS1) ]; // MPI send buffer for density and potential
+ real *RecvBuf = new real [ (long)NX0_TOT[0]*(1L<NPatchComma[lv][1]*PS1 ]; // MPI send buffer for 1D coordinate in slab
+ long *RecvBuf_SIdx = new long [ (long)NX0_TOT[0]*(1L<PotSg [lv] = SaveSg_Pot;
@@ -152,24 +161,57 @@ void Gra_AdvanceDt( const int lv, const double TimeNew, const double TimeOld, co
false, false ),
Timer_Gra_Advance[lv], Timing );
- amr->FluSg[0] = SaveSg_Flu;
+ amr->FluSg[lv] = SaveSg_Flu;
} // if ( Gravity )
- } // if ( lv == 0 )
+ } // if ( FullRefinedLv )
- else // lv > 0
+ else // if ( FullRefinedLv )
{
if ( Poisson && !Gravity )
+ {
+# if ( POT_SCHEME == HYPRE_POI )
+ Hypre_Solver( HYPRE_SOLVE_TYPE_POISSON, lv, TimeNew, TimeOld, dt, Poi_Coeff, SaveSg_Flu, NULL_INT, SaveSg_Pot );
+ TIMING_FUNC( Buf_GetBufferData( lv, NULL_INT, NULL_INT, SaveSg_Pot, POT_FOR_POISSON, _POTE, _NONE, Pot_ParaBuf, USELB_YES ),
+ Timer_GetBuf[lv][1], Timing );
+// must call Poi_StorePotWithGhostZone AFTER collecting potential for buffer patches
+# ifdef STORE_POT_GHOST
+ TIMING_FUNC( Poi_StorePotWithGhostZone( lv, SaveSg_Pot, true ),
+ Timer_Gra_Advance[lv], Timing );
+# endif
+# else
InvokeSolver( POISSON_SOLVER, lv, TimeNew, TimeOld, NULL_REAL, Poi_Coeff, NULL_INT, NULL_INT, SaveSg_Pot,
OverlapMPI, Overlap_Sync );
+# endif
+ }
else if ( !Poisson && Gravity )
+ {
InvokeSolver( GRAVITY_SOLVER, lv, TimeNew, TimeOld, dt, NULL_REAL, SaveSg_Flu, NULL_INT, NULL_INT,
OverlapMPI, Overlap_Sync );
-
+ }
else if ( Poisson && Gravity )
+ {
+# if ( POT_SCHEME == HYPRE_POI )
+ Hypre_Solver( HYPRE_SOLVE_TYPE_POISSON, lv, TimeNew, TimeOld, dt, Poi_Coeff, SaveSg_Flu, NULL_INT, SaveSg_Pot );
+ amr->PotSg [lv] = SaveSg_Pot;
+ amr->PotSgTime[lv][SaveSg_Pot] = TimeNew;
+ TIMING_FUNC( Buf_GetBufferData( lv, NULL_INT, NULL_INT, SaveSg_Pot, POT_FOR_POISSON, _POTE, _NONE, Pot_ParaBuf, USELB_YES ),
+ Timer_GetBuf[lv][1], Timing );
+// must call Poi_StorePotWithGhostZone AFTER collecting potential for buffer patches
+# ifdef STORE_POT_GHOST
+ TIMING_FUNC( Poi_StorePotWithGhostZone( lv, SaveSg_Pot, true ),
+ Timer_Gra_Advance[lv], Timing );
+# endif
+
+
+ InvokeSolver( GRAVITY_SOLVER, lv, TimeNew, TimeOld, dt, NULL_REAL, SaveSg_Flu, NULL_INT, NULL_INT,
+ OverlapMPI, Overlap_Sync );
+# else
InvokeSolver( POISSON_AND_GRAVITY_SOLVER, lv, TimeNew, TimeOld, dt, Poi_Coeff, SaveSg_Flu, NULL_INT, SaveSg_Pot,
OverlapMPI, Overlap_Sync );
+# endif
+ }
}
diff --git a/src/SelfGravity/Init_GreenFuncK.cpp b/src/SelfGravity/Init_GreenFuncK.cpp
index fdd8465c7f..f814341c22 100644
--- a/src/SelfGravity/Init_GreenFuncK.cpp
+++ b/src/SelfGravity/Init_GreenFuncK.cpp
@@ -2,7 +2,7 @@
#if ( defined GRAVITY && defined SUPPORT_FFTW )
-extern root_fftw::real_plan_nd FFTW_Plan_Poi;
+extern root_fftw::real_plan_nd FFTW_Plan_Poi[NLEVEL];
@@ -15,9 +15,9 @@ extern root_fftw::real_plan_nd FFTW_Plan_Poi;
// 2. The zero-padding method is implemented
// 3. Slab decomposition is assumed in FFTW
//
-// Parameter : None
+// Parameter : lv : Target level
//-------------------------------------------------------------------------------------------------------
-void Init_GreenFuncK()
+void Init_GreenFuncK( const int lv )
{
if ( MPI_Rank == 0 ) Aux_Message( stdout, "%s ...\n", __FUNCTION__ );
@@ -27,9 +27,11 @@ void Init_GreenFuncK()
if ( OPT__BC_POT != BC_POT_ISOLATED )
Aux_Message( stderr, "OPT__BC_POT != BC_POT_ISOLATED, why do you need to calculate the Green's function !?\n" );
+ if ( ! FFTW_Inited[lv] ) Aux_Error( ERROR_INFO, "FFTW not initialized lv = %02d !!\n", lv );
+ if ( GreenFuncK_Inited[lv] ) return;
// 1. get the array indices used by FFTW
- const int FFT_Size[3] = { 2*NX0_TOT[0], 2*NX0_TOT[1], 2*NX0_TOT[2] };
+ const int FFT_Size[3] = { 2*NX0_TOT[0]*(1L<dh[0];
- const double Coeff = -NEWTON_G*CUBE(dh0)/( (double)FFT_Size[0]*FFT_Size[1]*FFT_Size[2] );
+ const double dh = amr->dh[lv];
+ const double Coeff = -NEWTON_G*CUBE(dh)/( (double)FFT_Size[0]*FFT_Size[1]*FFT_Size[2] );
double x, y, z, r;
int kk;
long idx;
- GreenFuncK = (real*) root_fftw::fft_malloc(sizeof(real) * total_local_size);
+ GreenFuncK[lv] = (real*) root_fftw::fft_malloc(sizeof(real) * total_local_size);
for (int k=0; kFluSg[lv], NULL_INT, NULL_INT, DATA_GENERAL, _DENS, _NONE,
Rho_ParaBuf, USELB_YES );
@@ -388,8 +391,12 @@ void Aux_Record_Gravity()
if ( lv >= MinLv ) Timer_PoiPerf.Start();
for (int t=0; tPotSg[lv], Time[lv] );
+ if ( FullRefinedLv )
+ {
+ if ( ! FFTW_Inited[lv] ) Init_FFTW( lv );
+
+ CPU_PoissonSolver_FFT( Poi_Coeff, amr->PotSg[lv], Time[lv], lv );
+ }
else
InvokeSolver( POISSON_SOLVER, lv, Time[lv], NULL_REAL, NULL_REAL, Poi_Coeff,
diff --git a/src/TestProblem/Hydro/JeansInstability/Init_TestProb_Hydro_JeansInstability.cpp b/src/TestProblem/Hydro/JeansInstability/Init_TestProb_Hydro_JeansInstability.cpp
index 515e2b8a39..ea44bd00d8 100644
--- a/src/TestProblem/Hydro/JeansInstability/Init_TestProb_Hydro_JeansInstability.cpp
+++ b/src/TestProblem/Hydro/JeansInstability/Init_TestProb_Hydro_JeansInstability.cpp
@@ -23,6 +23,7 @@ static double Jeans_WaveSpeed; // propagation speed (sound speed in hydro o
static bool Jeans_Stable; // true/false --> Jeans stable/unstable
// =======================================================================================
+extern int Evolve_stage;
@@ -470,6 +471,59 @@ static void OutputError()
# endif
} // FUNCTION : OutputError
+
+
+
+//-------------------------------------------------------------------------------------------------------
+// Function : Mis_UserWorkBeforeNextSubstep_JeansInstability
+// Description : Template of user-specified work before proceeding to the next sub-step in EvolveLevel()
+// --> After fix-up and grid refinement on lv
+//
+// Note : 1. Invoked by EvolveLevel() using the function pointer "Mis_UserWorkBeforeNextSubstep_Ptr"
+//
+// Parameter : lv : Target refinement level
+// TimeNew : Target physical time to reach
+// TimeOld : Physical time before update
+// dt : Time interval to advance solution (can be different from TimeNew-TimeOld in COMOVING)
+//-------------------------------------------------------------------------------------------------------
+void Mis_UserWorkBeforeNextSubstep_JeansInstability( const int lv, const double TimeNew, const double TimeOld, const double dt )
+{
+
+ return;
+
+ for (int r=0; rNPatchComma[lv][1]; PID++)
+ {
+ for (int k=0; kpatch[0][lv][PID]->cornerL[0] + i,
+ amr->patch[0][lv][PID]->cornerL[1] + j,
+ amr->patch[0][lv][PID]->cornerL[2] + k,
+ amr->patch[0][lv][PID]->EdgeL[0],
+ amr->patch[0][lv][PID]->EdgeL[1],
+ amr->patch[0][lv][PID]->EdgeL[2],
+ amr->patch[ amr->FluSg[lv] ][lv][PID]->fluid[DENS][k][j][i],
+ amr->patch[ amr->FluSg[lv] ][lv][PID]->fluid[MOMX][k][j][i],
+ amr->patch[ amr->FluSg[lv] ][lv][PID]->fluid[MOMY][k][j][i],
+ amr->patch[ amr->FluSg[lv] ][lv][PID]->fluid[MOMZ][k][j][i],
+ amr->patch[ amr->FluSg[lv] ][lv][PID]->fluid[ENGY][k][j][i],
+ amr->patch[ amr->PotSg[lv] ][lv][PID]->pot[k][j][i]
+ );
+ }
+ }
+ }
+ MPI_Barrier( MPI_COMM_WORLD );
+ }
+
+} // FUNCTION : Mis_UserWorkBeforeNextSubstep_JeansInstability
#endif // #if ( MODEL == HYDRO && defined GRAVITY )
@@ -508,6 +562,7 @@ void Init_TestProb_Hydro_JeansInstability()
# ifdef SUPPORT_HDF5
Output_HDF5_InputTest_Ptr = LoadInputTestProb;
# endif
+ Mis_UserWorkBeforeNextSubstep_Ptr = Mis_UserWorkBeforeNextSubstep_JeansInstability;
# endif // if ( MODEL == HYDRO && defined GRAVITY )
diff --git a/src/configure.py b/src/configure.py
index 4a8a8e215e..66901e60d5 100755
--- a/src/configure.py
+++ b/src/configure.py
@@ -601,8 +601,9 @@ def load_arguments( sys_setting : SystemSetting ):
)
parser.add_argument( "--pot_scheme", type=str, metavar="TYPE", gamer_name="POT_SCHEME",
- default="SOR", choices=["SOR", "MG"],
+ default="SOR", choices=["SOR", "MG", "HYPRE_POI"],
depend={"gravity":True},
+ constraint={ "HYPRE_POI":{"hypre":True} },
help="Select the Poisson solver. SOR: successive-overrelaxation (recommended), MG: multigrid. "\
"Must be set when <--gravity> is enabled.\n"
)
@@ -612,6 +613,7 @@ def load_arguments( sys_setting : SystemSetting ):
depend={"gravity":True},
help="Store GRA_GHOST_SIZE ghost-zone potential for each patch on each side. "\
"Recommended when PARTICLE is enabled for improving accuaracy for particles around the patch boundaries. "\
+ "Useless for <--pot_scheme=HYPRE_POI>. "\
"Must be enabled for <--star_formation> + <--store_par_acc>.\n"
)
@@ -806,6 +808,11 @@ def load_arguments( sys_setting : SystemSetting ):
"Must compile libyt with JUPYTER_KERNEL. Must enable <--libyt>.\n"
)
+ parser.add_argument( "--hypre", type=str2bool, metavar="BOOLEAN", gamer_name="SUPPORT_HYPRE",
+ default=False,
+ help="Support HYPRE library.\n"
+ )
+
parser.add_argument( "--rng", type=str, metavar="TYPE", gamer_name="RANDOM_NUMBER",
default=None,
choices=["RNG_GNU_EXT", "RNG_CPP11"],
@@ -1109,7 +1116,7 @@ def warning( paths, **kwargs ):
# 3. Path
path_links = { "gpu":{True:"CUDA_PATH"}, "fftw":{"FFTW2":"FFTW2_PATH", "FFTW3":"FFTW3_PATH"},
"mpi":{True:"MPI_PATH"}, "hdf5":{True:"HDF5_PATH"}, "grackle":{True:"GRACKLE_PATH"},
- "gsl":{True:"GSL_PATH"}, "libyt":{True:"LIBYT_PATH"} }
+ "gsl":{True:"GSL_PATH"}, "libyt":{True:"LIBYT_PATH"}, "hypre":{True:"HYPRE_PATH"} }
for arg, links in path_links.items():
for val, p_name in links.items():
diff --git a/tool/vscode/README.md b/tool/vscode/README.md
new file mode 100644
index 0000000000..760afad070
--- /dev/null
+++ b/tool/vscode/README.md
@@ -0,0 +1,20 @@
+# Integrating GAMER with Visual Studio Code
+
+Run `sh tool/vscode/copy_to_vscode.sh` to copy necessary files to the `.vscode` directory to integrate GAMER with VS Code.
+
+Please check the [wiki](https://github.com/gamer-project/gamer/wiki/Developing-with-VS-Code) for more information.
+
+## Files in this directory
+- `c_cpp_properties.json`: Contains the IntelliSense configuration for VS Code.
+- `gamercpp.natvis`: Contains the visualization configuration for VS Code debugger.
+- `launch.json`: Contains the debug configuration for VS Code.
+- `settings.json`: Contains the settings for the editor in VS Code.
+- `tasks.json`: Contains the build configuration for VS Code.
+- `extract_macros.py`: Script to extract the macros from the `Makefile.log` to `c_cpp_properties.json`.
+- `copy_to_vscode.sh`: Script to copy the above files to `.vscode` directory.
+- `bin_working`: File for storing the name of the working directory under `bin/`.
+- `set_bin_working.sh`: Script to set up the path of the input files and the executable.
+- `build.sh`: Script for the build task.
+- `clean_work_dir.sh`: Script for the clean task.
+- `config.sh`: Script for the configuration task.
+- `README.md`: This file.
diff --git a/tool/vscode/bin_working b/tool/vscode/bin_working
new file mode 100644
index 0000000000..e69de29bb2
diff --git a/tool/vscode/build.sh b/tool/vscode/build.sh
new file mode 100644
index 0000000000..eb6a4d603d
--- /dev/null
+++ b/tool/vscode/build.sh
@@ -0,0 +1,19 @@
+PYTHON=python3
+
+bin_working=$(cat bin_working)
+if [ -z "${bin_working}" ]; then
+ echo "Error: Please set the working directory by the task 'set-working-bin' first."
+ exit 1
+fi
+
+cd ../src
+make clean
+
+cd ../.vscode
+${PYTHON} extract_macros.py
+
+cd ../src
+make -j 8
+
+cd ../bin/${bin_working}/
+cp ../gamer .
diff --git a/tool/vscode/c_cpp_properties.json b/tool/vscode/c_cpp_properties.json
new file mode 100644
index 0000000000..cfdaaa577f
--- /dev/null
+++ b/tool/vscode/c_cpp_properties.json
@@ -0,0 +1,15 @@
+{
+ "configurations": [
+ {
+ "name": "Linux",
+ "includePath": [
+ "${default}",
+ "${workspaceFolder}/include"
+ ],
+ "cStandard": "c17",
+ "cppStandard": "gnu++17",
+ "intelliSenseMode": "linux-gcc-x64"
+ }
+ ],
+ "version": 4
+}
\ No newline at end of file
diff --git a/tool/vscode/clean_work_dir.sh b/tool/vscode/clean_work_dir.sh
new file mode 100644
index 0000000000..83b71dd818
--- /dev/null
+++ b/tool/vscode/clean_work_dir.sh
@@ -0,0 +1,8 @@
+bin_working=$(cat bin_working)
+if [ -z "${bin_working}" ]; then
+ echo "Error: Please set the working directory by the task 'set-working-bin' first."
+ exit 1
+fi
+
+cd ../bin/${bin_working}/
+sh clean.sh
\ No newline at end of file
diff --git a/tool/vscode/config.sh b/tool/vscode/config.sh
new file mode 100644
index 0000000000..d3fcf06175
--- /dev/null
+++ b/tool/vscode/config.sh
@@ -0,0 +1,9 @@
+bin_working=$(cat bin_working)
+if [ -z "${bin_working}" ]; then
+ echo "Error: Please set the working directory by the task 'set-working-bin' first."
+ exit 1
+fi
+
+cd ../src
+cp ../bin/${bin_working}/generate_make.sh .
+sh generate_make.sh
\ No newline at end of file
diff --git a/tool/vscode/copy_to_vscode.sh b/tool/vscode/copy_to_vscode.sh
new file mode 100644
index 0000000000..8791289adc
--- /dev/null
+++ b/tool/vscode/copy_to_vscode.sh
@@ -0,0 +1,51 @@
+#!/bin/bash
+
+# Change to the directory where the script is located
+cd "$(dirname "$0")"
+
+# Set source and destination paths relative to the script's location
+SOURCE_DIR="."
+TARGET_DIR="../../.vscode"
+
+# Specify files to exclude (space-separated list, e.g., "file1 file2")
+EXCLUDE_FILES="copy_to_vscode.sh README.md"
+
+if [ ! -d "$TARGET_DIR" ]; then
+ echo "Target directory $TARGET_DIR does not exist."
+ mkdir -p "$TARGET_DIR"
+ echo "Target directory $TARGET_DIR created."
+fi
+
+# Function to check if a file is in the exclude list
+is_excluded() {
+ for excluded in $EXCLUDE_FILES; do
+ if [ "$excluded" = "$1" ]; then
+ return 0
+ fi
+ done
+ return 1
+}
+
+for file in "$SOURCE_DIR"/*; do
+ filename=$(basename "$file")
+
+ if is_excluded "$filename"; then
+ continue
+ fi
+
+ if [ -e "$TARGET_DIR/$filename" ]; then
+ echo "File $filename already exists in $TARGET_DIR. Overwrite? (y/n)"
+ read -r answer
+
+ if [ "$answer" = "y" ] || [ "$answer" = "Y" ]; then
+ cp "$file" "$TARGET_DIR/$filename"
+ echo "$filename overwritten."
+ else
+ echo "$filename not copied."
+ fi
+ else
+ cp "$file" "$TARGET_DIR/$filename"
+ echo "$filename copied."
+ fi
+
+done
diff --git a/tool/vscode/extract_macros.py b/tool/vscode/extract_macros.py
new file mode 100644
index 0000000000..196631ee25
--- /dev/null
+++ b/tool/vscode/extract_macros.py
@@ -0,0 +1,56 @@
+import re
+import json
+
+'''
+This script updates the "defines" section in the .vscode/c_cpp_properties.json file in Visual Studio Code
+by directly parsing the Makefile.log file to extract relevant compiler macros (definitions). It provides an automated
+way to synchronize project-specific preprocessor definitions with the VSCode configuration. This approach allows
+VSCode to recognize the defines directly from your build configuration, improving VSCode IntelliSense and error
+detection without manual intervention.
+
+Usage:
+1. Place this script under the `.vscode` folder of your project.
+2. Set the paths to `Makefile.log` and `c_cpp_properties.json` in the script as follows:
+ - `makefile_log_path`: Path to `Makefile.log`, which contains compiler settings output (default: "../src/Makefile.log").
+ - `c_cpp_properties_path`: Path to the VSCode C++ configuration file (default: "c_cpp_properties.json").
+3. Run the script each time `Makefile.log` is updated, or set it as a pre-build or post-build task in VSCode
+ to keep the configuration in sync.
+'''
+# Path to Makefile.log and c_cpp_properties.json
+makefile_log_path = "../src/Makefile.log"
+c_cpp_properties_path = "c_cpp_properties.json"
+
+# Pattern to match the setting in the format of " MODEL : HYDRO"
+pattern = re.compile(r":\s+(\w+)\s*:\s+(\w+)")
+defines = []
+
+# Read Makefile.log and extract macros
+with open(makefile_log_path, 'r') as log_file:
+ print(f"Reading {makefile_log_path} and extracting defines...")
+ for line in log_file:
+ match = pattern.search(line)
+ if match:
+ key, value = match.groups()
+ if value == 'False':
+ continue
+ elif value == 'True':
+ defines.append(f"{key}")
+ else:
+ defines.append(f"{key}={value}")
+print(f"Extracted {len(defines)} macros from {makefile_log_path}.")
+
+# Load c_cpp_properties.json
+with open(c_cpp_properties_path, 'r') as cpp_properties_file:
+ print(f"Loading {c_cpp_properties_path}...")
+ cpp_properties = json.load(cpp_properties_file)
+
+# Update the "defines" array in c_cpp_properties.json
+print("Updating defines in c_cpp_properties.json...")
+cpp_properties['configurations'][0]['defines'] = defines
+
+# Write the updated content back to c_cpp_properties.json
+with open(c_cpp_properties_path, 'w') as cpp_properties_file:
+ print(f"Writing updates to {c_cpp_properties_path}...")
+ json.dump(cpp_properties, cpp_properties_file, indent=4)
+
+print("Update completed successfully.")
diff --git a/tool/vscode/gamercpp.natvis b/tool/vscode/gamercpp.natvis
new file mode 100644
index 0000000000..81dd2d3fd4
--- /dev/null
+++ b/tool/vscode/gamercpp.natvis
@@ -0,0 +1,7 @@
+
+
+
+ Timer running
+ Timer stop: {Time*1.0e-6} sec
+
+
diff --git a/tool/vscode/launch.json b/tool/vscode/launch.json
new file mode 100644
index 0000000000..c77de7fd65
--- /dev/null
+++ b/tool/vscode/launch.json
@@ -0,0 +1,40 @@
+{
+ "version": "0.2.0",
+ "configurations": [
+ {
+ "name": "Debug GAMER",
+ "type": "cppdbg",
+ "request": "launch",
+ "program": "${workspaceFolder}/bin/${input:bin-working}/gamer",
+ "args": [],
+ "stopAtEntry": true,
+ "cwd": "${workspaceFolder}/bin/${input:bin-working}/",
+ "environment": [],
+ "externalConsole": false,
+ "MIMode": "gdb",
+ "setupCommands": [
+ {
+ "description": "Enable pretty-printing for gdb",
+ "text": "-enable-pretty-printing",
+ "ignoreFailures": true
+ }
+ ],
+ "miDebuggerPath": "/usr/bin/gdb",
+ "miDebuggerArgs": "-quiet",
+ "logging": {
+ "engineLogging": false
+ },
+ "internalConsoleOptions": "openOnSessionStart",
+ "preLaunchTask": "clean-work-dir",
+ "visualizerFile": "${workspaceFolder}/.vscode/gamercpp.natvis",
+ "showDisplayString": true
+ }
+ ],
+ "inputs": [
+ {
+ "id": "bin-working",
+ "type": "promptString",
+ "description": "Enter the working directory under bin/."
+ }
+ ]
+}
\ No newline at end of file
diff --git a/tool/vscode/set_bin_working.sh b/tool/vscode/set_bin_working.sh
new file mode 100644
index 0000000000..4e08330f81
--- /dev/null
+++ b/tool/vscode/set_bin_working.sh
@@ -0,0 +1,18 @@
+#!/bin/bash
+
+cd "$(dirname "$0")"
+
+# If $1 is given, use it as $chosen, else prompt the user
+if [ -n "$1" ]; then
+ chosen="$1"
+else
+ echo "Directories under bin/:"
+ for dir in ../bin/*/; do
+ echo "$(basename "$dir")"
+ done
+ echo -n "Enter a directory name from the above list: "
+ read chosen
+fi
+
+printf "%s" "$chosen" > bin_working
+echo "bin_working updated."
\ No newline at end of file
diff --git a/tool/vscode/settings.json b/tool/vscode/settings.json
new file mode 100644
index 0000000000..c80f4cee33
--- /dev/null
+++ b/tool/vscode/settings.json
@@ -0,0 +1,11 @@
+{
+ "C_Cpp.formatting": "disabled",
+ "files.associations": {
+ "Input__*": "ini",
+ "Record__*": "markdown"
+ },
+ "[cpp]": {
+ "editor.detectIndentation": false,
+ "editor.tabSize": 3
+ }
+}
\ No newline at end of file
diff --git a/tool/vscode/tasks.json b/tool/vscode/tasks.json
new file mode 100644
index 0000000000..d8ffbb0369
--- /dev/null
+++ b/tool/vscode/tasks.json
@@ -0,0 +1,92 @@
+{
+ "tasks": [
+ {
+ "label": "build-GAMER",
+ "type": "shell",
+ "command": "sh build.sh",
+ "group": {
+ "kind": "build",
+ "isDefault": false
+ },
+ "options": {
+ "cwd": "${workspaceFolder}/.vscode"
+ },
+ "detail": "After you configure the Makefile, run this task to build the project."
+ },
+ {
+ "label": "config-GAMER",
+ "type": "shell",
+ "command": "sh config.sh",
+ "options": {
+ "cwd": "${workspaceFolder}/.vscode"
+ },
+ "detail": "Run configure.py with generate_make.sh in the working directory.",
+ "problemMatcher": []
+ },
+ {
+ "label": "clean-work-dir",
+ "type": "shell",
+ "command": "sh clean_work_dir.sh",
+ "options": {
+ "cwd": "${workspaceFolder}/.vscode"
+ },
+ "detail": "Clean the working directory with clean.sh.",
+ "problemMatcher": []
+ },
+ {
+ "label": "set-working-bin",
+ "type": "shell",
+ "command": "sh",
+ "args": [
+ "set_bin_working.sh",
+ "${input:bin-working}"
+ ],
+ "options": {
+ "cwd": "${workspaceFolder}/.vscode"
+ },
+ "detail": "Choose the working directory for the binary files.",
+ "problemMatcher": []
+ },
+ {
+ "label": "updated_mac_launch",
+ "type": "shell",
+ "command": "sh",
+ "args": [
+ "updated_mac_launch.sh",
+ "${input:lldb-mi_path}"
+ ],
+ "options": {
+ "cwd": "${workspaceFolder}/.vscode"
+ },
+ "detail": "Update launch.json with the path to your lldb-mi executable for mac user.",
+ "problemMatcher": []
+ },
+ {
+ "label": "config-and-build",
+ "dependsOn": [
+ "config-GAMER",
+ "build-GAMER"
+ ],
+ "dependsOrder": "sequence",
+ "group": {
+ "kind": "build",
+ "isDefault": true
+ },
+ "detail": "Run config-GAMER and build-GAMER.",
+ "problemMatcher": []
+ }
+ ],
+ "inputs": [
+ {
+ "id": "bin-working",
+ "type": "promptString",
+ "description": "Enter the working directory under bin/."
+ },
+ {
+ "id": "lldb-mi_path",
+ "type": "promptString",
+ "description": "Enter the path to your lldb-mi executable."
+ }
+ ],
+ "version": "2.0.0"
+}
\ No newline at end of file
diff --git a/tool/vscode/updated_mac_launch.sh b/tool/vscode/updated_mac_launch.sh
new file mode 100644
index 0000000000..f610f95fb3
--- /dev/null
+++ b/tool/vscode/updated_mac_launch.sh
@@ -0,0 +1,50 @@
+#!/bin/bash
+
+# Check OS
+if [ "$(uname)" == "Darwin" ]; then
+ echo "OS: macOS. Start updating launch.json..."
+else
+ echo "OS: Linux. Do nothing."
+ exit 0
+fi
+
+cd "$(dirname "$0")"
+
+# If $1 is given, use it as $lldb_mi_path, else prompt the user
+if [ -n "$1" ]; then
+ lldb_mi_path="$1"
+else
+ echo -n "Enter your lldb-mi path. (e.g. /Users/user_namer/lldb-mi/build/src/lldb-mi): "
+ read lldb_mi_path
+fi
+
+# Parameters to change
+MIMODE="lldb"
+MIDEBUGGERPATH=$lldb_mi_path
+MIDEBUGGERARGS="-q"
+
+# Input file
+INPUT_FILE="launch.json"
+
+# Replace MIMode from "gdb" to "lldb"
+if grep -q "\"MIMode\": \"$MIMODE\"" "$INPUT_FILE"; then
+ echo "Warning: 'MIMode' is already set to '$MIMODE'."
+else
+ sed -i '' 's/"MIMode": "gdb"/"MIMode": "'"$MIMODE"'"/' "$INPUT_FILE"
+ echo "'MIMode' changed to '$MIMODE'."
+fi
+
+# Replace miDebuggerPath to the new path
+sed -i '' 's|"miDebuggerPath": .*|"miDebuggerPath": "'"$MIDEBUGGERPATH"'",|' "$INPUT_FILE"
+echo "'miDebuggerPath' changed to '$MIDEBUGGERPATH'."
+
+
+# Replace miDebuggerArgs to "-q"
+if grep -q "\"miDebuggerArgs\": \"$MIDEBUGGERARGS\"" "$INPUT_FILE"; then
+ echo "Warning: 'miDebuggerArgs' is already set to '$MIDEBUGGERARGS'."
+else
+ sed -i '' 's/"miDebuggerArgs": "-quiet"/"miDebuggerArgs": "'"$MIDEBUGGERARGS"'"/' "$INPUT_FILE"
+ echo "'miDebuggerArgs' changed to '$MIDEBUGGERARGS'."
+fi
+
+echo "launch.json updated successfully."
\ No newline at end of file
diff --git a/tool/wiki/sync_runtime_parameter.py b/tool/wiki/sync_runtime_parameter.py
index d5d7f3eb91..ffc8e3f6e9 100644
--- a/tool/wiki/sync_runtime_parameter.py
+++ b/tool/wiki/sync_runtime_parameter.py
@@ -50,6 +50,7 @@
"../../doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-General.md",
"../../doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Gravity.md",
"../../doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Hydro.md",
+ "../../doc/wiki/Runtime-Parameters-related/[Runtime-Parameters]-Hypre.md",
"../../doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Initial-Conditions.md",
"../../doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Interpolation.md",
"../../doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-MPI-and-OpenMP.md",