CompatHelper: bump compat for CairoMakie to 0.15, (keep existing compat)

.github/pull_request_template.md

@@ -0,0 +1,16 @@
+*We set up this PR template to confirm that you have considered all the necessary aspects of your contribution so that we are on the right track. Please fill out the template below to the best of your ability. If you have any questions, please feel free to get in touch.*
+
+## Proposed Changes
+*Give a brief overview of your contribution here in a few sentences. For bug fixes, please detail the bug and include a test case. For a new feature, please clearly describe the design, its rationale, the possible alternatives considered. Note that each PR should address only one issue*
+
+## Related Work
+*Resolve any issues (bug fix or feature request), note any related PRs, or mention interactions with the work of others, if any.*
+
+## Checklist
+*Put an x by all that apply. You can fill this out after submitting the PR.*
+
+- [ ] I am submitting my contribution to the development branch.
+- [ ] My contribution is commented and uses the code style of the project (see .JuliaFormatter.toml).
+- [ ] I have added a unit test to verify my contribution, if necessary.
+- [ ] I have updated appropriate documentation (README, documentation, tutorials, etc.), if applicable.
+- [ ] I have added an entry in NEWS.md, if necessary.

.github/workflows/ci.yml

@@ -32,8 +32,7 @@ jobs:
         arch:
           - x64
     steps:
-      - run: |
-          pip3 install meshio scipy --user
+      - run: pip3 install scipy --user
       - uses: actions/checkout@v2
       - uses: julia-actions/setup-julia@v1
         with:
@@ -56,6 +55,8 @@ jobs:
       - uses: julia-actions/julia-runtest@latest
         continue-on-error: ${{ matrix.version == 'nightly' }}
       - uses: julia-actions/julia-processcoverage@v1
-      - uses: codecov/codecov-action@v1
+      - uses: codecov/codecov-action@v4
         with:
-          file: lcov.info
+          file: lcov.info
+          token: ${{ secrets.CODECOV_TOKEN }}
+          fail_ci_if_error: true

.github/workflows/format.yml

@@ -2,6 +2,7 @@ name: FormatPR
 on:
   schedule:
     - cron: '0 0 * * *'
+permissions: write-all
 jobs:
   build:
     runs-on: ubuntu-latest
@@ -28,4 +29,4 @@ jobs:
       - name: Check outputs
         run: |
           echo "Pull Request Number - ${{ steps.cpr.outputs.pull-request-number }}"
-          echo "Pull Request URL - ${{ steps.cpr.outputs.pull-request-url }}"
+          echo "Pull Request URL - ${{ steps.cpr.outputs.pull-request-url }}"

Project.toml

@@ -7,24 +7,39 @@ version = "0.1.0"
 AbstractTrees = "1520ce14-60c1-5f80-bbc7-55ef81b5835c"
 CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0"
 DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
-HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f"
-JLD = "4138dd39-2aa7-5051-a626-17a0bb65d9c8"
 JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819"
-JuliaFormatter = "98e50ef6-434e-11e9-1051-2b60c6c9e899"
-Libdl = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
-MAT = "23992714-dd62-5051-b70f-ba57cb901cac"
-MKL_jll = "856f044c-d86e-5d09-b602-aeab76dc8ba7"
 MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195"
-MPICH_jll = "7cb0a576-ebde-5e09-9194-50597f1243b4"
 P4est = "7d669430-f675-4ae7-b43e-fab78ec5a902"
 Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a"
 PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"
-ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca"
 PythonCall = "6099a3de-0909-46bc-b1f4-468b9a2dfc0d"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
-SHA = "ea8e919c-243c-51af-8825-aaa63cd721ce"
+Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
-Tar = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e"
+
+[compat]
+AbstractTrees = "0.4"
+CairoMakie = "0.12, 0.13, 0.15"
+DocStringExtensions = "0.9"
+JLD2 = "0.4, 0.5"
+LinearAlgebra = "<0.0.1, 1"
+MPI = "0.20"
+P4est = "0.4"
+Parameters = "0.12"
+PrecompileTools = "^1"
+PythonCall = "0.9"
+RecipesBase = "^1"
+Reexport = "^1"
+SparseArrays = "1.10"
+SpecialFunctions = "^2"
+StaticArrays = "1.9"
+julia = "^1"
+
+[extras]
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+
+[targets]
+test = ["Test"]

README.md

@@ -1,2 +1,111 @@
 # KitAMR.jl
-Large-scale distributed computational fluid dynamics computing with adaptive mesh refinement
+
+![CI](https://img.shields.io/github/actions/workflow/status/vavrines/KitBase.jl/ci.yml?branch=main)
+[![codecov](https://img.shields.io/codecov/c/github/CFDML/KitAMR.jl)](https://codecov.io/gh/CFDML/KitAMR.jl)
+
+## Motivation
+**KitAMR.jl** is a distributed adaptive Cartesian grid solver for kinetic equations, developed based on [P4est.jl](https://github.com/trixi-framework/P4est.jl). Its goal is to achieve large-scale parallel solving of 2D and 3D flows across all regimes, leveraging GPUs, differentiable programming, and machine learning to enhance solving efficiency.
+
+Currently, the development is in its initial stages, with the following functionalities implemented:
+- 2D flow
+- Square domain grid generation
+- Geometric adaptation of physical space grids
+- Dynamic adaptation of physical space grids
+- Dynamic adaptation of velocity space grids
+- Load balancing based on physical space grid density
+- VanLeer reconstruction
+- UGKS fluxes
+- Post-processing
+
+## Theory
+### Kinetic Methods
+In kinetic theory, methods simulate the macroscopic motion of fluids by describing the evolution of distribution functions of fluid particles in velocity space over time. Such methods include direct simulation Monte Carlo (DSMC), discrete velocity method (DVM), lattice Boltzmann equation (LBE), gas-kinetic scheme, semi-Lagrangian method, implicit-explicit (IMEX) method, and others.
+### Distribution Function
+In kinetic methods, the state of the fluid is described by a distribution function $f(\mathbf{x},\mathbf{u},t)$. The distribution function represents the number of particles at a given moment in time at a specific physical space point in a particular velocity space element. Its normalization property is
+
+$$
+\rho(\mathbf{x},t)=\iiint_{-\infty}^{\infty}f(\mathbf{x},\mathbf{u},t)d\mathbf{u}
+$$
+
+, where $\rho(\mathbf{x},t)$ is the density of flows at $\mathbf{x}$ in time $t$. At higher temperatures, molecular internal degrees of freedom are excited, and the distribution function will depend on these degrees of freedom as well. For simplicity, we do not currently consider this situation.
+
+<div align="center">
+    <img src="https://i.postimg.cc/pTsjgW3T/20240703151108.jpg)](https://postimg.cc/ygSWRH14" alt="Image 1" width="300" title="distribution function of 2-D" style="margin-right: 1px;">
+</div>
+With the distribution function, macroscopic conservative variables can be caculated as follows:
+
+$$
+\begin{split}
+&\rho = \iiint_{-\infty}^{\infty}f(\mathbf{x},\mathbf{u},t)d\mathbf{u}\\
+&\rho\mathbf U = \iiint_{-\infty}^{\infty}\mathbf uf(\mathbf{x},\mathbf{u},t)d\mathbf{u}\\
+&\rho E = \iiint_{-\infty}^{\infty} \frac 12 \mathbf u^2f(\mathbf{x},\mathbf{u},t)d\mathbf{u}\\
+&\mathbf q=\iiint_{-\infty}^{\infty} \mathbf u\frac 12 \mathbf u^2f(\mathbf{x},\mathbf{u},t)d\mathbf{u}
+\end{split}
+$$
+
+### Boltzmann Equation
+The evolution of the distribution function follows the Boltzmann equation.
+
+$$
+\frac{\partial f}{\partial t}+\mathbf{u}\cdot\frac{\partial f}{\partial \mathbf{x}}+\mathbf{F}\cdot\frac{\partial f}{\partial \mathbf{u}}=\iiint_{-\infty}^{\infty}\int_0^{4\pi}(f^{\*}f_1^{\*}-ff_1)c_r\sigma d\Omega dc_1
+$$
+
+<!-- The Boltzmann equation holds a central position in kinetic theory. In the free molecular flow regime, the collisionless Boltzmann equation is applied, or its equilibrium solution, the Maxwellian distribution, is directly used. In the slip flow regime, the Boltzmann equation is expanded via the Chapman-Enskog method to obtain the first-order approximation, the Navier-Stokes equations, or the second-order approximation, the Burnett equations. Alternatively, a systematic asymptotic expansion yields a set of fluid mechanics equations. In the transition regime, solving flow problems involves approaches derived from the Boltzmann equation or its equivalent methods. -->
+The Boltzmann equation holds a central position in kinetic theory, which is an integro-differential equation, with a four-fold integral term on the right-hand side. The complexity introduced by this term makes solving the Boltzmann equation challenging. Therefore, as a simplification, kinetic methods often solve its model equations instead.
+### Shakhov Model Equation
+The simplified equation obtained after simplifying the collision term is known as the Shakhov model equation:
+
+$$
+\begin{cases}
+\frac{\partial f}{\partial t}+\mathbf{u}\cdot\frac{\partial f}{\partial \mathbf{x}}+\mathbf{F}\cdot\frac{\partial f}{\partial \mathbf{u}}=\frac{g^+-f}{\tau}\\
+g^+=g[1+\frac45(1-\mathrm{Pr})\lambda^2\frac{\mathbf{u}\cdot\mathbf{q}}{\rho}(2\lambda \mathbf{u}^2-5)]
+\end{cases}
+$$
+
+, Where $\mathrm {Pr}$ is the Prandtl number, $\mathbf{q}$ is the heat flux, and $\tau$ is the relaxation time typically given through phenomenological models to match experimental data. The Shakhov model can provide results close to the Boltzmann equation for mass, momentum, energy, and heat fluxes in situations where flow velocities are not very high and deviations from equilibrium are not very large. And the current solver is based on Shakhov model.
+
+## Computing Methods
+### Finite Volume Method (FVM)
+The underlying discretization of the Boltzmann equation is based on finite volume method, with its discrete form as follows:
+
+$$
+\{\bar U_j\Omega_j\}^{n+1} = 
+\{\bar U_j \Omega_j\}^n-\Delta t\sum_{\partial\Omega}\mathbf{F}^{\*}\cdot\Delta\mathbf{S}+\Delta t \bar{Q}_j\Omega_j
+$$
+
+, where $\bar U_j^n\coloneqq\frac{1}{\Omega_j}\int_{\Omega_j}Ud\Omega_j|^n$, $\bar Q_j^n\coloneqq\frac{1}{\Omega_j}\int_{\Omega_j}Qd\Omega_j$ are the cell-averaged conservative variable and source. And $\bar Q_j^{\*}$ and $\mathbf{F}^{\*}$ are respectively cell- and time-averaged sources and numerical flux, which are defined as $\bar Q_j^{\*}\coloneqq\frac{1}{\Delta t}\int_n^{n+1}\bar Q_jdt$ and $\mathbf F^{\*}\cdot \Delta \mathbf{S}\coloneqq \frac{1}{\Delta t}\int_n^{n+1}\mathbf F\cdot\Delta \mathbf{S}dt$.
+It is an exact relation for the time evolution of the space averaged conservative avriables over cell $j$ from time step $n$ to $n+1$. And the numerical flux $\mathbf F^{\*}$ identifies completely a scheme by the way it approximates the time-averaged physical flux along each cell face.
+### Numerical Flux
+The numerical scheme that the solver currently supports is UGKS, with interface variables obtained from VanLeer reconstruction. For specific implementation details, please refer to the relevant literature.
+
+### Discrete Velocity Method
+DVM is one of the popular approaches for solving rarefied flow problems. In this method, the particle velocity space is discretized into a finite set of points and the numerical quadrature rule is utilized to approximate the integration of moments. Since the particle velocity space is discretized into a finite set of points, the continuum Boltzmann equation is reduced to the corresponding discrete velocity Boltzmann equation (DVBE). A great variety of algorithms have been developed, including unified gas-kinetic scheme (UGKS), discrete unified gas-kinetic scheme (DUGKS), semi-Lagrangian method, etc. Overall, most of the above methods can be applied from free molecular regime to continuum regime. But in order to make the quadrature error to be small enough, a large number of discrete velocity points are usually required. In particular, for fluid flows near continuum regime, the computational cost of DVM is much larger than those traditional CFD methods based on the Navier-Stokes equation.
+
+### Adaptive Mesh Refinement (AMR)
+Considering the above limitations of DVM, we adopt Adaptive Mesh Refinement (AMR) to improve solving efficiency. AMR reallocates computational resources based on flow features, balancing efficiency and accuracy. 
+
+<div align="center">
+    <img src="https://i.postimg.cc/FzQx9cSD/PV-AMR.png" alt="Image 1" title="Simultaneous AMR" width="700" style="margin-right: 1px;">
+</div>
+
+<div align="center">
+    <img src="https://i.postimg.cc/hvp8Z2rk/ps.gif" alt="Image 1" title="physical space" width="300" style="margin-right: 1px;">
+    <img src="https://i.postimg.cc/rwVV31C8/d5ubx-sf3vl.gif" alt="Image 2" title="simultaneously, in velocity space" width = "300" style="margin-left: 1px;">
+</div>
+
+<!-- Generally speaking, the criteria for refining or coarsening physical space grids are based on the magnitude of macroscopic variable gradients. This approach aims to use finer grids to capture more pronounced variations, thereby modeling more complex flow fields. On the other hand, criteria for velocity space grids are based on the proportion of energy of particles within the grid relative to the total energy of all particles at that physical point. This is because particles with higher energy proportions exert a more significant influence on the macroscopic behavior of the flow, necessitating more accurate numerical integration.  -->
+The current solver discretizes with tree-based Cartesian grids, adapting the mesh in both physical and velocity space simultaneously, according to characteristics including velocity field gradients, boundary complexity, and energy proportion at phase points. Users can also customize the criteria for grid refinement or coarsening according to their own needs.
+
+### Cartesian Grids
+Cartesian grids offer advantages such as high grid quality and automation but struggle with complex boundary shapes. In the subsequent development of the solver, we will introduce Immersed Boundary Method (IBM) to address this challenge.
+
+## Installation
+The current solver does not yet provide a user interface. If you'd like to try it out, we recommend cloning the code repository from Git to your local machine and then using the Julia package manager to install its dependencies in the directory where the code resides:
+```julia
+julia> ]
+(v1.10) pkg> activate .
+(//your_path//) pkg> instantiate
+```
+Soon, we will provide a complete user interface with customizable boundaries, and we will register it in the [Julia package registry](https://github.com/JuliaRegistries/General). Stay tuned!
+## Current Usage
+Currently, flow properties, initial contitions and boundary conditions are defined in `Global_Data` in `/src/types.jl`. It should be noted specifically that the boundary follows the sequence of -x, +x, -y, +y. The maximum level of the refinement is defined in `/src/abstract.jl` by `DVM_PS_MAXLEVEL` and `DVM_VS_MAXLEVEL` respectively. Grid refinement and coarsening behaviors are defined in `adaptive.jl` and `vs_adaptive.jl`. Feel free to modify them according to your needs.

src/DVM.jl → deprecated/DVM.jl

@@ -1,12 +1,15 @@
-using Pkg
-Pkg.activate(@__DIR__)
+#using Pkg
+#Pkg.activate(@__DIR__)
 using MPI
 using P4est
-using P4estTypes
+#using P4estTypes
+include("../lib/P4estTypes/src/P4estTypes.jl")
+using .P4estTypes
+
 using StaticArrays
 using Parameters
 using SpecialFunctions
-using Makie
+#using Makie
 using CairoMakie
 using PythonCall
 scipy = pyimport("scipy")

example/cavity_sim.jl

@@ -0,0 +1,90 @@
+"""
+To execute the code (assuming 16 cores are going to be used):
+- from the main directory: mpiexecjl -n 16 --project=. julia example/cavity_sim.jl
+- from the file directory: mpiexecjl -n 16 --project=.. julia cavity_sim.jl
+"""
+
+using KitAMR, MPI, CairoMakie
+
+MPI.Init()
+ps4est, DVM_data = KitAMR.init()
+if MPI.Comm_rank(MPI.COMM_WORLD) == 0
+    KitAMR.save_set(DVM_data)
+end
+KitAMR.save_result(DVM_data)
+while true
+    if DVM_data.global_data.adapt_step == 10
+        KitAMR.update_slope!(DVM_data)
+        KitAMR.update_gradmax!(DVM_data)
+        KitAMR.PS_refine!(ps4est)
+        KitAMR.PS_coarsen!(ps4est)
+        KitAMR.PS_balance!(ps4est)
+        if DVM_data.global_data.partition_step == 30
+            KitAMR.PS_partition!(ps4est, DVM_data)
+            DVM_data.global_data.partition_step = 0
+        end
+        if DVM_data.global_data.vs_adapt_step == 50
+            KitAMR.vs_refine!(DVM_data)
+            KitAMR.vs_coarsen!(DVM_data)
+            DVM_data.global_data.vs_adapt_step = 0
+        end
+        KitAMR.update_ghost!(ps4est, DVM_data)
+        KitAMR.update_neighbor!(ps4est, DVM_data)
+        KitAMR.update_faces!(ps4est, DVM_data)
+        DVM_data.global_data.adapt_step = 0
+    end
+    KitAMR.update_Δt!(DVM_data) # 798.125 μs (41515 allocations: 806.91 KiB)
+    KitAMR.update_slope!(DVM_data) # 1.304 s (22645678 allocations: 1.43 GiB) -> final: 336.005 ms (32818 allocations: 1.35 MiB)
+    KitAMR.DVM_data_exchange!(ps4est, DVM_data, Val(2)) # 0.121942
+    KitAMR.update_flux!(DVM_data) # 4.437 s (174354747 allocations: 5.30 GiB)->3.268 s (75782427 allocations: 4.84 GiB) -> final: 540.574 ms (1301847 allocations: 1.42 GiB)->411.685 ms (1262009 allocations: 1.51 GiB)
+    KitAMR.update_volume!(DVM_data) # 1.262 s (59874178 allocations: 1.55 GiB) ->469.264 ms (9070978 allocations: 623.28 MiB)-> final 50.184 ms (248053 allocations: 158.25 MiB)
+    KitAMR.DVM_data_exchange!(ps4est, DVM_data, Val(1))
+    DVM_data.global_data.adapt_step += 1
+    DVM_data.global_data.vs_adapt_step += 1
+    DVM_data.global_data.partition_step += 1
+    KitAMR.save_result(DVM_data)
+    if MPI.Comm_rank(MPI.COMM_WORLD) == 0
+        @show DVM_data.global_data.gas.sim_time
+        pp = PointerWrapper(ps4est)
+        @show pp.local_num_quadrants[]
+    end
+    DVM_data.global_data.gas.sim_time > DVM_data.data_set.end_time && break
+    MPI.Barrier(MPI.COMM_WORLD)
+    if MPI.Comm_rank(MPI.COMM_WORLD) == MPI.Comm_size(MPI.COMM_WORLD) - 1
+        ps_data = DVM_data.trees.data[end][end]
+        @show ps_data.vs_data.vs_num
+    end
+end
+KitAMR.save_VS_final(DVM_data)
+meshes = KitAMR.collect_mesh(DVM_data)
+solutions = KitAMR.collect_solution(ps4est, DVM_data)
+
+if MPI.Comm_rank(MPI.COMM_WORLD) == 0
+    x, y = KitAMR.reshape_mesh(meshes)
+    f = Figure()
+    ax = Axis(f[1, 1])
+    KitAMR.mesh_plot(x, y, ax)
+    save("mesh.png", f)
+
+    x, y, variable = KitAMR.reshape_solutions(solutions, DVM_data.global_data, :prim, 4)
+    variable = 1 ./ variable
+    f = Figure()
+    ax = Axis(f[1, 1])
+    co = contourf!(x, y, variable)
+    Colorbar(f[1, 2], co)
+    save("test_vs_adaptive.png", f)
+end
+
+KitAMR.finalize!(ps4est, DVM_data)
+
+phase_cells = 0
+for i in eachindex(DVM_data.trees.data)
+    for j in eachindex(DVM_data.trees.data[i])
+        global phase_cells
+        phase_cells += DVM_data.trees.data[i][j].vs_data.vs_num
+    end
+end
+phase_cells = MPI.Allreduce(phase_cells, +, MPI.COMM_WORLD)
+if MPI.Comm_rank(MPI.COMM_WORLD) == 0
+    @show phase_cells
+end